Dear Prof Spencer and dear Prof Moxey,
Thank you very much for your kind and prompt reply. I am glad to know about this implementation. I successfully compiled the development branch of the code on Cray-XC40 and tried it for the example test case ('ChanFlow_m8_Flowrate.xml') and then channel flow-2D and it works.
As soon as I try to use it for case-3DH1D, the simulation starts but does not proceed, it stuck just after writing the initial file at 0th timestep, I even ran the simulation for ~40 CPU-hrs. The simulation proceeds only when processors are less than or equal to 2 but after that (>2), it stuck. I tried to debug it and found out that it stuck at "MeasureFlowrate" function in ./solvers/IncNavierStokesSolver/EquationSystemsVelocityCorrectionScheme.cpp specifically at:
m_comm->AllReduce(flowrate, LibUtilities::ReduceSum);
however, this implementation looks correct, but I don't know whats wrong here.
In addition, the 3DH1D simulation with <=2 processors crashes, therefore, I checked the outputs of few variables and found out that "m_flowrateArea" is wrong in this case. The Lx*Ly*Lz= 4Pi*2*1.333Pi, therefore "m_flowrateArea" for the 3D case, it should be 8.37, but code's output is 128 while for 2D case, its correct (i.e. 2). This, "m_flowrateArea"=128 is the result of Ly*HomModesZ (should be Ly*Lz?). Moreover, the "flowrate" from the function turns out to be negative for the 3DH1D case during the very first step and I guess it also causes the crashing. While it is positive and remains constant in 2D case during the entire run.
Do you have any suggestions which I can try and also please correct me if I am doing or understanding it in a wrong way.
Snippet of parameter.xml file:
<!-- xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx -->
<CONDITIONS>
<SOLVERINFO>
<I PROPERTY="SolverType" VALUE="VelocityCorrectionScheme"/>
<I PROPERTY="EQTYPE" VALUE="UnsteadyNavierStokes"/>
<I PROPERTY="AdvectionForm" VALUE="Convective"/>
<I PROPERTY="Projection" VALUE="Galerkin"/>
<I PROPERTY="TimeIntegrationMethod" VALUE="IMEXOrder1"/>
<I PROPERTY="HOMOGENEOUS" VALUE="1D"/>
<I PROPERTY="DEALIASING" VALUE="ON"/>
<I PROPERTY="USEFFT" VALUE="FFTW" />
<I PROPERTY="SPECTRALHPDEALIASING" VALUE="True" />
</SOLVERINFO>
<PARAMETERS>
<P> TimeStep = 0.00001 </P>
<P> NumSteps = 1000000 </P>
<P> IO_CheckSteps = 500 </P>
<P> IO_InfoSteps = 1 </P>
<P> Kinvis = 1.0/180.0 </P>
<P> Flowrate = 15.0 </P> <!--For Re_tau=180 and Re_b=2700 (based on H) -->
<P> IO_FlowSteps = 1 </P>
<P> HomModesZ = 64 </P>
<P> LZ = 4.0*PI/3.0 </P>
</PARAMETERS>
<VARIABLES>
<V ID="0"> u </V>
<V ID="1"> v </V>
<V ID="2"> w </V>
<V ID="3"> p </V>
</VARIABLES>
<BOUNDARYREGIONS>
<B ID="0"> C[1] </B>
<B ID="1"> C[2] </B>
<B ID="2"> C[3] </B>
</BOUNDARYREGIONS>
<BOUNDARYCONDITIONS>
<REGION REF="0">
<D VAR="u" VALUE="0" />
<D VAR="v" VALUE="0" />
<D VAR="w" VALUE="0" />
<N VAR="p" USERDEFINEDTYPE="H" VALUE="0" />
</REGION>
<REGION REF="1">
<P VAR="u" VALUE="[2]" />
<P VAR="v" VALUE="[2]" />
<P VAR="w" VALUE="[2]" />
<P VAR="p" VALUE="[2]" />
</REGION>
<REGION REF="2">
<P VAR="u" VALUE="[1]" USERDEFINEDTYPE="Flowrate" /> <!--Cross section for Q = outlet. -->
<P VAR="v" VALUE="[1]" />
<P VAR="w" VALUE="[1]" />
<P VAR="p" VALUE="[1]" />
</REGION>
</BOUNDARYCONDITIONS>
<FUNCTION NAME="FlowrateForce">
<E VAR="ForceX" VALUE="1" /> <!--Flow is in the X-Direction, so explicitly defined -->
<E VAR="ForceY" VALUE="0" />
<E VAR="ForceZ" VALUE="0" />
</FUNCTION>
<FUNCTION NAME="InitialConditions">
<F VAR="u" FILE="TurbChFl_3DH1D_83.rst" /> <!--Initialize with fully developed turbulent flow -->
<F VAR="v" FILE="TurbChFl_3DH1D_83.rst" />
<F VAR="w" FILE="TurbChFl_3DH1D_83.rst" />
<F VAR="p" FILE="TurbChFl_3DH1D_83.rst" />
</FUNCTION>
</CONDITIONS>
<!-- xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx -->
Thanks you very much!
Sandeep