hPolynomialThermo.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "hPolynomialThermo.H"
27 #include "IOstreams.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class EquationOfState, int PolySize>
33 (
34  const word& name,
35  const dictionary& dict
36 )
37 :
38  EquationOfState(name, dict),
39  hf_
40  (
41  dict
42  .subDict("thermodynamics")
43  .lookupBackwardsCompatible<scalar>({"hf", "Hf"})
44  ),
45  sf_
46  (
47  dict
48  .subDict("thermodynamics")
49  .lookupBackwardsCompatible<scalar>({"sf", "Sf"})
50  ),
51  CpCoeffs_
52  (
53  dict.subDict("thermodynamics").lookup
54  (
55  "CpCoeffs<" + Foam::name(PolySize) + '>'
56  )
57  ),
58  hCoeffs_(),
59  sCoeffs_()
60 {
61  hCoeffs_ = CpCoeffs_.integral();
62  sCoeffs_ = CpCoeffs_.integralMinus1();
63 
64  // Offset h poly so that it is relative to the enthalpy at Tstd
65  hCoeffs_[0] += hf_ - hCoeffs_.value(Tstd);
66 
67  // Offset s poly so that it is relative to the entropy at Tstd
68  sCoeffs_[0] += sf_ - sCoeffs_.value(Tstd);
69 }
70 
71 
72 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
73 
74 template<class EquationOfState, int PolySize>
76 (
77  Ostream& os
78 ) const
79 {
81 
82  dictionary dict("thermodynamics");
83  dict.add("hf", hf_);
84  dict.add("sf", sf_);
85  dict.add
86  (
87  word("CpCoeffs<" + Foam::name(PolySize) + '>'),
88  CpCoeffs_
89  );
90  os << indent << dict.dictName() << dict;
91 }
92 
93 
94 // * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
95 
96 template<class EquationOfState, int PolySize>
97 Foam::Ostream& Foam::operator<<
98 (
99  Ostream& os,
101 )
102 {
103  pt.write(os);
104  return os;
105 }
106 
107 
108 // ************************************************************************* //
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
virtual Ostream & write(const char)=0
Write character.
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:123
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:843
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:1014
Enthalpy based thermodynamics package using a polynomial function of temperature for the constant hea...
hPolynomialThermo(const EquationOfState &pt, const scalar hf, const scalar sf, const Polynomial< PolySize > &CpCoeffs, const typename Polynomial< PolySize >::intPolyType &hCoeffs, const Polynomial< PolySize > &sCoeffs)
Construct from components.
void write(Ostream &os) const
Write to Ostream.
A class for handling words, derived from string.
Definition: word.H:62
const dimensionedScalar Tstd
Standard temperature.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:227
dictionary dict