thermophysicalModels LIB in OpenFOAM

This post analyze classes about species, mixture and thermo in thermophysicalModels LIB.

UML 类图

下图是热物理库部分类的 UML 类图 。欲获取更佳阅读体验,可以右键在新标签页打开该图。

求解器中创建对象

solver 中的 creatFields.H 中:

autoPtr<psiReactionThermo> pThermo(psiReactionThermo::New(mesh));
psiReactionThermo& thermo = pThermo();

创建 thermo,在 reactingFoam 自带算例中,创建的是 hePsiThermo 的对象,但是这里 thermo 是基类 psiReactionThermo 的引用,即只能使用 psiReactionThermo 这一继承链的方法。

basicSpecieMixture& composition = thermo.composition();

创建composition,是指系统中的混合物。可以调用继承链中 mixture 这一分支的方法。

燃烧类的求解器中,一般都没有对单个组分进行操作,这些操作都在热物理库里边进行。

multiComponentMixture 类中的私有数据 speciesData_ 的类型是 PtrList<sutherlandTransport<species::thermo<janafThermo<perfectGas<specie>>,sensibleEnthalpy>>>,它是组分的 list,可以用它获取组分的所有信息。solver 中可以通过 multiComponentMixturesolvercomposition 的类型是 basicSpecieMixture,需要强制转换成 multiComponentMixture 才行 ) 的 speciesData() 或者 getLocalThermo(speciei) 这两个函数来调用。如:

dynamic_cast<const multiComponentMixture<gasHThermoPhysics>&>(composition).speciesData()[speciei].rho(p,T)

也可以通过 multiComponentMixture 的派生类 reactingMixture 来调用。

dynamic_cast<const reactingMixture<gasHThermoPhysics>&>(composition).speciesData()[speciei].rho(p,T)

也可以通过派生类的派生类 SpecieMixture 来调用。

dynamic_cast<const SpecieMixture<reactingMixture<gasHThermoPhysics>>&>(composition).speciesData()[speciei].rho(p,T)

上述操作说明可以使用组分列表调用单个组分的函数(即 sutherlandTransport<species::thermo<janafThermo<perfectGas<specie>>,sensibleEnthalpy>> 继承链中的任何函数),如 speciesData_[i].rho(p, T),这里的 rho(p,T) 来自于 perfectGas 类。

也可以在求解器中重新构造组分列表:

autoPtr<psiReactionThermo> pThermo(psiReactionThermo::New(mesh));
psiReactionThermo& thermo = pThermo();

basicSpecieMixture& composition = thermo.composition();
PtrList<volScalarField>& Y = composition.Y();

PtrList<gasHThermoPhysics> specieData(Y.size());
forAll(specieData, i)
{
//对于每个组分,通过 thermo 调用 Lib 中的 speciesData_ 来构造 solver 中的 specieData
specieData.set
(
i,
new gasHThermoPhysics
(
dynamic_cast<const reactingMixture<gasHThermoPhysics>&>
(thermo).speciesData()[i]
)
);
}

// 每个网格构造一个 mixture,这个 mixture 相当于 multiComponentMixture 中的私有数据 mixture_。
// 它是混合物等效成的一个组分,和单个组分同一类型,即 `gasHThermoPhysics`。
forAll(T, celli)
{
gasHThermoPhysics mixture(0.*specieData[0]);
forAll(Y, speciei)
{
mixture += Y[speciei][celli]*specieData[speciei];
}
}

gasHThermoPhysicssutherlandTransport<species::thermo<janafThermo<perfectGas<specie>>,sensibleEnthalpy>> 的别名,它是单个组分的类型。
事实上,multiComponentMixture 类中的函数 cellMixture 函数可以用于获取 mixture_

template<class ThermoType>
const ThermoType& Foam::multiComponentMixture<ThermoType>::cellMixture
(
const label celli
) const
{
mixture_ = Y_[0][celli]*speciesData_[0];

for (label n=1; n<Y_.size(); n++)
{
mixture_ += Y_[n][celli]*speciesData_[n];
}

return mixture_;
}

使用方法参考 hePsiThermo.C:

const scalarField& hCells = this->he_;
const scalarField& pCells = this->p_;

scalarField& TCells = this->T_.primitiveFieldRef();
scalarField& psiCells = this->psi_.primitiveFieldRef();
scalarField& muCells = this->mu_.primitiveFieldRef();
scalarField& alphaCells = this->alpha_.primitiveFieldRef();

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]);
}

RTS 选择热物理模型

从字典文件中的设置到最终的对象构建

下面介绍我们是如何从 thermophysicalProperties 字典文件中的下列设置构建对象的:

thermoType
{
type hePsiThermo;
mixture reactingMixture;
transport sutherland;
thermo janaf;
energy sensibleEnthalpy;
equationOfState perfectGas;
specie specie;
}

首先需要一个函数读取这些信息:

basicThermoTemplates.C 中的 lookupThermo 函数:

template<class Thermo, class Table>
typename Table::iterator Foam::basicThermo::lookupThermo
(
const dictionary& thermoDict,
Table* tablePtr //即 fvMeshConstructorTablePtr_
)
{
// Construct the name of the thermo package from the components
const word thermoTypeName
(
word(thermoTypeDict.lookup("type")) + '<'
+ word(thermoTypeDict.lookup("mixture")) + '<'
+ word(thermoTypeDict.lookup("transport")) + '<'
+ word(thermoTypeDict.lookup("thermo")) + '<'
+ word(thermoTypeDict.lookup("equationOfState")) + '<'
+ word(thermoTypeDict.lookup("specie")) + ">>,"
+ word(thermoTypeDict.lookup("energy")) + ">>>"
);
//通过字典文件读入的数据,构造 thermoTypeName: //hePsiThermo<reactingMixture<sutherland<janaf<perfectGas<specie>>,sensibleEnthalpy>>>

return lookupThermo<Thermo, Table>
(
thermoTypeDict,
tablePtr,
nCmpt,
cmptNames,
thermoTypeName
);
}

最后返回值是调用 5 个参数的 lookupThermo 函数:

template<class Thermo, class Table>
typename Table::iterator Foam::basicThermo::lookupThermo
(
const dictionary& thermoTypeDict,
Table* tablePtr,
const int nCmpt,
const char* cmptNames[],
const word& thermoTypeName
)
{
// Lookup the thermo package
typename Table::iterator cstrIter = tablePtr->find(thermoTypeName);
return cstrIter;
}
//上述步骤是根据字典文件在备选库中选择我们需要的组合,那么备选库是怎么形成的呢?

第一个 lookupThermo 函数的作用是,将thermophysicalProperties 字典文件中的设置翻译成:

hePsiThermo<reactingMixture<sutherland<janaf<perfectGas<specie>>,sensibleEnthalpy>>>

但该字符串并不是最终的 class,比如 janaf,实际上类名是 janafThermo

第二个 lookupThermo 函数的作用就是,根据该字符串,在 RTS 库中找到我们选择的类。

找到之后,在下面这个 New 函数中新建对象:

template<class Thermo>
Foam::autoPtr<Thermo> Foam::basicThermo::New(const fvMesh& mesh,const word& phaseName)
{
typename Thermo::fvMeshConstructorTable::iterator cstrIter =
lookupThermo<Thermo, typename Thermo::fvMeshConstructorTable>
(
thermoDict,
Thermo::fvMeshConstructorTablePtr_
);
//fvMeshConstructorTablePtr_
//这个是在 runTimeSelectionTables.H 的 declareRunTimeSelectionTable 函数中定义的。
return autoPtr<Thermo>(cstrIter()(mesh, phaseName));
}

那么上面提到的字符串,怎么和 RTS 库中的组合匹配的呢?下面来看 RTS 库的生成过程。

按照Giskard的博客提到的 RTS 流程:

  1. 基类类体里调用 TypeNamedeclareRunTimeSelectionTable 两个函数,类体外面调用 defineTypeNameAndDebugdefineRunTimeSelectionTableaddToRunTimeSelectionTable 三个函数;
  2. 基类中需要一个静态 New 函数作为 selector
  3. 派生类类体中需要调用 TypeName 函数,类体外调用 defineRunTimeSelectionTableaddToRunTimeSelectionTable 两个宏函数。

makeReactionThermo.H 中有两个个重要函数:defineThermoPhysicsReactionThermomakeThermoPhysicsReactionThermos

下面一一介绍:

defineThermoPhysicsReactionThermo

其中的 defineThermoPhysicsReactionThermo, 这个函数通过一组形参列表得到 CThermo##Mixture##ThermoPhys,然后再转换成一串简化的字符串,用于 RTS 查询。这里的形参列表就是实例化所需要提供的各种类。(BaseReactionThermo,CThermo,Mixture,ThermoPhys)对应(psiReactionThermo,hePsiThermo,reactingMixture,sutherlandTransport)

#define defineThermoPhysicsReactionThermo(BaseReactionThermo,CThermo,Mixture,ThermoPhys) \
typedef CThermo//hePsiThermo \
< \
BaseReactionThermo,//psiReactionThermo \
SpecieMixture \
< \
Mixture//reactingMixture \
< \
ThermoPhys//sutherlandTransport \
> \
> \
> CThermo##Mixture##ThermoPhys; \
\
defineTemplateTypeNameAndDebugWithName \
( \
CThermo##Mixture##ThermoPhys, \
(#CThermo"<" + Mixture<ThermoPhys>::typeName() + ">").c_str(), \
0 \
)
//作用:把尖括号表示的类变成字符串 CThermo##Mixture##ThermoPhys,然后再变成简化的尖括号。

如:

hePsiThermo
<
psiReactionThermo,
SpecieMixture
<
reactingMixture
<
sutherlandTransport
<
species::thermo
<
janafThermo<perfectGas<specie>>,sensibleEnthalpy
>
>
>
>
>

变成 RTS 库中保存的字符串:
hePsiThermo<reactingMixture<sutherland<janaf<perfectGas<specie>>,sensibleEnthalpy>>>
隐藏了一些信息,如 psiReactionThermoSpecieMixture

具体过程:
typeName() 依次调用 reactingMixturesutherlandTransport, species::thermo, janafThermo, sensibleEnthalpy, perfectGas, specietypeName()。返回值依次是:

"reactingMixture<" + ThermoType::typeName() + '>'
sutherland<" + Thermo::typeName() + '>'
Thermo::typeName() + ',' + Type<thermo<Thermo, Type>>::typeName()
"janaf<" + EquationOfState::typeName() + '>'
"sensibleEnthalpy"
"perfectGas<" + word(Specie::typeName_()) + '>'

最底层的 specie 类只有一个 ClassName("specie");
我们找到一个文件className.H 中:

#define ClassName(TypeNameString)                             \
ClassNameNoDebug(TypeNameString); \
static int debug

#define ClassNameNoDebug(TypeNameString) \
static const char* typeName_() { return TypeNameString; } \
static const ::Foam::word typeName

makeThermoPhysicsReactionThermos

调用了上边的 defineThermoPhysicsReactionThermo 用于创建 RTS 中的 typeName,调用了 addThermoPhysicsThermo(位于 makeThermo.H), 用于写入到 RTS 库。

#define makeThermoPhysicsReactionThermos(BaseThermo,BaseReactionThermo,CThermo,Mixture,ThermoPhys)  \                                                                          
defineThermoPhysicsReactionThermo(BaseReactionThermo,CThermo,Mixture,ThermoPhys); \
\
addThermoPhysicsThermo(basicThermo, CThermo##Mixture##ThermoPhys); \
addThermoPhysicsThermo(fluidThermo, CThermo##Mixture##ThermoPhys); \
addThermoPhysicsThermo(BaseThermo, CThermo##Mixture##ThermoPhys); \
addThermoPhysicsThermo(BaseReactionThermo, CThermo##Mixture##ThermoPhys)

函数的形参:(BaseThermo,BaseReactionThermo,CThermo,Mixture,ThermoPhys),这是实例化需要提供的所有信息,比如其中的一种组合:(psiThermo,psiReactionThermo,hePsiThermo,reactingMixture,gasHThermoPhysics)(所有可能的组合在 psiReactionThermos.C 文件中列出)。

#define addThermoPhysicsThermo(BaseThermo,CThermoMixtureThermoPhys) //定义于makeThermo.H  
addToRunTimeSelectionTable \
( \
BaseThermo, \
CThermoMixtureThermoPhys, \
fvMesh \
);

这里到底往 RTS Tableadd 了什么?
BaseThermopsiThermoCThermoMixtureThermoPhyshePsiThermo##reactingMixture##gasHThermoPhysics

那前边提到的保存在 RTS 中的简化的尖括号字符串是干嘛的?

hePsiThermo<reactingMixture<sutherland<janaf<perfectGas<specie>>,sensibleEnthalpy>>>

其实这两个的含义是一样的,因为

typedef
sutherlandTransport
<
species::thermo
<
janafThermo
<
perfectGas<specie>
>,
sensibleEnthalpy
>
> gasHThermoPhysics;

typedef sutherlandTransport<species::thermo<janafThermo<perfectGas<specie>>,sensibleEnthalpy>> gasHThermoPhysics;

但是此二者在 RTS 里边有什么联系呢?

Author: Yan Zhang
Link: https://openfoam.top/thermodynamicLIB/
Copyright Notice: All articles in this blog are licensed under CC BY-NC-SA 4.0 unless stating additionally.
微信打赏