[SylvainWang]

hello, foamers,
I am now running simulation of a fluidic osciillator to predict its performance, especially its oscillation frequency. Before, I utilized FLUENT without any problem. I have turned to OF since about 5 months ago because of the parallel license problem of Fluent.
The fluidic oscillator is working after the coanda effect, and has one inlet, two outlets and two feedback loops. The direction of the main jet from the inlet will switch from one outlet to the other automatically. the velocity in the inlet throat would always be sonic.

I am now doing the numerical sensitivity study concerning my case. The solver I chose is sonicFoam since the flow inside the oscillator is transonic. I have tried different ddt scheme, Maximum Courant number, convection scheme with kEpsilon family turbulent model without any problem.
However, when I turned to the kOmegaSST turbulent model, I can never reach a converged results anymore, even though I have tried different discretization scheme(1st order, 2nd orde upwind, TVD), decrease the time step and dense the mesh(average Y+ 5.7, max Y+ 27, and average y+ 1.6, Max y+ 20 (not precise)).

the case setups are attached in the following link, including all the numerical schemes and the initial fields obtained from KEpsilon model.
https://www.dropbox.com/s/bflljy4b7e...or.tar.gz?dl=0

Hope someone can help. Thanks in advance.
Sylvain.

[RodriguezFatz]

Sylvain, in your /0 directory I see you are using high-Re boundary conditions and you write about y+ average 5.7 or 1.6. This means low-Re mesh. Are you aware of this discrepancy?
The fvSchemes looks pretty safe.
A log-file would be interesting.

[SylvainWang]

Originally Posted by RodriguezFatz

Sylvain, in your /0 directory I see you are using high-Re boundary conditions and you write about y+ average 5.7 or 1.6. This means low-Re mesh. Are you aware of this discrepancy?
The fvSchemes looks pretty safe.
A log-file would be interesting.

hi, Rodri,

Thanks very much for you reply. In fact, to the best of my knowledge, kOmegaSST is a blended kOmega and kEpsilon turbulent model, which means a blended high Re and Low Re turbulent model. And in OmegaWallFunction or EpsilonWallFunction, the value of turbulent factors are calculated depend on the Y+ according to the law of wall. (in the code, the separation value between viscous sublayer and log layer is about 11). However, in general, if the y+ is smaller than 5 which means the first grid lies in the viscous sublayer, the results should be better. That's why I'am always working on refine the mesh to solve this problem. I've never thought that it is caused by the wall condition modeling choice since with the kEpsilon turbulent model, it works fine.


https://www.dropbox.com/s/qzqhqpaelt...sh.tar.gz?dl=0
follow this link , you can find a new finer mesh file and schèmes setup, and a running log file. In this case, I've just chosen the simplest scheme to see if it works. unfortunately, it didn't work.

To tell the truth, I am not clear about the low-Re wall condition. I'll try it and let you know the results.

Sylvain

[RodriguezFatz]

(basic) k-omega and k-epsilon are both high-Re models. Only with some extensions close to the wall they are applicable for low-Re meshes. In OpenFoam these extensions are implemented by your boundary conditions. So if you don't use low-Re boundary conditions the models might even not converge, or they might converge (such as your k-epsilon run) but probably give bad results. You can not run a case with high-Re boundary conditions and low-Re mesh and expect convergence or good results.
I would try a coarser wall mesh and see if it runs.

Edit: You say "it didn't work"... What didn't work?

[SylvainWang]

Originally Posted by RodriguezFatz

(basic) k-omega and k-epsilon are both high-Re models. Only with some extensions close to the wall they are applicable for low-Re meshes. In OpenFoam these extensions are implemented by your boundary conditions. So if you don't use low-Re boundary conditions the models might even not converge, or they might converge (such as your k-epsilon run) but probably give bad results. You can not run a case with high-Re boundary conditions and low-Re mesh and expect convergence or good results.
I would try a coarser wall mesh and see if it runs.

Edit: You say "it didn't work"... What didn't work?

I just said that just with finer mesh, it didn't work. I have modified the wall condition of k to kLowReWallFunction, and the others kept the same since alphat and mut are calculated depend on k, and omegaWallFunction has already combined the Low AND High Re effects. Are my words right?
========
I just got the results of new setting: didn't work....
here are the setting with coarser mesh and log file.
https://www.dropbox.com/s/66zh5uk6bu...sh.tar.gz?dl=0

[RodriguezFatz]

I just know it for the incompressible case, but there you actually need a b.c. for nut for low-Re meshes that implements some damping for nut.
My question was: If you write "it didn't work" what do you mean by that. What did not work? Did the solution diverge or what?

[SylvainWang]

Originally Posted by RodriguezFatz

I just know it for the incompressible case, but there you actually need a b.c. for nut for low-Re meshes that implements some damping for nut.
My question was: If you write "it didn't work" what do you mean by that. What did not work? Did the solution diverge or what?

ok, it means the solution diverges, just crashed

[RodriguezFatz]

1) You changed the numerical schemes between the two cases. Why did you do that?
2) What happend to the boundary conditions? In the 0.0075 directory I see some values lists.

[SylvainWang]

Originally Posted by RodriguezFatz

1) You changed the numerical schemes between the two cases. Why did you do that?
2) What happend to the boundary conditions? In the 0.0075 directory I see some values lists.

1) In pursuing a converged solution, not a precise one, I think a first order scheme is easier.
2) the 0.0075 directory is the initial results obtained from the kEpsilon model, I just calculate the omega by epsilon/k/0.09. There are only three boundary in each variety, inlet, outlet and wall_all. You can locate them easily.

The 2nd order scheme indeed has some effect, however, with my previous setup, it just postpone the time of crash. Nevertheless, I'll try the 2nd order upwind scheme with kLowReWallFunction

[RodriguezFatz]

1) But in the new case you use limitedLinear and not upwind. I would use upwind if I have trouble with convergence. You should also use "Gauss linear uncorrected" as a lapacianScheme instead of the one you use right now.
2) See your "wall_all" b.c. for "k". I think your intention is to use kqrwallfunction which means a zero gradient b.c. for k. In your case you set k in the 0.0075 directory at the "wall_all" boundary to fixed values. This is not the same.

[SylvainWang]

Originally Posted by RodriguezFatz

1) But in the new case you use limitedLinear and not upwind. I would use upwind if I have trouble with convergence. You should also use "Gauss linear uncorrected" as a lapacianScheme instead of the one you use right now.
2) See your "wall_all" b.c. for "k". I think your intention is to use kqrwallfunction which means a zero gradient b.c. for k. In your case you set k in the 0.0075 directory at the "wall_all" boundary to fixed values. This is not the same.

1)yes, i prefer the limitedLinear since in my mind, this TVD scheme is just like a blend of upwind and central linear, with limited number 1.0 and 0, so I appoint the limited number to 1.0 to stand for the upwind scheme. (correct me if I am wrong). I will try the Gauss linear uncorrected lapacianScheme.
2)in this case: oscillatorNormalMesh, the wall condition of k is kLowReWallFunction as you recommend. I just look up the code of it. Indeed, it is a fixedValue B.C depend on the y+ value. While the kqRWallFunction is a zeroGradient B.C. Which one would be better in my case? I prefer the kLowReWallFunction.
3) in the code of kLowReWallFuntion, the calculation of y+ value which used to separate the viscous and log region is shown below. But I don't understand why there needs a loop, just offer the constant should not be better?

scalar kLowReWallFunctionFvPatchScalarField::yPlusLam
59 (
60 const scalar kappa,
61 const scalar E
62 )
63 {
64 scalar ypl = 11.0;
65
66 for (int i=0; i<10; i++)
67 {
68 ypl = log(max(E*ypl, 1))/kappa;
69 }
70
71 return ypl;
72 }

[SylvainWang]

Originally Posted by RodriguezFatz

1) But in the new case you use limitedLinear and not upwind. I would use upwind if I have trouble with convergence. You should also use "Gauss linear uncorrected" as a lapacianScheme instead of the one you use right now.
2) See your "wall_all" b.c. for "k". I think your intention is to use kqrwallfunction which means a zero gradient b.c. for k. In your case you set k in the 0.0075 directory at the "wall_all" boundary to fixed values. This is not the same.

I am trying to modify the inlet and outlet b.c of k, omega, U , P. In fact, I am wondering how to decide the parameters's value of reversed flow in the inletOutlet condition. How to calculate the mixingLength?(in my case, 5%of the inlet width is OK?)

In my case, the inlet width is 0.01m=10mm, with total pressure inlet condition 2.5e5 Pa. the outlet is a region directly open to the atmosphere with the totoal pressure outlet condtion 1e5 Pa.
p:

outlet
{
type totalPressure;
rho rho;
psi none;
gamma 1.4;
p0 uniform 100000;
value uniform 1e5;
}
wall_all
{
type zeroGradient;
}
inlet
{
type totalPressure;
rho rho;
psi none;
gamma 1.4;
p0 uniform 250000;
value uniform 2.5e5;
}

U:

outlet
{
type inletOutlet; //zeroGradient;
inletValue $internalField;
value $internalField;
}
wall_all
{
type fixedValue;
value uniform (0 0 0);
}
inlet
{
type pressureInletVelocity;
value uniform (0 10 0);
}

k:

outlet
{
type inletOutlet;
inletValue uniform 0.25;
value $internalField; }
inlet
{
type turbulentIntensityKineticEnergyInlet;
intensity 0.05;
value uniform 1;
}
wall_all
{
type compressible::kLowReWallFunction;
value uniform 1;
}

omega:

outlet
{
type inletOutlet;
inletValue uniform 600 ;
value $internalField;
}
inlet
{
type compressible::turbulentMixingLengthFrequencyInlet;
mixingLength 0.0005;
value uniform 3000;
}
wall_all
{
type compressible:megaWallFunction;
value uniform 5000;
}

[RodriguezFatz]

Ok, I didn't know that kLowReWallFunction works for both high-Re and low-Re meshes. Inappropriate name sorry. Probably your b.c. are ok.
If you look at the code of limitedlinear, you find
max(min(twoByk_*r, 1), 0)
twoByk_ = 2/k and k is the 1.0 in your case.
So your limiter function is
max(min(r, 1), 0)
This is not the same as upwind as I understand it.

Edit: But again, I would try the uncorrected laplacian. For me it looks like all the crap in your log starts with the k/omega bounding errors. In my experience they often stem from the laplacian.

[SylvainWang]

Originally Posted by RodriguezFatz

Ok, I didn't know that kLowReWallFunction works for both high-Re and low-Re meshes. Inappropriate name sorry. Probably your b.c. are ok.
If you look at the code of limitedlinear, you find
max(min(twoByk_*r, 1), 0)
twoByk_ = 2/k and k is the 1.0 in your case.
So your limiter function is
max(min(r, 1), 0)
This is not the same as upwind as I understand it.

Edit: But again, I would try the uncorrected laplacian. For me it looks like all the crap in your log starts with the k/omega bounding errors. In my experience they often stem from the laplacian.

I didn't know either before. Never mind.
Yes, you r right. they are not the same. I just feel they are similar. Well, TVD scheme seems alway more complex, but if you treat them as upwind or linear or similar, they will become quite simple, right?
I'll try both the kqr and klowRe.
Until now, I havn't got the converged results. This problem has puzzled me more than two weeks. It crashed quickly. How do you think about the value I set of the inlet and outlet boundary conditions, k, omega, p and U?

[RodriguezFatz]

I don't run compressible cases, so I am the wrong one to judge the b.c.

[SylvainWang]

Originally Posted by RodriguezFatz

I don't run compressible cases, so I am the wrong one to judge the b.c.

good news! I have got the converged results with new boundary condition. I just modified the inletOutlet condition to zeroGradient condition with Euler time scheme and upwind omega advection scheme (the others 2nd order upwind).
However, this leads to another question about the inletOutlet b.c. In which cases we need to asign a fixed value to reversed flow which meets the inletOutlet b.c. ? I've seen that in many cases in the OF tutorial the outlet conditions are inletOutlet. But I just can not understand why.
In my case, the outlet in fact is the open space in atmosphere. It indeed has reverse flow.