View Issue Details

Jump to NotesJump to History
IDProjectCategoryView StatusDate SubmittedLast Update
0000669OpenFOAMBugpublic2012-10-26 17:052012-10-29 12:37
Reporteruser528Assigned Tohenry  
PrioritynormalSeverityminorReproducibilityalways
Status resolvedResolutionfixed 
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.

Activities

henry

2012-10-29 12:37

manager   ~0001757

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

Issue History

Date Modified Username Field Change
2012-10-26 17:05 user528 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