CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Running, Solving & CFD

Negative temperature and pressure in sonicFoam :O

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree7Likes
  • 1 Post By KarenRei
  • 1 Post By KarenRei
  • 1 Post By gallon
  • 3 Post By Vyssion
  • 1 Post By Icemaan

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   August 22, 2016, 06:22
Red face Negative temperature and pressure in sonicFoam :O
  #1
New Member
 
Vladislav Galas
Join Date: Jun 2016
Posts: 6
Rep Power: 8
gallon is on a distinguished road
Hello everyone, here's a problem I've been trying to solve for a long time.
I have an stl model of a bullet and I want to calculate its supersonic movement under normal conditions.
However at some moment the solver yields negative pressure in the trail => negative temperature => fails.


Now for the details.
1. Case setup.
compressible/sonicFoam/ras/nacaAirfoil was used as an example/template for BC.
For time step I define ~1e-9 by default, but I tried steps down to 1e-12, the result is the same. With these steps, inlet velocity 850 m/s and cell size down to 1e-4 m the Courant number is about ~0.001 which should be more than sufficient.

2. Mesh
I tried sHM, cfMesh and ICEM CFD. Layers or no layers, the general behaviour is the same: the larger the cell size and the smaller the velocity, the farther the computation goes before failing. At velocity ~400 m/s and cell size ~1-2 mm the computations manages to come to a steady solution but it is wrong because the shock wave is not resolved.

3. The error
Here's the output:
Code:
[185] --> FOAM FATAL ERROR:
[185] Maximum number of iterations exceeded
[185]
[185]     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) const [with Thermo = Foam::hConstThermo<Foam::perfectGas<Foam::specie> >; Type = Foam::sensibleInternalEnergy; Foam::scalar = double; Foam::species::thermo<Thermo, Type> = Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleInternalEnergy>]
[185]     in file /common/software/OpenFOAM/OpenFOAM-4.0/src/thermophysicalModels/specie/lnInclude/thermoI.H at line 69.
I opened thermoI.H, added some printfs, recompiled and rerun. The T is negative at that point, so it is not the root of the problem.

4. The search for success
What else I tried?.. Switching schemes in fvSchemes, a hundred different mesh setups, adjusting timestep

Yet I took `nacaairfoil` example, made it 3D by changing empty->zeroGradient and refineMesh command and directed the inlet normally to the wing plane to test. Well, it worked fine! Even at 7 Machs...

Here's the full case setup attached. The bullet.stl is changed to a simple handmade model for understandable reasons. I checked, the result is the same. PS controlDict is for OF 4. You can disable postprocessing for compatibility.
To run:
Code:
blockMesh
surfaceFeatureExtract
mesh
sonicFoam
Decomposition is strongly advised for perfomance.
bullet.zip

PPS I tried this on another model but there's a good chance that if you flip the wind (exchange inlet with outlet and swap sign in 0/defines/U), everything will go well.
Attached Images
File Type: png T.PNG (15.4 KB, 1472 views)
gallon is offline   Reply With Quote

Old   September 13, 2016, 06:08
Default
  #2
Member
 
Anonymouse
Join Date: Dec 2015
Posts: 98
Rep Power: 9
KarenRei is on a distinguished road
I'm having the same problem with sonicFoam I'm not as far along as you at finding solutions as it's taken me days to compile the latest version of OpenFOAM (I'm on Fedora, so building it is always a pain, the makefiles are missing libs, some libs get put in the wrong place, etc)

Temperature = -0.22K on my last run for some reason... low courant number... no luck.
KarenRei is offline   Reply With Quote

Old   September 15, 2016, 07:36
Default
  #3
Member
 
Anonymouse
Join Date: Dec 2015
Posts: 98
Rep Power: 9
KarenRei is on a distinguished road
Found the solution in my case I dropped the p tolerance in fvsolutions by two orders of magnitude, that totally stabilized the run (and let me dramatically increase my timestep!)
student666 likes this.
KarenRei is offline   Reply With Quote

Old   September 25, 2016, 16:48
Default
  #4
New Member
 
Vladislav Galas
Join Date: Jun 2016
Posts: 6
Rep Power: 8
gallon is on a distinguished road
Quote:
Originally Posted by KarenRei View Post
Found the solution in my case I dropped the p tolerance in fvsolutions by two orders of magnitude, that totally stabilized the run (and let me dramatically increase my timestep!)
Thanks for the info. Didn't work for me, unfortunately.
gallon is offline   Reply With Quote

Old   September 26, 2016, 00:40
Default
  #5
Member
 
Anonymouse
Join Date: Dec 2015
Posts: 98
Rep Power: 9
KarenRei is on a distinguished road
Quote:
Originally Posted by gallon View Post
Thanks for the info. Didn't work for me, unfortunately.
Unfortunate. I've been doing a lot of different runs with different solvers, and randomly running into such negative temperature excursions... and each time it really seems like it's just guesswork to find things that work. It's almost always by some combination of tweaking fvSchemes, fvSolutions, and increasing the mesh resolution in sensitive areas where oscillatory or turbulent behavior spirals out of control. But I've not yet found some sort of easy, always-works solution. As you noted, shrinking the timestep rarely seem to be the solution if the courant number is already sufficiently low.

What's your fvSchemes that you're using right now? And what does "checkMesh -constant -allGeometry -allTopology" say? No problems with the mesh? Orthogonality? Skewness?

Oh, one more thing: one "excursion" case I had was due to me setting up the geometry wrong; I thought I had linked some faces up proper boundary conditions but it had been the wrong faces, leaving some faces with boundaries that they shouldn't have had and some faces with no boundary conditions set whatsoever. And it led to some really weird behavior. So probably good to double check your geometry - sometimes geometry can "pass" all checks but still be wrong.
Paul Caicedo likes this.
KarenRei is offline   Reply With Quote

Old   September 26, 2016, 17:36
Default
  #6
New Member
 
Vladislav Galas
Join Date: Jun 2016
Posts: 6
Rep Power: 8
gallon is on a distinguished road
Quote:
Originally Posted by KarenRei View Post
It's almost always by some combination of tweaking fvSchemes, fvSolutions, and increasing the mesh resolution in sensitive areas where oscillatory or turbulent behavior spirals out of control. But I've not yet found some sort of easy, always-works solution. As you noted, shrinking the timestep rarely seem to be the solution if the courant number is already sufficiently low.
The weird thing is that the finer I make the mesh the quicker negative p/T is achieved, even if timestep is reduced accordingly... while the solution is generally expected to converge better on a finer grid.

Quote:
Originally Posted by KarenRei View Post
What's your fvSchemes that you're using right now? And what does "checkMesh -constant -allGeometry -allTopology" say? No problems with the mesh? Orthogonality? Skewness?
Oh, one more thing: one "excursion" case I had was due to me setting up the geometry wrong; I thought I had linked some faces up proper boundary conditions but it had been the wrong faces, leaving some faces with boundaries that they shouldn't have had and some faces with no boundary conditions set whatsoever. And it led to some really weird behavior. So probably good to double check your geometry - sometimes geometry can "pass" all checks but still be wrong.
I guess I tried almost every combination of numeric schemes
checkMesh usually tells everything is ok. I also tried both real STL from CAD and a simple self-made "projectile" formed as a half of an ellipsoid.
The geometry looks rather simple: a cylinder (one circle for inlet, other parts for outlet) with a wall inside...

If only someone reported the reason of this weird behavior!
Paul Caicedo likes this.
gallon is offline   Reply With Quote

Old   September 28, 2016, 12:53
Default
  #7
New Member
 
Vyssion's Avatar
 
Vyssion
Join Date: Mar 2016
Posts: 12
Rep Power: 8
Vyssion is on a distinguished road
One thing that I have used successfully in the past is to set up a temperature limit within an fvOptions file:

Code:
TemperatureLimit1
{
	type		temperatureLimitsContraint;
	active		true;
	selectionMode	all;

	temperatureLimitsConstraintCoeffs
	{
		Tmin	 200;
		Tmax	 2000;
	}
}
student666, KarenRei and Ramzy1990 like this.
__________________
Optimism, pessimism, f*ck that: We're going make it happen. As God is my bloody witness, I'm hell-bent on making it work.
- Elon Musk in 2008, after three unsuccessful Falcon launches.
Vyssion is offline   Reply With Quote

Old   September 28, 2016, 20:39
Default
  #8
Member
 
Anonymouse
Join Date: Dec 2015
Posts: 98
Rep Power: 9
KarenRei is on a distinguished road
Ooh neat - didn't know that was possible!
KarenRei is offline   Reply With Quote

Old   October 6, 2016, 18:54
Default
  #9
Member
 
Anonymouse
Join Date: Dec 2015
Posts: 98
Rep Power: 9
KarenRei is on a distinguished road
Quote:
Originally Posted by Vyssion View Post
One thing that I have used successfully in the past is to set up a temperature limit within an fvOptions file:

Code:
TemperatureLimit1
{
	type		temperatureLimitsContraint;
	active		true;
	selectionMode	all;

	temperatureLimitsConstraintCoeffs
	{
		Tmin	 200;
		Tmax	 2000;
	}
}
Hmm. This isn't working for me.

First off, apparently in newer versions of OpenFoam it's now "limitTemperature", not "temperatureLimitsConstraint". To get it to start running, I ended up having to go with (in system/fvOptions):

temperatureLimit
{
type limitTemperature;
active true;

limitTemperatureCoeffs
{
selectionMode all;
min 100;
max 5000;
}
}

(I don't know what's supposed to go on that "temperatureLimit line, I think it's just a name?)

It's clearly parsing this file because if I have any sort of typo it doesn't start. But the limits don't actually seem to limit anything:

--> FOAM Warning :
From function Foam::scalar Foam::janafThermo<EquationOfState>::limit(Foam::sc alar) const [with EquationOfState = Foam:erfectGas<Foam::specie>; Foam::scalar = double]
in file /opt/OpenFOAM/OpenFOAM-dev/src/thermophysicalModels/specie/lnInclude/janafThermoI.H at line 105
attempt to use janafThermo<EquationOfState> out of temperature range 1e-06 -> 6000; T = -1.3418404
--> FOAM Warning :
From function Foam::scalar Foam::janafThermo<EquationOfState>::limit(Foam::sc alar) const [with EquationOfState = Foam:erfectGas<Foam::specie>; Foam::scalar = double]
in file /opt/OpenFOAM/OpenFOAM-dev/src/thermophysicalModels/specie/lnInclude/janafThermoI.H at line 105
attempt to use janafThermo<EquationOfState> out of temperature range 1e-06 -> 6000; T = -1.3418404
--> FOAM Warning :
From function Foam::scalar Foam::janafThermo<EquationOfState>::limit(Foam::sc alar) const [with EquationOfState = Foam:erfectGas<Foam::specie>; Foam::scalar = double]
in file /opt/OpenFOAM/OpenFOAM-dev/src/thermophysicalModels/specie/lnInclude/janafThermoI.H at line 105
attempt to use janafThermo<EquationOfState> out of temperature range 1e-06 -> 6000; T = -1.3418404


--> FOAM FATAL ERROR:
Maximum number of iterations exceeded

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) const [with Thermo = Foam::janafThermo<Foam:erfectGas<Foam::specie> >; Type = Foam::sensibleInternalEnergy; Foam::scalar = double; Foam::species::thermo<Thermo, Type> = Foam::species::thermo<Foam::janafThermo<Foam:erf ectGas<Foam::specie> >, Foam::sensibleInternalEnergy>]
in file /opt/OpenFOAM/OpenFOAM-dev/src/thermophysicalModels/specie/lnInclude/thermoI.H at line 66.

FOAM aborting

Am I doing something wrong?
KarenRei is offline   Reply With Quote

Old   October 7, 2016, 07:38
Default
  #10
New Member
 
Vyssion's Avatar
 
Vyssion
Join Date: Mar 2016
Posts: 12
Rep Power: 8
Vyssion is on a distinguished road
Quote:
Originally Posted by KarenRei View Post
First off, apparently in newer versions of OpenFoam it's now "limitTemperature", not "temperatureLimitsConstraint"
I don't use the latest version of OpenFOAM so it may be different yes.

Quote:
Originally Posted by KarenRei View Post
--> FOAM Warning :
From function Foam::scalar Foam::janafThermo<EquationOfState>
This janafThermo error is a nightmare to deal with sometimes - your boundary conditions and/or thermalProperties script may not have the correct thermal properties for your setup.

There isn't too much help that I can offer with this issue tbh... I did come across this thread which someone was whinging about it aimlessly - but it did have some links to other threads which apparently helped other people, so I wish you luck... sorry
__________________
Optimism, pessimism, f*ck that: We're going make it happen. As God is my bloody witness, I'm hell-bent on making it work.
- Elon Musk in 2008, after three unsuccessful Falcon launches.
Vyssion is offline   Reply With Quote

Old   October 7, 2016, 11:01
Default
  #11
Member
 
Join Date: Jul 2013
Posts: 39
Rep Power: 11
cfdsolver1 is on a distinguished road
I am having kind of same problem with chtMultiRegionSimpleFoam. Negative temperature is a problem in my case.
cfdsolver1 is offline   Reply With Quote

Old   October 7, 2016, 12:22
Default
  #12
New Member
 
Vladislav Galas
Join Date: Jun 2016
Posts: 6
Rep Power: 8
gallon is on a distinguished road
Quote:
Originally Posted by KarenRei View Post
Hmm. This isn't working for me.
True, for OF 4.0 syntax has to be adjusted, but it's done easily.

I track min/max p/U/T using postprocessing, so it worked for me like this (p/T/U respectively):


Looks like the constraint tries to work but is applied too seldom or something.

You can try to add this to controlDict to get minmax postprocessing:
Code:
functions
{
    minMax
    {
        type            fieldMinMax;
        libs ("libfieldFunctionObjects.so");
        enabled         true;
        log             true;
        write           true;
        // Fields to be monitored - runTime modifiable
        fields(T U p);
    }
}
This is the script I wrote to view data with gnuplot in realtime:
Code:
#!/usr/bin/bash
set -e

if [ $# != 1 ]; then
  echo "Usage: gp.sh <FILE>"
  exit -1
fi

glaunch () {
cat << EOF > gpscript."$2"
  plot "$1" every 3::$2 using 1:3 title 'x' with points, "$1" every 3::$2 using 1:8 title 'x' with points
  pause 1
  reread
EOF

gnuplot ./gpscript."$2"
}

glaunch $1 0 &
glaunch $1 1 &
glaunch $1 2 &
Simply issue something like:
Code:
$ ./this_script.sh ./postProcessing/minMax/0/fieldMinMax.dat
having launched the solver or after it finishes.
Attached Images
File Type: jpg tLimit.jpg (48.0 KB, 1425 views)
gallon is offline   Reply With Quote

Old   November 1, 2017, 05:25
Default Good news
  #13
New Member
 
Vladislav Galas
Join Date: Jun 2016
Posts: 6
Rep Power: 8
gallon is on a distinguished road
Good news everybody. I tweaked system/fvSchemes so that the solver doesn't fail.
Code:
gradSchemes
{
    default cellLimited Gauss linear 1; //default: Gauss linear;
}

divSchemes
{
    div(phi,U) Gauss linearUpwind cellLimited Gauss linear 1; // default: Gauss limitedLinearV 1;
}

laplacianSchemes
{
    default Gauss linear limited 0.333; // default: Gauss linear limited corrected 0.5;
}

snGradSchemes
{
    default limited 0.333; // default: corrected;
}
Honestly, I have no idea how and why this works but it works. Probably some of these changes are redundant but right now I am not ready to try it. Anyway, it works as is.

Perhaps fellow members could explain why.
gallon is offline   Reply With Quote

Old   January 19, 2018, 03:35
Post
  #14
New Member
 
Tarandeep
Join Date: Dec 2015
Posts: 3
Rep Power: 8
Icemaan is on a distinguished road
Quote:
Originally Posted by gallon View Post
Good news everybody. I tweaked system/fvSchemes so that the solver doesn't fail.
Code:
gradSchemes
{
    default cellLimited Gauss linear 1; //default: Gauss linear;
}

divSchemes
{
    div(phi,U) Gauss linearUpwind cellLimited Gauss linear 1; // default: Gauss limitedLinearV 1;
}

laplacianSchemes
{
    default Gauss linear limited 0.333; // default: Gauss linear limited corrected 0.5;
}

snGradSchemes
{
    default limited 0.333; // default: corrected;
}
Honestly, I have no idea how and why this works but it works. Probably some of these changes are redundant but right now I am not ready to try it. Anyway, it works as is.

Perhaps fellow members could explain why.
Hi,

I had the same negative temperature issue while solving for Mach 2 external flow using sonicFoam. I did a similar thing. I guess one of the reason this tweak works is because of the order of the divergence schemes. I only changed the divergence scheme from 'llimitedLinear V 1' (which is second order) to 'Gauss Upwind' (which is first order). This took care of this problem for my simulation. Also, after obtaining a good enough solution with first order scheme, I was able to change it to second order and continue, although at a smaller timestep. I did get the same negative temperature problem again with the second order simulation, but this time it was only due to trying to up the timestep too fast.
Ramzy1990 likes this.
Icemaan is offline   Reply With Quote

Old   April 12, 2018, 00:45
Default
  #15
New Member
 
Paul
Join Date: Nov 2017
Posts: 15
Rep Power: 7
Paul Caicedo is on a distinguished road
Quote:
Originally Posted by gallon View Post
True, for OF 4.0 syntax has to be adjusted, but it's done easily.


Looks like the constraint tries to work but is applied too seldom or something.

You can try to add this to controlDict to get minmax postprocessing:
Code:
functions
{
    minMax
    {
        type            fieldMinMax;
        libs ("libfieldFunctionObjects.so");
        enabled         true;
        log             true;
        write           true;
        // Fields to be monitored - runTime modifiable
        fields(T U p);
    }
}
Hello Vladislav,

Perhaps, do you know how to set this in version 5.x? As an example, in previous versions to use MRF you needed to include the data in the system/fvOptions file, but now it is located in constant/MRFProperties file.

Regards,
Paul
Paul Caicedo is offline   Reply With Quote

Old   August 27, 2018, 05:59
Default
  #16
Member
 
Foad
Join Date: Aug 2017
Posts: 55
Rep Power: 7
foadsf is on a distinguished road
Adding a fvOptions file in system folder with :

Code:
temperatureLimit
{
    type limitTemperature;
    active true;

    limitTemperatureCoeffs
    {
        selectionMode all;
        min <Tmin>;
        max <Tmax>;
    }
}
seems to resolve the issue.
foadsf is offline   Reply With Quote

Old   January 21, 2019, 06:56
Default Temperature problem
  #17
New Member
 
Athanasios Niotis
Join Date: Aug 2018
Posts: 12
Rep Power: 6
thanosniotis is on a distinguished road
Quote:
Originally Posted by foadsf View Post
Adding a fvOptions file in system folder with :

Code:
temperatureLimit
{
    type limitTemperature;
    active true;

    limitTemperatureCoeffs
    {
        selectionMode all;
        min <Tmin>;
        max <Tmax>;
    }
}
seems to resolve the issue.
That is good! I will implement it and we will see what happens.
But I suspect that the problem with the negative temperature exists when I use turbulence flow (Kepsilon), in case of laminar flow everything works fine. I dont know if any fellow foamer has noticed something relevant.
thanosniotis is offline   Reply With Quote

Reply

Tags
sonicfoam, sonicfoam convergence, supersonic

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Wind tunnel Boundary Conditions in Fluent metmet FLUENT 6 October 30, 2019 13:23
Ansys CFX problem: unexpected very high temperatures in premix laminar combustion faizan_habib7 CFX 4 February 1, 2016 18:00
Calculation of the Governing Equations Mihail CFX 7 September 7, 2014 07:27
Water subcooled boiling Attesz CFX 7 January 5, 2013 04:32
what the result is negatif pressure at inlet chong chee nan FLUENT 0 December 29, 2001 06:13


All times are GMT -4. The time now is 22:25.