|
[Sponsors] |
Temperature calculation from total enthalpy in OpenFOAM (XiFOAM) |
|
LinkBack | Thread Tools | Search this Thread | Display Modes |
August 25, 2017, 12:10 |
Temperature calculation from total enthalpy in OpenFOAM (XiFOAM)
|
#1 |
New Member
Jim KIT
Join Date: Aug 2012
Location: Germany
Posts: 25
Rep Power: 14 |
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; } Code:
template<class BasicPsiThermo, class MixtureType> void Foam::heheuPsiThermo<BasicPsiThermo, MixtureType>::correct() { if (debug) { Info<< "entering heheuPsiThermo<BasicPsiThermo, MixtureType>::correct()" << endl; } // force the saving of the old-time values this->psi_.oldTime(); calculate(); if (debug) { Info<< "exiting heheuPsiThermo<BasicPsiThermo, MixtureType>::correct()" << endl; } } Code:
template<class BasicPsiThermo, class MixtureType> void Foam::heheuPsiThermo<BasicPsiThermo, MixtureType>::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]); } } } } 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 Thermo, template<class> class Type> inline Foam::scalar Foam::species::thermo<Thermo, Type>::THE ( const scalar he, const scalar p, const scalar T0 ) const { return Type<thermo<Thermo, 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 Thermo, template<class> class Type> inline Foam::scalar Foam::species::thermo<Thermo, Type>::THa ( const scalar ha, const scalar p, const scalar T0 ) const { return T ( ha, p, T0, &thermo<Thermo, Type>::Ha, &thermo<Thermo, Type>::Cp, &thermo<Thermo, Type>::limit ); } and from the same h-data the temperature would be calculated using Newton-Raphson method: Code:
template<class Thermo, template<class> class Type> inline Foam::scalar Foam::species::thermo<Thermo, Type>::T ( scalar f, scalar p, scalar T0, scalar (thermo<Thermo, Type>::*F)(const scalar, const scalar) const, scalar (thermo<Thermo, Type>::*dFdT)(const scalar, const scalar) const, scalar (thermo<Thermo, Type>::*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<Thermo, Type>::T(scalar f, scalar T0, " "scalar (thermo<Thermo, Type>::*F)" "(const scalar) const, " "scalar (thermo<Thermo, Type>::*dFdT)" "(const scalar) const, " "scalar (thermo<Thermo, Type>::*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<class EquationOfState> inline Foam::scalar Foam::janafThermo<EquationOfState>::ha ( const scalar p, const scalar T ) const { const coeffArray& a = coeffs(T); return this->RR* ( ((((a[4]/5.0*T + a[3]/4.0)*T + a[2]/3.0)*T + a[1]/2.0)*T + a[0])*T + a[5] ); } 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 03:48. |
|
October 8, 2020, 09:16 |
|
#2 |
New Member
Join Date: Sep 2020
Posts: 11
Rep Power: 6 |
Did you get the solution for this? Can anyone give a comment on this topic?
|
|
Thread Tools | Search this Thread |
Display Modes | |
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
whats the cause of error? | immortality | OpenFOAM Running, Solving & CFD | 13 | March 24, 2021 07:15 |
Sensible Enthalpy vs. Total enthalpy | Sethi | OpenFOAM Programming & Development | 3 | March 23, 2017 21:48 |
How does fluent calculate the total enthalpy? | bap | FLUENT | 0 | June 22, 2016 10:55 |
Error in run Batch file | saba1366 | CFX | 4 | February 10, 2013 01:15 |
temperature / enthalpy fields depending on type of fvPatchField | astein | OpenFOAM Programming & Development | 0 | June 28, 2010 07:10 |