# Temperature calculation from total enthalpy in OpenFOAM (XiFOAM)

 Register Blogs Members List Search Today's Posts Mark Forums Read August 25, 2017, 13:10 Temperature calculation from total enthalpy in OpenFOAM (XiFOAM) #1 New Member   Jim KIT Join Date: Aug 2012 Location: Germany Posts: 25 Rep Power: 10 Hi all, I'm currently working with XiFoam-Solver (OF 2.4) and I'm trying to understand how exactly the temperature is calculated from the total enthalpy. Actually I'm looking for a reference where it is explained clearly. It would be great If you could share such a reference with me if you are aware of any. As far as I understood, in XiFoam we are dealing with one step chemistry reactants --> products So it is possible to calculate the Nasa polynomial constants using the data for let's say methan-air combustion using equivalence ratio and "thermoData" of OpenFOAM. The temperature is calculated by calling thermo.correct(); from EaEqn.H Because the XiFoam uses: Code: ```thermoType { type heheuPsiThermo; mixture homogeneousMixture; transport sutherland; thermo janaf; equationOfState perfectGas; specie specie; energy absoluteEnthalpy; }``` The function correct(); can be found in "heheuPsiThermo.C" as: Code: ```template void Foam::heheuPsiThermo::correct() { if (debug) { Info<< "entering heheuPsiThermo::correct()" << endl; } // force the saving of the old-time values this->psi_.oldTime(); calculate(); if (debug) { Info<< "exiting heheuPsiThermo::correct()" << endl; } }``` The function calculate can be found in the same h-data: Code: ```template void Foam::heheuPsiThermo::calculate() { const scalarField& hCells = this->he_.internalField(); const scalarField& heuCells = this->heu_.internalField(); const scalarField& pCells = this->p_.internalField(); scalarField& TCells = this->T_.internalField(); scalarField& TuCells = this->Tu_.internalField(); scalarField& psiCells = this->psi_.internalField(); scalarField& muCells = this->mu_.internalField(); scalarField& alphaCells = this->alpha_.internalField(); forAll(TCells, celli) { const typename MixtureType::thermoType& mixture_ = this->cellMixture(celli); TCells[celli] = mixture_.THE ( hCells[celli], pCells[celli], TCells[celli] ); psiCells[celli] = mixture_.psi(pCells[celli], TCells[celli]); muCells[celli] = mixture_.mu(pCells[celli], TCells[celli]); alphaCells[celli] = mixture_.alphah(pCells[celli], TCells[celli]); TuCells[celli] = this->cellReactants(celli).THE ( heuCells[celli], pCells[celli], TuCells[celli] ); } forAll(this->T_.boundaryField(), patchi) { fvPatchScalarField& pp = this->p_.boundaryField()[patchi]; fvPatchScalarField& pT = this->T_.boundaryField()[patchi]; fvPatchScalarField& pTu = this->Tu_.boundaryField()[patchi]; fvPatchScalarField& ppsi = this->psi_.boundaryField()[patchi]; fvPatchScalarField& ph = this->he_.boundaryField()[patchi]; fvPatchScalarField& pheu = this->heu_.boundaryField()[patchi]; fvPatchScalarField& pmu_ = this->mu_.boundaryField()[patchi]; fvPatchScalarField& palpha_ = this->alpha_.boundaryField()[patchi]; if (pT.fixesValue()) { forAll(pT, facei) { const typename MixtureType::thermoType& mixture_ = this->patchFaceMixture(patchi, facei); ph[facei] = mixture_.HE(pp[facei], pT[facei]); ppsi[facei] = mixture_.psi(pp[facei], pT[facei]); pmu_[facei] = mixture_.mu(pp[facei], pT[facei]); palpha_[facei] = mixture_.alphah(pp[facei], pT[facei]); } } else { forAll(pT, facei) { const typename MixtureType::thermoType& mixture_ = this->patchFaceMixture(patchi, facei); pT[facei] = mixture_.THE(ph[facei], pp[facei], pT[facei]); ppsi[facei] = mixture_.psi(pp[facei], pT[facei]); pmu_[facei] = mixture_.mu(pp[facei], pT[facei]); palpha_[facei] = mixture_.alphah(pp[facei], pT[facei]); pTu[facei] = this->patchFaceReactants(patchi, facei) .THE(pheu[facei], pp[facei], pTu[facei]); } } } }``` Now, to calculate the temperature the function THE is called from mixture_ as mixture_.THE After some search I think I have to go to "specie/thermo/thermo/thermoI.H". I'm not sure why. It would be great if you could tell me if I am right? and if yes why?. So in "specie/thermo/thermo/thermoI.H": Code: Code: ``` template class Type> inline Foam::scalar Foam::species::thermo::THE ( const scalar he, const scalar p, const scalar T0 ) const { return Type >::THE(*this, he, p, T0); }``` Code: ``` ``` Then it seems it goes to "specie/thermo/absoluteEnthalpy/absoluteEnthalpy.H". But I don’t know, where it is defined for the code that we are using "absolute enthalpy". With other word, where in the code the energy form will be read? Code: ```//- Temperature from absolute enthalpy // given an initial temperature T0 scalar THE ( const Thermo& thermo, const scalar h, const scalar p, const scalar T0 ) const { return thermo.THa(h, p, T0); }``` Code: and again back to "specie/thermo/thermo/thermoI.H" to calculate THa Code: ```template class Type> inline Foam::scalar Foam::species::thermo::THa ( const scalar ha, const scalar p, const scalar T0 ) const { return T ( ha, p, T0, &thermo::Ha, &thermo::Cp, &thermo::limit ); }``` and from the same h-data the temperature would be calculated using Newton-Raphson method: Code: ```template class Type> inline Foam::scalar Foam::species::thermo::T ( scalar f, scalar p, scalar T0, scalar (thermo::*F)(const scalar, const scalar) const, scalar (thermo::*dFdT)(const scalar, const scalar) const, scalar (thermo::*limit)(const scalar) const ) const { scalar Test = T0; scalar Tnew = T0; scalar Ttol = T0*tol_; int iter = 0; do { Test = Tnew; Tnew = (this->*limit) (Test - ((this->*F)(p, Test) - f)/(this->*dFdT)(p, Test)); if (iter++ > maxIter_) { FatalErrorIn ( "thermo::T(scalar f, scalar T0, " "scalar (thermo::*F)" "(const scalar) const, " "scalar (thermo::*dFdT)" "(const scalar) const, " "scalar (thermo::*limit)" "(const scalar) const" ") const" ) << "Maximum number of iterations exceeded" << abort(FatalError); } } while (mag(Tnew - Test) > Ttol); return Tnew; }``` So, F is probably the absolute Enthalpy as a function of temperature and the constants of Nasa polynomial from "src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H" Code: ```template inline Foam::scalar Foam::janafThermo::ha ( const scalar p, const scalar T ) const { const coeffArray& a = coeffs(T); return this->RR* ( ((((a/5.0*T + a/4.0)*T + a/3.0)*T + a/2.0)*T + a)*T + a ); }``` and dFdT is its derivative with respect to T. So my questions: 1: How are F and dFdT defined? 2: Do we need pressure (p) for absolute enthalpy? 3: Don't we need velocity to calculate the kinetic energy and substract it from the absoulute enthalpy to calculate temperature? 4:Is the Nasa polynomial (see http://combustion.berkeley.edu/gri_m...nasa_plnm.html) written for absolute or sensible enthalpy? Any help will be appreciated. Last edited by sharifi; August 26, 2017 at 04:48.   October 8, 2020, 10:16 #2 New Member   Join Date: Sep 2020 Posts: 11 Rep Power: 2 Did you get the solution for this? Can anyone give a comment on this topic?  Thread Tools Search this Thread Show Printable Version Email this Page Search this Thread: Advanced Search Display Modes Linear Mode Switch to Hybrid Mode Switch to Threaded Mode 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 OffTrackbacks are Off Pingbacks are On Refbacks are On Forum Rules Similar Threads Thread Thread Starter Forum Replies Last Post immortality OpenFOAM Running, Solving & CFD 12 July 11, 2019 18:33 Sethi OpenFOAM Programming & Development 3 March 23, 2017 22:48 bap FLUENT 0 June 22, 2016 11:55 saba1366 CFX 4 February 10, 2013 02:15 astein OpenFOAM Programming & Development 0 June 28, 2010 08:10

All times are GMT -4. The time now is 11:18.