[peob]

I've already submitted this to the OpenCFD Mantis bug report, but I thought I would post this here also, to see if anyone else has experienced this issue...

I'm running OpenFOAM-2.1.1 (but this also happens with version 2.1.0).
I can run the sprayFoam/aachenBomb tutorial just fine.

Then I turn on chemistry (in constant/chemistryProperities I change "chemistry off;" to "chemistry on;"), and run the case from scratch, and it also runs just fine.

Then I turn on chemistry (as above) and combustion ( in constant/combustionProperities I change "active false;" to "active true;"), and the case experiences a "FOAM FATAL ERROR" in speciesThermoI.H after about 0.000769 seconds of simulated time.

I've attached a copy of the log file (log.txt) with most of the intermediate steps deleted so it isn't so large.

If anyone else has experienced this problem, please comment.
If anyone knows what the problem is and how to fix it, please post a reply.

Thanks.
Phil

[kalle]

I have no resolution, but I can just say that I see exactly the same in several setups. I could not pinpoint the issue yet. Only clue so far is that is happens as soon as T gets in range 1350-1400 Kin all cases. therm.dat differs only from other version (1.6-ext) on species not in use for this tutorial.

I'll get back if I find something.

K

[kalle]

Ok... did some digging. Attached are plots generated with this code:

            for(label T=300; T<2400; T++)
            {
                Info << T << "   " << (this->*dFdT)(T) << "   " << (this->*F)(T) << endl;
            }

Which I put into the Foam::specieThermo<Thermo>::T
It seems there is a discontinuity in the enthalpy at ~1390 K, which likely cause the non-convergence. The fault is likely in therm.dat which have diffenrent values for the thermo compared to for instance

ftp://ftp.technion.ac.il/pub/support...mics/THERM.DAT

Changing therm.dat might help here

K

[ed_teller]

Hi

I placed the therm.dat file in the etc/thermoData folder, after renaming the old therm.dat, and it did not solve the problem... in which folder should i place it... ? in aachenBomb/chemkin...?


BR

Alex

[wyldckat]

Hi Alex,

Originally Posted by ed_teller

I placed the therm.dat file in the etc/thermoData folder, after renaming the old therm.dat, and it did not solve the problem... in which folder should i place it... ? in aachenBomb/chemkin...?

Did you remove the header text that is in the new file?

As for the file that should be replaced, you need to check the file "constant/thermophysicalProperties" that is in your case folder!

Best regards,
Bruno

[ed_teller]

Hi wyldckat

I changed the reference in thermophysicalProperties, pasted the file in chemkin folder, after removing the header... and launched the app... but it seesm that the therm.dat file is ill formed anyway... even after removing the header... :

Reading g
Creating combustion model

Selecting combustion model PaSR<psiChemistryCombustionModel>
Selecting psiChemistryModel ODEChemistryModel<gasThermoPhysics>
Selecting thermodynamics package hsPsiMixtureThermo<reactingMixture<gasThermoPhysics>>
Selecting chemistryReader chemkinReader
[0] 
[1] 
[1] 
[1] --> FOAM FATAL ERROR: 
[1] while reading thermodynamic specie name on line 3
'[2] r REF ELEMENT T 6/12AG 1.   0.   0.   0.S   200.000  1235.080 1000.        1
[2] 
[2] --> FOAM FATAL ERROR: 
[2] while reading thermodynamic specie name on line 3
'Ag cr REF ELEMENT T 6/12AG 1.   0.   0.   0.S   200.000  1235.080 1000.        1
[2] 
[3] 
[3] 
[3] --> FOAM FATAL ERROR: 
[3] while reading thermodynamic specie name on line 3
'Ag cr REF ELEMENT T 6/12AG 1.   0.   0.   0.S   200.000  1235.080 1000.        1
[3] 
[3]     From function chemkinReader::lex()
[3]     in file chemistryReaders/chemkinReader/chemkinLexer.L at line [2]     From function chemkinReader::lex()
[2]     in file chemistryReaders/chemkinReader/chemkinLexer.L at line 1497.
[3] 
FOAM parallel run exiting
[3] 
[0] 
[0] --> FOAM FATAL ERROR: 
[0] while reading thermodynamic specie name on line 3
'Ag cr REF ELEMENT T 6/12AG 1.   0.   0.   0.S   200.000  1235.080 1000.        1
[0] 
[0]     From function chemkinReader::lex()
[0]     in file chemistryReaders/chemkinReader/chemkinLexer.L at line 1497.
[0] 
FOAM parallel run exiting
[0] 

[1] 
[1]     From function chemkinReader::lex()
[1]     in file chemistryReaders/chemkinReader/chemkinLexer.L at line 1497.
[1] 
FOAM parallel run exiting
[1] 
1497.
[2] 
FOAM parallel run exiting
[2] 
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD 
with errorcode 1.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.
--------------------------------------------------------------------------
--------------------------------------------------------------------------
mpirun has exited due to process rank 0 with PID 4846 on
node iskandhar exiting improperly. There are two reasons this could occur:

1. this process did not call "init" before exiting, but others in
the job did. This can cause a job to hang indefinitely while it waits
for all processes to call "init". By rule, if one process calls "init",
then ALL processes must call "init" prior to termination.

2. this process called "init", but exited without calling "finalize".
By rule, all processes that call "init" MUST call "finalize" prior to
exiting or it will be considered an "abnormal termination"

This may have caused other processes in the application to be
terminated by signals sent by mpirun (as reported here).
--------------------------------------------------------------------------
[iskandhar:04845] 3 more processes have sent help message help-mpi-api.txt / mpi-abort
[iskandhar:04845] Set MCA parameter "orte_base_help_aggregate" to 0 to see all help / error messages

the refered line in therm.dat file is as follows... :

Ag cr REF ELEMENT T 6/12AG 1. 0. 0. 0.S 200.000 1235.080 1000. 1


the corresponding line in etc/thermoData/them.dat is as follows :

(CH2O)3 70590C 3H 6O 3 G 0200.00 4000.00 1500.00 1

A slightly different format....

BR

Alex

[wyldckat]

Mmm... I guess you're meant only to replace the lines that you actual need...

[ris]

Hello,

I have also experienced the problem mentioned above when turning on combustion and chemistry.
I modified the original C7H16 in therm.dat with the one linked here: ftp://ftp.technion.ac.il/pub/support...mics/THERM.DAT

and rerun the simulation only to receive the same error as before at around 0.00077 seconds. I am not sure whether I should change the coefficients for all the species within the default reaction scheme.
I am using OpenFOAM-2.1.1, please could someone confirm whether this problem is present in the GIT version of OpenFOAM.

Many Thanks,

ris

[wyldckat]

Greetings Ris,

Originally Posted by ris

I am using OpenFOAM-2.1.1, please could someone confirm whether this problem is present in the GIT version of OpenFOAM.

If you can provide an example case, or give instructions based on one of the existing tutorials, I can test it in OpenFOAM 2.1.x for you.

Best regards,
Bruno

[ris]

Hello Bruno,

The case in question is the tutorial aachenbomb located in /tutorials/lagrangian/sprayFoam/.
The only changes made to the case were:

Turning chemistry on in constant/chemistryProperties

and

Turning combustion to true in contant/combustionProperties.

I have attached an example case with the changes above.

Thank you for your assistance,

Regards,

ris

[wyldckat]

Hi Ris,

OK, so after dwelling a bit on this I finally figured out what Kalle was talking about.

The following changes seem to work:

1. Copy the "therm.dat" file from OpenFOAM's global location:
   

   cp $WM_PROJECT_DIR/etc/thermoData/therm.dat chemkin/

2. Edit the file "constant/thermophysicalProperties" and change the location of the "therm.dat" file:
   

   //CHEMKINThermoFile "~OpenFOAM/thermoData/therm.dat";
   CHEMKINThermoFile "$FOAM_CASE/chemkin/therm.dat";

3. Edit this "chemkin/therm.dat" file and locate "C7H16". It should have this:
   

   C7H16             P10/95C   7H  16    0    0G   200.000  5000.000  1391.000    1
    2.22148969e+01 3.47675750e-02-1.18407129e-05 1.83298478e-09-1.06130266e-13    2
   -3.42760081e+04-9.23040196e+01-1.26836187e+00 8.54355820e-02-5.25346786e-05    3
    1.62945721e-08-2.02394925e-12-2.56586565e+04 3.53732912e+01                   4

4. Locate the respective entry in the aforementioned "THERM.DAT" from that ftp location:
   

   C7H16,n-heptane   P10/85C  7.H 16.   0.   0.G   200.000  6000.000 1000.        1 
    2.04565203E+01 3.48575357E-02-1.09226846E-05 1.67201776E-09-9.81024850E-14    2 
   -3.25556365E+04-8.04405017E+01 1.11532994E+01-9.49419773E-03 1.95572075E-04    3 
   -2.49753662E-07 9.84877715E-11-2.67688904E+04-1.59096837E+01-2.25846141E+04    4

   I'm assuming this is the relevant entry, based on the searches I've made, where the "C7H16" that is defined in OpenFOAM is "NC7H16" in other databases.
5. On our "chemkin/therm.dat" file, we should merge the two and get this:
   

   C7H16             P10/85C  7.H 16.   0.   0.G   200.000  6000.000 1000.        1
    2.04565203E+01 3.48575357E-02-1.09226846E-05 1.67201776E-09-9.81024850E-14    2
   -3.25556365E+04-8.04405017E+01 1.11532994E+01-9.49419773E-03 1.95572075E-04    3
   -2.49753662E-07 9.84877715E-11-2.67688904E+04-1.59096837E+01-2.25846141E+04    4

With these changes, I was able to go beyond 0.00088 seconds and reached a max temperature of 2448.29 K. I didn't bother continuing the simulation, since it takes sooooo long to run


Now we come to the big problem: how on Earth is this going to be validated? The values present in OpenFOAM were possibly validated with some project, problem is with which one?

Best regards,
Bruno

PS: the mantis bug report mentioned in the first post is this one: http://www.openfoam.org/mantisbt/view.php?id=574

[ris]

Hi Bruno,

Thank you very much for this workaround, I will make the changes you have specified and update this thread with the outcome.
With respect to validating this workaround, I should be able to get access to experimental data regarding diesel injection and will try to see how the simulation will compare. Its not ideal but it will give us an idea on what is happening.

Many Thanks,

ris

[kalle]

Sorry Bruno for not clarifying enough what I meant with my post. You've got my idea anyhow.

I just got inspired so I decided to calculate Cp and H @Tcommon manually using the data in therm.dat retrieved July 14 from OpenCFD. You can paste the following into python to see for yourself:

T=1391
#Cp/R
2.22148969e+01 +  3.47675750e-02*T  -1.18407129e-05*pow(T,2) +  1.83298478e-09*pow(T,3) -1.06130266e-13*pow(T,4)        # -3.42760081e+04 -9.23040196e+01

-1.26836187e+00 + 8.54355820e-02*T -5.25346786e-05*pow(T,2) + 1.62945721e-08*pow(T,3) -2.02394925e-12*pow(T,4)            #-2.56586565e+04 3.53732912e+01        

#H/RT
2.22148969e+01 +  3.47675750e-02*T/2  -1.18407129e-05*pow(T,2)/3 +  1.83298478e-09*pow(T,3)/4 -1.06130266e-13*pow(T,4)/5 -3.42760081e+04/T #-9.23040196e+01

-1.26836187e+00 + 8.54355820e-02*T/2 -5.25346786e-05*pow(T,2)/3 + 1.62945721e-08*pow(T,3)/4 -2.02394925e-12*pow(T,4)/5 -2.56586565e+04/T # 3.53732912e+01

... and surprise is: The two polynomials give very close values at Tcommon = 1391! Now, I am confused... The figures I generated in the earlier post clearly shows that the crash-causing cell has corrupt thermodata with a discontinuity at 1391 K. Since C7H16 had it's Tcommon at that temperature I assumed the therm.dat was the problem. Likely I sign error on a coefficient or something like that.

How to use thermodata: http://www.me.berkeley.edu/gri_mech/data/nasa_plnm.html

I have no clue right now what is wrong here... fishy.

K

[kalle]

The way janafThermo defines enthalpy for mixtures (weighting together their coefficients, see operators+= and operator+ in janafThermoI.H) does not allow for using different Tcommon for different species in the same mixture. Janaf should crash upon this in the case debug switches are activated. I did not manage to trigger that however. In this case; maybe C7H16 data is OK, but Tcommon for that specie (1391) is used also for mixtures containing also other species (in this case all the other species have Tcommon = 1000K). As those species have inconsistent Cp and H at 1391K for their two sets of coefficient, which causes non-convergence there. Indeed, the Cp-figure in my previous post shows reasonable Cp up to 1000K.

Conclusion: be careful not to mix thermodata's with different Tcommons. I know Cantera for instance, will complain when you use Chemkin thermo data with several Tcommons. OpenFOAM does not check this by default, but it causes obviously confusing problems. The other thermo data for C7H16 in this thread have Tcommon=1000K and the issue goes away.

K

[wyldckat]

Very nice explanation Kalle!

But now that I think of it, this looks familiar... I think this is still an incomplete feature in OpenFOAM which was discussed sometime ago here: http://www.openfoam.org/mantisbt/view.php?id=327 - with any hope, it'll be properly fixed in the next version of OpenFOAM...

[kalle]

I did not know about that bug report. Thanks for point that out, Bruno! I have updated the other bug report about our findings here.

For now on, I'd guess we should stick to thermo data with only one Tcommon.

K

[immortality]

hello
I used janaf coefficients in openfoam but an error in thermo model displayed.i draw that in matlab and it showed a large discontinuity in 1000K.i used coefficients that come in openfoamwiki and don't know what is wrong.how you used the coeffs?i need air or nitropen and oxygen as a mixture.
Thanks for helping.

[wyldckat]

For more information about Ehsan's problem, check his more complete thread here: http://www.cfd-online.com/Forums/ope...incorrect.html

Best regards,
Bruno

[CFDUser_]

thankyou wyldckat

[u22]

Dear all,

Sorry, if a little bit off topic.
Why etc/thermoData/therm.dat includes 15 coeffs for each specie? Why not 14 (7 for low + 7for high temperature)?

And why in etc/thermoData/thermData imported only 14 from 15? Is it correct to use 15 coeffs if the OF accepts 14 only?