2016-10-28 11:12 BST

View Issue Details
IDProjectCategoryView StatusLast Update
0000669OpenFOAM[All Projects] Bugpublic2012-10-29 12:37
Reporterpkalinin
Assigned Tohenry
PrioritynormalSeverityminorReproducibilityalways
StatusresolvedResolutionfixed
Product Version
Target VersionFixed in Version
Summary0000669: Wrong calculation of calculating Cp in addition operator in hConstThermo
DescriptionIt seems to me that operator+ in hConstThermo calculates the Cp in a wrong way (and this applies also to Hf and to other arithmetical operators).

Consider the following code:

int main(int argc, char *argv[])
{
typedef constTransport<specieThermo<hConstThermo<perfectGas> > > ThermoType;

IFstream f("thermo.dict");
dictionary dict(f);

ThermoType t1(dict.subDict("specie1"));
ThermoType t2(dict.subDict("specie2"));

Info << "W= " << t1.W() << " " << t2.W() << " " << (t1+t2).W() << endl;
Info << "Cp= " << t1.cp(1) << " " << t2.cp(1) << " " << (t1+t2).cp(1) << endl;

return(0);
}
Run it on the following thermo.dict file:

specie1{
specie {
nMoles 1;
molWeight 1;
}
thermodynamics {
Cp 1;
Hf 0;
}
transport {
mu 1;
Pr 1;
}
}
specie2{
specie {
nMoles 1;
molWeight 0.5;
}
thermodynamics {
Cp 2;
Hf 0;
}
transport {
mu 1;
Pr 1;
}
}
The program will output
W= 1 0.5 0.75
Cp= 1 1 1.125
However, the I expect the following output:
W= 1 0.5 0.75
Cp= 1 1 1

Indeed, here two species are defined having different molar masses (1 and 0.5 kg/kmol) and different heat capacitance per kg (1 and 2), but note that they have the same heat capacitance per mole --- 1 J/(kmol K).

The program reads both specie data, adds them, and outputs the resulting heat capacitance. As expected, cp() for both species separately is 1 (cp() returns J/(kmol K)). The heat capacitance per kilomole of mixture of these two species should therefore be 1 too, independently of the mixture composition, therefore the program should output 1, not 1.125.

-----

As I understand it, the problems comes from the hConstThermo addition operator that calculates the thermal capacity by the following code (hConstThermoI.H, lines 209-210):

ct1.nMoles()/eofs.nMoles()*ct1.Cp_
+ ct2.nMoles()/eofs.nMoles()*ct2.Cp_
This would be correct if hConstThermo::Cp_ was the thermal capacity per mole (J/(kmol K) or J/(mol K)). However, Cp_ is the thermal capacity per kg (J/(kg K)) --- see, for example, hConstThermo::cp(), which returns J/(kmol K) and is implemented as Cp_*this->W().

So, the correct code should be something like

(ct1.nMoles()/eofs.nMoles()*ct1.Cp_*ct1.W()
+ ct2.nMoles()/eofs.nMoles()*ct2.Cp_*ct2.W())
/eofs.W()
The same applies to Hf calculations and to other arithmetical operators in hConstThermo (operator-, operator+=, operator -=).

Note also that a similar code exists in janafThermo:
highCpCoeffs[coefLabel] =
molr1*jt1.highCpCoeffs_[coefLabel]
+ molr2*jt2.highCpCoeffs_[coefLabel];
However, here highCpCpeffs_ are per mole, not per kg (and janafThermo::cp() does not refer to W() ), so this equation is a correct one.
TagsNo tags attached.
Attached Files

 Relationships
 Relationships

 Notes ~0001757 henry (manager) 2012-10-29 12:37 Thanks for the detailed bug-report. On further consideration we have concluded that it would be better to convert the Cp_ and Hf_ to mole-based on input and back again on output in the same manner as in hPolynomialThermo which is more consistent and also resolves the issue you reported. Resolved by commit 74eb38ab125c177165c9704f405bcbfe12cf26c8
 Notes

 Date Modified Username Field Change Issue History 2012-10-26 17:05 pkalinin New Issue 2012-10-29 12:37 henry Note Added: 0001757 2012-10-29 12:37 henry Status new => resolved 2012-10-29 12:37 henry Resolution open => fixed 2012-10-29 12:37 henry Assigned To => henry
 Issue History