[jvinuales]

Hi everyone,

I am OpenFOAM newbie, so bear with me!

I am a bit confused about cyclicAMI and how to set it up... I have not been able to find some tutorials that replicate what I am trying to do.

So, I am studying 3D periodic porous structures. Previously, I had access to commercial software (STAR CCM+) and I will have a set up with a single unit cell with periodic conditions in all 3 directions (in the flow direction, a special type of periodic BC which is called "fully developed flow" that allows specifying certain mass flow rate through the boundary). This set up keeps the domain (and thus cell count) very low, which is important as I have to characterize many different geometries (100-1000) every time.

The objective is to replicate this set up in OpenFOAM (if possible). I was thinking cyclicAMI may be the best option, as creating periodic meshes can be challenging for complex geometries (as my case). I will be using simpleFOAM at the moment (with a modified version to account for heat transfer later on)

- Is it cyclicAMI the right choice?
- For the special case of the periodicity in the flow direction, can cyclicAMI be used? how can I specify certain flow through the boundary?
- Once I created the mesh, I just have to modify the following files:
- constant/polyMesh/boundary
- 0/U , 0/p , 0/T ...
Is this workflow correct? what will be correct set up in each file?

Thanks!

Javier

[ybapat]

You can use cyclicAMI. But I think translation periodic will need additional implementation. It can be either pressure gradient specification or mass flow rate specification. You can use it for rotational periodic directly.

[.bastian]

Originally Posted by jvinuales

[...]
- Is it cyclicAMI the right choice?
- For the special case of the periodicity in the flow direction, can cyclicAMI be used? how can I specify certain flow through the boundary?
[...]

Hello Javier.

Lets begin with some of your questions.
Generally there are two main options to set periodicity.
If you are able to generate a conformal mesh, cyclic would be the choice: User Guide: Cyclic
In case of non conformal meshes you may choose CyclicAMI oder GGI (foam-extend): User Guide: cyclicAMI
Conformal mesh means, that the nodes of the faces you want to connect are matching together perfectly.

As far i know, you can't specify a flow through an cyclic/cyclicAMI interface. Its only a connection of both sides.

It would be nice if you can give us a few more informations like:
- which OpenFOAM version do you want to use?
- how do you generate the mesh?

Greetings
bastian

[jvinuales]

Thank you for your answers.

I thought I would go about this step by step:
- first, will try to set up cyclicAMI in the 2 directions perpendicular to the flow.
- second, I will try to implement periodicity in the flow direction.

The first case should be more straight-forward. I managed to set it up in 1 direction, and it ran without problems. However, the results seemed odd as not as much flow as expected was getting through the cyclicAMI boundary. Furthermore, if I added the second pair it will not run (file attached).

Regarding the second case, I saw that some flow rate can be forced through the boundary with fvOptions (momentumsource:meanVelocityForce), and hopefully works well with cyclicAMI. For the pressure, I saw jumpCyclicAMI, that could work for representing a jump in pressure from inlet to outlet. However, I am not sure about other fields such temperature and turbulence variables.

Javier

[jvinuales]

Originally Posted by .bastian

Hello Javier.


It would be nice if you can give us a few more informations like:
- which OpenFOAM version do you want to use?
- how do you generate the mesh?

I am running openFOAM in a windows machine. I am currently working with 2 versions (have not decided yet what works better for me): BlueCFD with OpenFOAM 5.x and v1912 with an Ubuntu Bash.

In terms of meshing, currently using cf-mesh+ (with a trial)... also still trying to figure out what works better for me.

[jvinuales]

I made some progress with this case.

I have the periodicity in all three directions, and the case run without problems.

I used cyclicAMI in all the directions, and add the following to fvOptions:

momentumSource
{
    type            meanVelocityForce;

    selectionMode   all;

    fields          (U);
    Ubar            (-0.6 0 0);
}

pRefCell	0;
pRefValue	0;

The problem is that the results are not very physical or realistic, as there is no pressure drop along the flow (obviously, as the pressure is mapped from one face to another).

Is there any boundary condition I can apply for the pressure that results in a physically realistic condition? Like an undefined constant offset?? or any other set up recommendations?

case:
https://drive.google.com/file/d/1R85...ew?usp=sharing

gyroid_p.JPG
gyroid_U.JPG
gyroid_U_zx.JPG

[jvinuales]

I found the fixedJumpAMI boundary condition that could be applied in the flow direction:

0/p

    x1
	{
		type            fixedJumpAMI;
		patchType       cyclicAMI;
		jump            uniform 1;
		value           $internalField;
	}
    x2
    {
		type            fixedJumpAMI;
		patchType       cyclicAMI;
		value           $internalField;
    }

0/U

    x1
    {

        type		cyclicAMI;

    }
    x2
    {

        type		cyclicAMI;

    }

This, in theory, achieves the expected behaviour:

Gyroid_p.JPG
Gyroid_U.JPG

I will start trying to set up turbulence and temperature. Regarding turbulence, can I assume that the turbulence parameters are periodic (as U)? hence I can use cyclic?

[pete20r2]

I've done a few simulations almost exactly like this meanvelocityForce is a closed loop control method that applies a pressure gradient to the flow to sustain the desired flow rate. If you watch the logs there should be a printout of the current pressure drop the solver is using.

I've found this to be the most reliable and convergent method. Especially if you want to study the flow at particular flow rates rather than applying pressure drops and hoping to find a flow rate.

I used cyclic turbulence parameters and U in my study and p ends up being cyclic with a jump. T can be difficult because if it is cyclic it will blow up. Additionally, if the flow is compressible then the physical properties are dependant on T so a cyclic study is non-physical.

The other way I have achieved good results is using mappedPatch.
Inlet is set to mappedPatch type and condition is mapped.
Turbulence is mapped from outlet to inlet, outlets are zeroGradient.
U is mapped the same but with setAverage true; as is T.
p is the only field that does not get mapped. outlet uses the fixedMean condition and inlet uses zeroGradient.
Start the flow overdamped and then lower the damping after it's starred and you should see convergence.
You may also need to extend your domain to be 2,3, or 4 repetitions of the geometry.

[jvinuales]

Thanks for your answer! It really helps.

Originally Posted by pete20r2

I've found this to be the most reliable and convergent method. Especially if you want to study the flow at particular flow rates rather than applying pressure drops and hoping to find a flow rate.

Yes, I found convergence is much better with meanVelocityForce than using a jump in pressure. As you point out, I also need to control the flow rate through the domain rather than the pressure jump, so it seems more logical as well. However, I do not understand what boundary condition should apply in p as the cyclicAMi will result in no gradient in pressure, meaning that the results are unphysical... am I missing something?

Originally Posted by pete20r2

I used cyclic turbulence parameters and U in my study and p ends up being cyclic with a jump.

do you use cyclicAMI in 0/p? If I do so I get no gradient/jump

Originally Posted by pete20r2

T can be difficult because if it is cyclic it will blow up. Additionally, if the flow is compressible then the physical properties are dependant on T so a cyclic study is non-physical.

I am using incompressible fluid, so I am not concern about compressibility. Also, I can live with non-temperature-dependant properties for this study. I wonder if a fixedJump in T, together with a fixed gradient in the wall would result in a const heat flux conditions... I need to investigate more.

Originally Posted by pete20r2

The other way I have achieved good results is using mappedPatch.

I will give a try to this... do you still use meanVelocityForce to keep the flow going? or just the initial conditions?

Originally Posted by pete20r2

Start the flow overdamped and then lower the damping after it's starred and you should see convergence.

I am not sure I fully understand what you mean with this statement?

[pete20r2]

However, I do not understand what boundary condition should apply in p as the cyclicAMi will result in no gradient in pressure, meaning that the results are unphysical... am I missing something?

meanVelocityForce will apply a gradient to the entire p field.
The meanVelocityForce fvOption will cause these outputs:

Start of execution

Creating finite volume options from "constant/fvOptions"

Selecting finite volume options type meanVelocityForce
    Source: momentumSource
    - selecting all cells
    - selected 60000 cell(s) with volume 16
    Initial pressure gradient = 0

Per iteration

Time = 0.8

PIMPLE: iteration 1
smoothSolver:  Solving for Ux, Initial residual = 0.0112044, Final residual = 5.90986e-06, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 0.0606708, Final residual = 3.37513e-06, No Iterations 3
smoothSolver:  Solving for Uz, Initial residual = 0.0579524, Final residual = 2.73287e-06, No Iterations 3
Pressure gradient source: uncorrected Ubar = 0.1335, pressure gradient = 2.67342e-05
GAMG:  Solving for p, Initial residual = 0.13367, Final residual = 0.00203606, No Iterations 2
time step continuity errors : sum local = 1.22471e-06, global = 1.74992e-20, cumulative = -2.29145e-19
Pressure gradient source: uncorrected Ubar = 0.1335, pressure gradient = 2.6722e-05
GAMG:  Solving for p, Initial residual = 0.0107536, Final residual = 7.18508e-07, No Iterations 5
time step continuity errors : sum local = 4.79701e-10, global = -5.57251e-20, cumulative = -2.84871e-19
Pressure gradient source: uncorrected Ubar = 0.1335, pressure gradient = 2.66637e-05
ExecutionTime = 0.52 s  ClockTime = 2 s

fieldAverage fieldAverage1 write:
    Calculating averages

Courant Number mean: 0.302291 max: 0.522438

This is a kinematic pressure because we are in an incompressible solver. If you take that gradient and multiply by the length of your channel, you will get the kinematic pressure drop between inlet and outlet.
Please look at the channel395 tutorial in "openfoam/tutorials/incompressible/pimpleFoam/LES/channel395". It works just fine for RAS cases, but that tutorial just happens to be where the function is demonstrated.

T can be added as a passive scalar using a function in controlDict called scalarTransport. I have a rough example but I'm not sure it works properly.

Or you can duplicate the case and create a new case for the scalarTransportFoam solver to run after your original case on the frozen field.

About mapped patches:

I will give a try to this... do you still use meanVelocityForce to keep the flow going? or just the initial conditions?

No, mappedPatch provides for mapping an outlet to an inlet and ensuring a specified flow rate with the setAverage option.

[jvinuales]

ok, I finally understand!! Thanks!!

My problem was that the pressure field is not updated... I was just checking the results with Paraview and with functions as pressureDifferencePatch and these do not reflect the pressure gradient.

Is there a way to update the p field with the pressure gradient source so the results in Paraview and the postprocessing functions reflect the actual p?

[pete20r2]

Sure you can, but not automatically without using more functions.
I do this by using the calculator function in paraview.

Hit ctrl+space to bring up the filter search and look for "calculator".
Grab the pressure gradient from the last iteration in your log file.
Multiply the pressure gradient by the point position in the direction of your flow, for me that is the X coordinate. The gradient may need to be negative to be in the correct direction but that should be obvious.

You now have a new field called realP that includes the pressure gradient.
While you're in there you could also multiply by density to calculate the real static pressure instead of kinematic pressure.

[jvinuales]

Thanks again.

Originally Posted by pete20r2

T can be added as a passive scalar using a function in controlDict called scalarTransport. I have a rough example but I'm not sure it works properly.

Or you can duplicate the case and create a new case for the scalarTransportFoam solver to run after your original case on the frozen field.

I figured out that in terms of temperature, for a constant heat flux case, there is an analogy with pressure i.e. there is a constant jump in the flow for subsequent repeatable structures. (source: Lattice Boltzmann simulations of heat transfer in fully developed periodic incompressible flows, Zimeng Wang et al).

If that is the case, jumpAMI should work correctly and lay physical results, always making sure that gradT and dT (jump) are consistent with each other for the fluid properties. However, the temperature increases without control, as there is no temperature reference or any ways to specify the inlet temperature (for example).

In the case of the pressure this is avoided through pRefCell and pRefValue... is there any equivalent formulation for temperature?

Also, I saw the possibility to try to force the inlet average temperature to be maintained as here:

Maintaining average temperature at Inlet in cyclic Boundary condition

did you have any success implementing this option?

[bany]

Hi, peter, i am interested in your method. But when using these in a pipe to get a fully developed pipe flow, i am always failed.
For your mapped method, did you add a initial perturbation to the flow field such as turbulentInlet or other? Because i found if there is no initial perturbation , it is difficult to produce a fully developed pipe flow.
And for my case with a pipe of 20D length, U is mapped with setAverage ture (110 0 0). k is also mapped with setAverage no. p and nut is not mapped. After about 100 flow times, i cannot see a fully developed flow structure and the Umean is gradually increasing.
By the way, i found add a wallFunction to nut only change the U profile near the wall having few effect on the entire U proflie.

[pete20r2]

Hi bany,
Sorry for the delayed reply. What are your pressure conditions?
Also every quantity except for pressure should be mapped, this includes the turbulent fields.
I have not used many LES simulations but from my understanding LES is much more sensitive to initial conditions than RAS. Obviously in a periodic flow we should see convergence to the same solution over a wide range of initial conditions. I think you are right to try and introduce some initial perturbations.

[bany]

Originally Posted by pete20r2

Hi bany,
Sorry for the delayed reply. What are your pressure conditions?
Also every quantity except for pressure should be mapped, this includes the turbulent fields.
I have not used many LES simulations but from my understanding LES is much more sensitive to initial conditions than RAS. Obviously in a periodic flow we should see convergence to the same solution over a wide range of initial conditions. I think you are right to try and introduce some initial perturbations.

Hi, thanks for your replying. The following is my setup:

[pete20r2]

Your zip file is broken, do you want to try again?

[bany]

Originally Posted by pete20r2

Your zip file is broken, do you want to try again?

OK.
The y+~5, and the solver pisoFoam is used.
Much thanks!

[pete20r2]

Is this not what you are expecting?


I changed a few things:
Pipe made shorter
Mapping taken from exit patch instead of offset
Mapped nut also
p field uses fixedFluxPressure (probably has no effect in this context)
Have a look at how blockMeshDict is defined for the patches
Added minIter 1 to U and p solutions

I had to turn some of your functions off to get it to run but shouldn't effect solution. See attached

pipe-mapped-short2.zip

[bany]

Originally Posted by pete20r2

Is this not what you are expecting?

I changed a few things:
Pipe made shorter
Mapping taken from exit patch instead of offset
Mapped nut also
p field uses fixedFluxPressure (probably has no effect in this context)
Have a look at how blockMeshDict is defined for the patches
Added minIter 1 to U and p solutions

I had to turn some of your functions off to get it to run but shouldn't effect solution. See attached

What i am expecting is that a fully developed turbulent circular pipe flow with mapped BC. So the wall functions are used. Now, according to your results, there is no turbulence and the U profile is also parabolic.However, i also cannot get a correct turbulent flow in the pipe.