[maike151]

Hello everyone,


I know there are several threads about this topic, but none of them helped me.


I am trying to simulate a rather simple air flow with two inlets with two different temperatures and I want to see what the temperature field at the outlet looks like.


In the second timestep I get the following error:


[6]
[6]
[6] --> FOAM FATAL ERROR:
[6] Negative initial temperature T0: -24430.2
[6]
[6] From function Foam::scalar Foam::species::thermo<Thermo, Type>::T(Foam::scalar, Foam::scalar, Foam::scalar, Foam::scalar (Foam::species::thermo<Thermo, Type>::*)(Foam::scalar, Foam::scalar) const, Foam::scalar (Foam::species::thermo<Thermo, Type>::*)(Foam::scalar, Foam::scalar) const, Foam::scalar (Foam::species::thermo<Thermo, Type>::*)(Foam::scalar) const, bool) const [with Thermo = Foam::hPolynomialThermo<Foam::icoPolynomial<Foam:: specie> >; Type = Foam::sensibleEnthalpy; Foam::scalar = double; Foam::species::thermo<Thermo, Type> = Foam::species::thermo<Foam::hPolynomialThermo<Foam ::icoPolynomial<Foam::specie> >, Foam::sensibleEnthalpy>]
[6] in file /opt/OpenFOAM/OpenFOAM-7/src/thermophysicalModels/specie/lnInclude/thermoI.H at line 56.
[6]


The mesh is only a block mesh, so nothing fancy and can (probably) be excluded as a cause of error.


Maybe I am totally mistaken with my BCs or thermophysical properties and just don't notice. I've tried tons of things and would appreciate any help.


The case is attached here: https://www.dropbox.com/sh/28zkn83n2...y3DEyriFa?dl=0




Thanks in advance,


Maike

[dlahaye]

1/ on short term: negative temperatures might arise due to (too) high temperatures in others parts of the domain. In that case it might be valuable to use a limiter on the temperature as suggested in post #7 by HPE in this thread: Why is it difficult to get rhoSimpleFoam solver running correctly ? ;

2/ on long term: in might be valuable to list all the options you have tried along with the outcome. This overview might be valuable to you. It might be valuable to post this overview, allowing others to think along. I will be interested in giving this overview a look.

Good luck.

[maike151]

Hello dlahaye,


thank you for your answer!


I tried limitTemperature (as in constant/fvOptions.backup) but the error still occurs. It's like this is just ignored even though it is recognized in initialization:


"Selecting finite volume options model type limitTemperature
Source: limitT
- selecting all cells
- selected 13248 cell(s) with volume 6.4e-05".


What else did I try..


1. Simplifying the mesh/make it coarser/finer


2. In fvSolution SIMPLE options put
-rhoMin/rhoMax or rhoMin [ 1 -3 0 0 0 ]/rhoMax [ 1 -3 0 0 0 ] gives this:


--> FOAM Warning :
From function Foam:ressureControl:ressureControl(const volScalarField&, const volScalarField&, const Foam::dictionary&, bool)
in file cfdTools/general/pressureControl/pressureControl.C at line 125
Reading "/home/schick/OpenFOAM/schick-7/run/2020-10-05_Schlebusch_SolarInstitutJuelich_Mischelement/TestCase4/skaliert_vollaufgeloest/schneller/ohneGewebe/system/fvSolution.SIMPLE" from line 83 to line 98
'rhoMax' specified rather than 'pMax' or 'pMaxFactor'
This is supported for backward-compatibility but 'pMax' or 'pMaxFactor' are more reliable.
--> FOAM Warning :
From function Foam:ressureControl:ressureControl(const volScalarField&, const volScalarField&, const Foam::dictionary&, bool)
in file cfdTools/general/pressureControl/pressureControl.C at line 182
Reading "/home/schick/OpenFOAM/schick-7/run/2020-10-05_Schlebusch_SolarInstitutJuelich_Mischelement/TestCase4/skaliert_vollaufgeloest/schneller/ohneGewebe/system/fvSolution.SIMPLE" from line 83 to line 98
'rhoMin' specified rather than 'pMin' or 'pMinFactor'
This is supported for backward-compatibility but'pMin' or 'pMinFactor' are more reliable.
pressureControl
pMax 319989
pMin 24614.5



-pMin/pMax or pMinFactor/pMaxFactor: still same error


3. Try thousands of combinations of relaxationFactors, even about 0.0001 for every field and equation: still same error


4. Try all upwind schemes: still same error


5. Try different BCs:
-For U: Try smaller mass flow rates
-For p: Try fixedValue and totalPressure at the outlet
-For T: Try zeroGradient and inletOutlet at the outlet
-For k and epsilon: Try different possible values


That's all I can think of right now.

[dlahaye]

It might be the case that OpenFoam prints T values before limiting them. Does limiting T have no influence whatsoever on the convergence?

Does it make sense to try an isothermal case using simpleFoam first?

[maike151]

No, there is no difference at all.


Sure, I can try that.

[dlahaye]

edit:

The fact that your T explodes after one iteration is likely to explain why limitT does not work. What is your domain like? What are your boundary conditions? What is your initial condition for temperature do you impose?

Let us try simpleFoam first.

[maike151]

So I did some changes in order to use simpleFoam


-delete thermopyhsical properties
-delete everything related to T, rho, h,...

-changed U inlet BC back to fixedValue
-changed p outlet BC back to fixedValue 0 and dimension of p back to [0 2 -2 0 0 0 0];

-changed kinematic velocity in transport properties to the average of the kinematic velocities the air at the two inlet temperatures had (not sure about this)


Then, in order to achieve my convergence criteria, refined the mesh by factor 1.5. Started running with

relaxationFactors
{

    fields
    {
    p        0.1;
    }    

    equations
    {
        U               0.7;
        k               0.1;
        omega           0.1;
        epsilon        0.1;
    }
}

Then, after about 1000 timesteps, changed to

relaxationFactors
{

    fields
    {
    p        0.5;
    }    

    equations
    {
        U               0.8;
        k               0.3;
        omega           0.3;
        epsilon        0.3;
    }
}

and the solution converged in about 1900 iterations. You can find the case here: https://www.dropbox.com/sh/vfj9fm0g4...1Hc-cSg9a?dl=0

[maike151]

Originally Posted by dlahaye

edit:

The fact that your T explodes after one iteration is likely to explain why limitT does not work. What is your domain like? What are your boundary conditions? What is your initial condition for temperature do you impose?

Let us try simpleFoam first.

You may find all the files in the first post

[dlahaye]

Looks wonderful.

Do I understand it correctly that your mesh does not include boundary layer refinements?

What happens when increasing the inlet velocities? Would you need boundary layer refinement to capture the gradient in the velocity profile on the boundary?

[maike151]

No, it does not include boundary layer refinements.


How could I tell that I needed it? (I'm sorry, I'm obviously not a CFD expert)

[shock77]

Hi,


if the simulation crashes after the first iteration, in most cases the BCs are wrong. I would try running a laminar simulation with slip walls to reduce the potential error sources. After that try using a fixedValue for U that matches your flowrate and see whether it works with that.



What are your min and max values for T?


I am not familiar with rhoSimpleFoam, but usually the temperature is basically calculated by subtracting the kinitec energie from the total internal energy. That means, that you are getting really high velocites in some locations.

[dlahaye]

No worries.

You can make line plots of the main velocity component on line perpendicular to the flow.

You can compute yplus values in postprocessing stage using -postprocess yplus

You can include insert boundary layers in the mesh and see what the effect is on the convergence.

Your choice ;-)

[maike151]

Hi shock77, thanks for your input!


I will try out your suggestions.



Tmin is 274 K and Tmax 1400 K, because T at one inlet is 1223.15 K.

I am not familiar with rhoSimpleFoam, but usually the temperature is basically calculated by subtracting the kinitec energie from the total internal energy. That means, that you are getting really high velocites in some locations.

That makes sense. Maybe I should try limitVelocity instead of limitTemperature then?

[maike151]

Thank you dlahaye! I will try all of that!

[shock77]

When I have some troubles with negative Temperatures, I usually set Tmax to like 10000 and Tmin to 10 to ensure the correctness of the result. I would not try to limit the velocity, because I think your BCs are the problem. Ofc you can limit everything do that point, where you make it work, but what benefit does that give you? You cant say wheter your results are trustworthy or not.

[maike151]

Sorry that it took me so long to answer.


I first tried out dlahayes suggestions.



The yPlus value for my mesh without boundary layer refinement was:
min 27.57 max 173.92 average 89.6


Then I used the refineWallLayer function first with factor 0.5 and I got
min 16.19 max 157.08 average 46.13


With factor 0.4 I got
min 11.58 max 150.69 average 37.61


When I got any smaller I got convergence problems. I also attached the pictures of the velocity profile line plot. It does not really change from 0.5 to 0.4.


So then I tried my rhoSimpleFoam case with refinement 0.5 but still got the error. What yPlus value is desirable?

[dlahaye]

Not sure what good values of yplus are. Possibly you want to reduce the max value of yplus. Possibly it helps to visualize yplus in paraview.

In case that non-convergence persists, it might be valuable to look into the residual field and refine the mesh where the residual remains large.

[maike151]

I visualized yPlus in Paraview and yPlus was the largest along the edges parallel to flow direction. So in addition to refineWallLayer I refined the Mesh at those edges by defining those cells as a cellSet with topoSet and use refineMesh there. It reduced yPlus: min = 8.84904, max = 85.9036, average = 35.7021

A picture of the resulting mesh is attached.



For the simpleFoam case the residuals got down to 1e-4, but then it started diverging. A picture is also attached. I don't know if I am happy with that because something must still be wrong, right?


How can I visualize the residual field in paraView? I tried using postProcess -func residuals, but this only generates the dat-files in the postPorcessing folder.



Tried the rhoSimpleFoam case with this mesh, still same error :-(

[maike151]

Originally Posted by shock77

Hi,


if the simulation crashes after the first iteration, in most cases the BCs are wrong. I would try running a laminar simulation with slip walls to reduce the potential error sources. After that try using a fixedValue for U that matches your flowrate and see whether it works with that.



What are your min and max values for T?


I am not familiar with rhoSimpleFoam, but usually the temperature is basically calculated by subtracting the kinitec energie from the total internal energy. That means, that you are getting really high velocites in some locations.

Hello shock77,


I finally tried out your suggestions. What I did was:
-use the mesh I did with dlahayes help

-change turbulenceProperties back to laminar
-delete everything related to k, epsilon,...
-change wall BC for U to slip


But the error still occured. The negative initial temperature didn't even change a bit.


Then I put fixedValue for U inlet BC with the velocity that matches my flow rate - still same error. I also tried really small mass flow rates and velocities - no change at all.

[NaiveYoung]

hi,maike
I 've just met the same problem as you do.
I drew my mesh by NUMECA, which uses millimetre as scale, but OpenFoam uses meter as scale.

And I changed my scale to meter, and it works well.