热物理类的RTS机制
OpenFOAM利用RTS机制来动态创建模型实例，核心思想在于利用C++中静态成员会在类的实例化之前创建，
以buoyantSimpleFoam.C对于生成热物理类的代码的为例，其中
rhoThermo::New(mesh)
作者网址：https://blog.csdn.net/dsfsdffgfd
调用rhoThermo类中静态成员函数:
Foam::autoPtr<Foam::rhoThermo> Foam::rhoThermo::New
(
const fvMesh& mesh,
const word& phaseName
)
{
return basicThermo::New<rhoThermo>(mesh, phaseName);
}
返回了基类basicThermo中名为New的静态模板类，在RTS机制中也被称为Selector,其函数的定义如下:
template<class Thermo>
Foam::autoPtr<Thermo> Foam::basicThermo::New
(
const fvMesh& mesh,
const word& phaseName
)
{
IOdictionary thermoDict
(
IOobject
(
phasePropertyName(dictName, phaseName),
mesh.time().constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE,
false
)
);
auto cstrIter =
lookupThermo<Thermo, typename Thermo::fvMeshConstructorTable>
(
thermoDict,
Thermo::fvMeshConstructorTablePtr_
);
return autoPtr<Thermo>(cstrIter()(mesh, phaseName));
}
template<class Thermo, class Table>
typename Table::iterator Foam::basicThermo::lookupThermo
(
const dictionary& thermoDict,
Table* tablePtr
)
{
if (thermoDict.isDict("thermoType"))
{
const dictionary& thermoTypeDict = thermoDict.subDict("thermoType");
Info<< "Selecting thermodynamics package " << thermoTypeDict << endl;
if (thermoTypeDict.found("properties"))
{
std::initializer_list<const char*> cmptNames
{
"type",
"mixture",
"properties",
"energy"
};
// Construct the name of the thermo package from the components
const word thermoTypeName
(
thermoTypeDict.get<word>("type") + '<'
+ thermoTypeDict.get<word>("mixture") + '<'
+ thermoTypeDict.get<word>("properties") + ','
+ thermoTypeDict.get<word>("energy") + ">>"
);
return lookupThermo<Thermo, Table>
(
thermoTypeDict,
tablePtr,
cmptNames,
thermoTypeName
);
}
else
{
std::initializer_list<const char*> cmptNames
{
"type",
"mixture",
"transport",
"thermo",
"equationOfState",
"specie",
"energy"
};
// Construct the name of the thermo package from the components
const word thermoTypeName
(
thermoTypeDict.get<word>("type") + '<'
+ thermoTypeDict.get<word>("mixture") + '<'
+ thermoTypeDict.get<word>("transport") + '<'
+ thermoTypeDict.get<word>("thermo") + '<'
+ thermoTypeDict.get<word>("equationOfState") + '<'
+ thermoTypeDict.get<word>("specie") + ">>,"
+ thermoTypeDict.get<word>("energy") + ">>>"
);
return lookupThermo<Thermo, Table>
(
thermoTypeDict,
tablePtr,
cmptNames,
thermoTypeName
);
}
}
else
{
const word thermoTypeName(thermoDict.get<word>("thermoType"));
Info<< "Selecting thermodynamics package " << thermoTypeName << endl;
// Table iterator, not const_iterator
auto cstrIter = tablePtr->find(thermoTypeName);
if (!cstrIter.found())
{
FatalIOErrorInLookup
(
thermoDict,
Thermo::typeName,
thermoTypeName,
*tablePtr
) << exit(FatalIOError);
}
return cstrIter;
}
}
首先定义了读写的字典IOdictionary thermoDict，其中dictName为静态变量定义为:const Foam::word Foam::basicThermo::dictName("thermophysicalProperties"); ,下一步是调用了静态模板函数lookupThermo，其模板参数为rhoThermo::fvMeshConstructorTable，也就是返回值rhoThermo::fvMeshConstructorTable::iterator，参数为Thermo::fvMeshConstructorTablePtr_。lookupThermo函数先读取thermophysicalProperties文件，三种字典形式通过上述代码合为一体，例如字典经过读入写成:
hePsiThermo<pureMixture<const<hConst<perfectGas<specie>>,sensibleEnthalpy>>>，或者直接读入这种形式
thermoType hePsiThermo<pureMixture<const<hConst<perfectGas<specie>>,sensibleEnthalpy>>>;，把这个值赋给thermoTypeName后返回rhoThermo::fvMeshConstructorTablePtr_->find(thermoTypeName)。这段代码的HashTable是通过如下宏定义:
declareRunTimeSelectionTable
(
autoPtr,
rhoThermo,
fvMesh,
(const fvMesh& mesh, const word& phaseName),
(mesh, phaseName)
)
#define declareRunTimeSelectionTable(autoPtr,baseType,argNames,argList,parList)\
\
/* Construct from argList function pointer type */ \
typedef autoPtr<baseType> (*argNames##ConstructorPtr)argList; \
\
/* Construct from argList function table type */ \
typedef HashTable<argNames##ConstructorPtr, word, string::hash> \
argNames##ConstructorTable; \
\
/* Construct from argList function pointer table pointer */ \
static argNames##ConstructorTable* argNames##ConstructorTablePtr_; \
\
/* Table constructor called from the table add function */ \
static void construct##argNames##ConstructorTables(); \
\
/* Table destructor called from the table add function destructor */ \
static void destroy##argNames##ConstructorTables(); \
\
/* Class to add constructor from argList to table */ \
template<class baseType##Type> \
class add##argNames##ConstructorToTable \
{ \
public: \
\
static autoPtr<baseType> New argList \
{ \
return autoPtr<baseType>(new baseType##Type parList); \
} \
\
explicit add##argNames##ConstructorToTable \
( \
const word& lookup = baseType##Type::typeName \
) \
{ \
construct##argNames##ConstructorTables(); \
if (!argNames##ConstructorTablePtr_->insert(lookup, New)) \
{ \
std::cerr<< "Duplicate entry " << lookup \
<< " in runtime selection table " << #baseType \
<< std::endl; \
error::safePrintStack(std::cerr); \
} \
} \
\
~add##argNames##ConstructorToTable() \
{ \
destroy##argNames##ConstructorTables(); \
} \
\
add##argNames##ConstructorToTable \
(const add##argNames##ConstructorToTable&) = delete; \
\
void operator= \
(const add##argNames##ConstructorToTable&) = delete; \
}; \
\
/* Class to add constructor from argList to table */ \
/* Remove only the entry (not the table) upon destruction */ \
template<class baseType##Type> \
class addRemovable##argNames##ConstructorToTable \
{ \
public: \
\
const word name; /* Lookup name for later removal */ \
\
static autoPtr<baseType> New argList \
{ \
return autoPtr<baseType>(new baseType##Type parList); \
} \
\
explicit addRemovable##argNames##ConstructorToTable \
( \
const word& lookup = baseType##Type::typeName \
) \
: \
name(lookup) \
{ \
construct##argNames##ConstructorTables(); \
argNames##ConstructorTablePtr_->set(lookup, New); \
} \
\
~addRemovable##argNames##ConstructorToTable() \
{ \
if (argNames##ConstructorTablePtr_) \
{ \
argNames##ConstructorTablePtr_->erase(name); \
} \
} \
\
addRemovable##argNames##ConstructorToTable \
(const addRemovable##argNames##ConstructorToTable&) = delete; \
\
void operator= \
(const addRemovable##argNames##ConstructorToTable&) = delete; \
};
上面代码翻译一下:
/* Construct from argList function pointer type */
typedef autoPtr<rhoThermo> (*fvMeshConstructorPtr)(const fvMesh& mesh, const word& phaseName);
/* Construct from argList function table type */
typedef HashTable<fvMeshConstructorPtr, word, string::hash> fvMeshConstructorTable;
/* Construct from argList function pointer table pointer */
static fvMeshConstructorTable* fvMeshConstructorTablePtr_;
/* Table constructor called from the table add function */
static void constructfvMeshConstructorTables();
/* Table destructor called from the table add function destructor */
static void destroyfvMeshConstructorTables();
/* Class to add constructor from argList to table */
template<class rhoThermoType>
class addfvMeshConstructorToTable
{
public:
static autoPtr<rhoThermo> New (const fvMesh& mesh, const word phaseName)
{
return autoPtr<rhoThermo>(new rhoThermoType (mesh, phaseName));
}
explicit addfvMeshConstructorToTable
(
const word& lookup = rhoThermoType::typeName
)
{
constructfvMeshConstructorTables();
if (!fvMeshConstructorTablePtr_->insert(lookup, New))
{
std::cerr<< "Duplicate entry " << lookup
<< " in runtime selection table " << "rhoThermo"
<< std::endl;
error::safePrintStack(std::cerr);
}
}
~addfvMeshConstructorToTable()
{
destroyfvMeshConstructorTables();
}
addfvMeshConstructorToTable
(const addfvMeshConstructorToTable&) = delete;
void operator=
(const addfvMeshConstructorToTable&) = delete;
};
/* Class to add constructor from argList to table */
/* Remove only the entry (not the table) upon destruction */
template<class rhoThermoType>
class addRemovablefvMeshConstructorToTable
{
public:
const word name; /* Lookup name for later removal */
static autoPtr<rhoThermo> New (const fvMesh& mesh, const word phaseName)
{
return autoPtr<rhoThermo>(new rhoThermoType (mesh, phaseName));
}
explicit addRemovablefvMeshConstructorToTable
(
const word& lookup = rhoThermoType::typeName
)
:
name(lookup)
{
constructfvMeshConstructorTables();
fvMeshConstructorTablePtr_->set(lookup, New);
}
~addRemovablefvMeshConstructorToTable()
{
if (fvMeshConstructorTablePtr_)
{
fvMeshConstructorTablePtr_->erase(name);
}
}
addRemovablefvMeshConstructorToTable
(const addRemovablefvMeshConstructorToTable&) = delete;
void operator=
(const addRemovablefvMeshConstructorToTable&) = delete;
};
上端代码主要是创建了HashTable，热物理类的独特性在于HashTable的元素不是在派生类中添加，而是通过下面宏进行定义:
makeThermos
(
rhoThermo,
heRhoThermo,
pureMixture,
constTransport,
sensibleEnthalpy,
hConstThermo,
perfectGas,
specie
);
这里直接对上面的宏进行翻译得到如下:
// **************************************步骤一*****************************************************
typedef
constTransport
<
species::thermo
<
hConstThermo
<
perfectGas
<
Specie
>
>,
sensibleEnthalpy
>
>
constTransportsensibleEnthalpyhConstThermoperfectGasSpecie
typedef
heRhoThermo
<
rhoThermo,
pureMixture<constTransportsensibleEnthalpyhConstThermoperfectGasSpecie>
>
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie;
// **************************************步骤二*****************************************************
template<>
const ::Foam::word heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie::typeName("heRhoThermo<pureMixture<constTransportsensibleEnthalpyhConstThermoperfectGasSpecie>>")
// **************************************步骤三*****************************************************
addToRunTimeSelectionTable
(
basicThermo,
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie,
fvMesh
);
addToRunTimeSelectionTable
(
basicThermo,
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie,
fvMeshDictPhase
);
addToRunTimeSelectionTable
(
fluidThermo,
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie,
fvMesh
);
addToRunTimeSelectionTable
(
fluidThermo,
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie,
fvMeshDictPhase
);
addToRunTimeSelectionTable
(
rhoThermo,
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie,
fvMesh
);
addToRunTimeSelectionTable
(
rhoThermo,
heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie,
fvMeshDictPhase
);
上面的代码步骤一先用typedef给一定模板参数下constTransport和heRhoThermo的别称，再在步骤二中对这一参数下的模板类特例话其静态变量，最后的步骤三利用addToRunTimeSelectionTable宏命令，在rhoThermo,fluidThermo,basicThermo的HashTable中添加上面参数的模板类，再次把上述宏命令展开:
rhoThermo::addfvMeshConstructorToTable<heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie>
addthisTypeheRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecieConstructorTorhoThermoTable_;
终于，再次回到了rhoThermo中通过宏定义模板类 addfvMeshConstructorToTable ，实例化了一个addfvMeshConstructorToTable<heRhoThermopureMixtureconstTransportsensibleEnthalpyhConstThermoperfectGasSpecie>对象，在其构造函数中，首先调用了constructfvMeshConstructorTables函数，这一函数通过如下的宏进行定义:
#define defineRunTimeSelectionTableConstructor(baseType,argNames) \
\
/* Table constructor called from the table add function */ \
void baseType::construct##argNames##ConstructorTables() \
{ \
static bool constructed = false; \
if (!constructed) \
{ \
constructed = true; \
baseType::argNames##ConstructorTablePtr_ \
= new baseType::argNames##ConstructorTable; \
} \
}
#define defineRunTimeSelectionTableDestructor(baseType,argNames) \
\
/* Table destructor called from the table add function destructor */ \
void baseType::destroy##argNames##ConstructorTables() \
{ \
if (baseType::argNames##ConstructorTablePtr_) \
{ \
delete baseType::argNames##ConstructorTablePtr_; \
baseType::argNames##ConstructorTablePtr_ = nullptr; \
}
}
\
// Create pointer to hash-table of functions
#define defineRunTimeSelectionTablePtr(baseType,argNames) \
\
/* Define the constructor function table */ \
baseType::argNames##ConstructorTable* \
baseType::argNames##ConstructorTablePtr_(nullptr)
展开后写为:
/* Table constructor called from the table add function */
void rhoThermo::constructfvMeshConstructorTables()
{
static bool constructed = false;
if (!constructed)
{
constructed = true;
rhoThermo::fvMeshConstructorTablePtr_
= new rhoThermo::fvMeshConstructorTable;
}
}
void rhoThermo::destroyfvMeshConstructorTables()
{
if (rhoThermo::fvMeshConstructorTablePtr_)
{
delete rhoThermo::fvMeshConstructorTablePtr_;
rhoThermo::fvMeshConstructorTablePtr_ = nullptr;
}
}
rhoThermo::fvMeshConstructorTable* \
rhoThermo::fvMeshConstructorTablePtr_(nullptr)
上面的constructed为局部静态变量，在计算中只会初始化一次，不会重复调用，也就是对于rhoThermo只会生成一个fvMeshConstructorTablePtr_，接下来调用了
fvMeshConstructorTablePtr_->insert(lookup, New)，回到HashTable的定义:
typedef autoPtr<rhoThermo> (*fvMeshConstructorPtr)(const fvMesh& mesh, const word& phaseName);
typedef HashTable<fvMeshConstructorPtr, word, string::hash> fvMeshConstructorTable;
其key为’lookup’即之前特例化的模板名heRhoThermo<pureMixture<constTransportsensibleEnthalpyhConstThermoperfectGasSpecie>>，value为返回值autoPtr<rhoThermo>参数const fvMesh& mesh, const word& phaseName的函数指针，其值为类中的静态函数static autoPtr<rhoThermo> New (const fvMesh& mesh, const word phaseName)，具体调动关系可以详见表驱动模式。层层代入最后返回到buoyantSimple.C的代码就是:
autoPtr<rhoThermo>(autoPtr<rhoThermo>(
new heRhoThermo<rhoThermo,
pureMixture<constTransport<species::thermo<hConstThermo<perfectGas<Specie>>,
sensibleEnthalpy>>>> (mesh, phaseName)))
也就是:
autoPtr<rhoThermo> pThermo(
new heRhoThermo<rhoThermo,
pureMixture<constTransport<species::thermo<hConstThermo<perfectGas<Specie>>,sensibleEnthalpy>>>>
(mesh, null));
热物理类的代码实现
通过热物理类的RTS机制代码，同样以字典为例最终代码得到了指向heRhoThermo<rhoThermo,pureMixture<constTransport<species::thermo<hConstThermo<perfectGas<Specie>>,sensibleEnthalpy>>>>实例的指针，首先看看heRhoThermo的代码实现上的区别:
Appendix I
能量方程可以表示为:
ρ
D
E
D
t
⏞
总能
=
−
∇
⋅
(
p
U
)
⏞
压力功
−
∇
⋅
q
˙
′
′
⏞
热流
+
∇
⋅
(
τ
⋅
U
)
⏞
粘性力功
+
ρ
g
⋅
U
⏞
重力功
+
ρ
q
˙
′
′
′
⏞
热源
(1)
\overbrace{\rho\frac{DE}{Dt}}^{\textbf{总能}}= \overbrace{-\nabla\cdot(p\bold{U})}^{\textbf{压力功}} -\overbrace{\nabla\cdot\bold{\dot{q}''}}^{\textbf{热流}} +\overbrace{\nabla\cdot(\bold{\tau}\cdot\bold{U})}^{\textbf{粘性力功}} +\overbrace{\rho\bold{g}\cdot\bold{U}}^{\textbf{重力功}} +\overbrace{\rho \dot{q}'''}^{\textbf{热源}} \tag{1}
ρDtDE
总能=−∇⋅(pU)
压力功−∇⋅q˙′′
热流+∇⋅(τ⋅U)
粘性力功+ρg⋅U
重力功+ρq˙′′′
热源(1)
(1)为能量方程的非守恒形式，通过随体导数的展开公式
D
D
t
≡
∂
∂
t
+
U
⋅
∇
\frac{D}{Dt}\equiv\frac{\partial}{\partial t}+\bold{U}\cdot\nabla
DtD≡∂t∂+U⋅∇公式及质量守恒可以把并把任一张量
Q
\bold{Q}
Q的非守恒形式展开为守恒形式:
ρ
D
Q
D
t
≡
∂
ρ
Q
∂
t
+
∇
⋅
(
ρ
Q
U
)
\rho\frac{D\bold{Q}}{Dt}\equiv\frac{\partial\rho \bold{Q}}{\partial t}+\nabla\cdot(\rho\bold{Q}\bold{U})
ρDtDQ≡∂t∂ρQ+∇⋅(ρQU)
把(1)写为守恒形式并把总能量 E E E展开为内能 e e e和机械能 K K K之和得到:
∂ ρ e ∂ t + ∇ ⋅ ( ρ e U ) + ∂ ρ K ∂ t + ∇ ⋅ ( ρ K U ) = − ∇ ⋅ ( p U ) − ∇ ⋅ q ˙ ′ ′ + ∇ ⋅ ( τ ⋅ U ) + ρ g ⋅ U + ρ q ˙ ′ ′ ′ (2) \frac{\partial \rho e}{\partial t}+\nabla\cdot(\rho e\bold{U}) +\frac{\partial \rho K}{\partial t}+\nabla\cdot(\rho K\bold{U}) = -\nabla\cdot(p\bold{U}) -\nabla\cdot\bold{\dot{q}''} +\nabla\cdot(\bold{\tau}\cdot\bold{U}) +\rho\bold{g}\cdot\bold{U} +\rho \dot{q}''' \tag{2} ∂t∂ρe+∇⋅(ρeU)+∂t∂ρK+∇⋅(ρKU)=−∇⋅(pU)−∇⋅q˙′′+∇⋅(τ⋅U)+ρg⋅U+ρq˙′′′(2)
公式(2)即为基于内能形式的能量方程。根据焓的定义:
h = e + p ρ h = e+\frac{p}{\rho} h=e+ρp
把:
∂ ρ e ∂ t + ∇ ⋅ ( ρ e U ) = ∂ ρ ( h − p ρ ) ∂ t + ∇ ⋅ ( ρ ( h − p ρ ) U ) = ∂ ρ h ∂ t + ∇ ⋅ ( ρ h U ) − ∂ p ∂ t − ∇ ⋅ ( p U ) \frac{\partial \rho e}{\partial t}+\nabla\cdot(\rho e\bold{U}) = \frac{\partial \rho (h-\frac{p}{\rho})}{\partial t}+\nabla\cdot(\rho (h - \frac{p}{\rho})\bold{U}) =\frac{\partial \rho h}{\partial t}+\nabla\cdot(\rho h\bold{U})- \frac{\partial p}{\partial t} - \nabla\cdot(p\bold{U}) ∂t∂ρe+∇⋅(ρeU)=∂t∂ρ(h−ρp)+∇⋅(ρ(h−ρp)U)=∂t∂ρh+∇⋅(ρhU)−∂t∂p−∇⋅(pU)
代入公式(2)可以得到基于焓形式的能量方程:
∂ ρ h ∂ t + ∇ ⋅ ( ρ h U ) + ∂ ρ K ∂ t + ∇ ⋅ ( ρ K U ) = ∂ p ∂ t − ∇ ⋅ q ˙ ′ ′ + ∇ ⋅ ( τ ⋅ U ) + ρ g ⋅ U + ρ q ˙ ′ ′ ′ (3) \frac{\partial \rho h}{\partial t}+\nabla\cdot(\rho h\bold{U}) +\frac{\partial \rho K}{\partial t}+\nabla\cdot(\rho K\bold{U}) = \frac{\partial p}{\partial t} -\nabla\cdot\bold{\dot{q}''} +\nabla\cdot(\bold{\tau}\cdot\bold{U}) +\rho\bold{g}\cdot\bold{U} +\rho \dot{q}''' \tag{3} ∂t∂ρh+∇⋅(ρhU)+∂t∂ρK+∇⋅(ρKU)=∂t∂p−∇⋅q˙′′+∇⋅(τ⋅U)+ρg⋅U+ρq˙′′′(3)
内能形式与焓形式的能量方程是完全一致的，只是形式不同，因此在理论上两者没有差别。
对动量方程:
∂
ρ
U
∂
t
+
∇
⋅
(
ρ
U
U
)
=
∇
⋅
τ
+
ρ
g
−
∇
p
\frac{\partial \rho \bold{U}}{\partial t}+\nabla\cdot(\rho \bold{U}\bold{U}) =\nabla\cdot\bold{\tau}+\rho\bold{g}-\nabla p
∂t∂ρU+∇⋅(ρUU)=∇⋅τ+ρg−∇p
两边同时点乘速度
U
\bold{U}
U可以得到机械能守恒方程:
∂
ρ
K
∂
t
+
∇
⋅
(
ρ
K
U
)
=
(
∇
⋅
τ
)
⋅
U
+
ρ
g
⋅
U
−
∇
p
⋅
U
(4)
\frac{\partial \rho K}{\partial t}+\nabla\cdot(\rho K\bold{U}) = (\nabla\cdot\bold{\tau})\cdot\bold{U}+\rho\bold{g}\cdot\bold{U}-\nabla p\cdot\bold{U} \tag{4}
∂t∂ρK+∇⋅(ρKU)=(∇⋅τ)⋅U+ρg⋅U−∇p⋅U(4)
把公式(4)分别代入(2)(3)并且去除粘性力影响，内能形式的能量方程的源项为
−
p
(
∇
⋅
U
)
-p(\nabla\cdot\bold{U})
−p(∇⋅U):
∂
ρ
e
∂
t
+
∇
⋅
(
ρ
e
U
)
=
−
p
(
∇
⋅
U
)
−
∇
⋅
q
˙
′
′
+
ρ
q
˙
′
′
′
(5)
\frac{\partial \rho e}{\partial t}+\nabla\cdot(\rho e\bold{U}) = -p(\nabla\cdot\bold{U}) -\nabla\cdot\bold{\dot{q}''} +\rho \dot{q}''' \tag{5}
∂t∂ρe+∇⋅(ρeU)=−p(∇⋅U)−∇⋅q˙′′+ρq˙′′′(5)
而焓形式的能量方程则替换为
∂
p
∂
t
+
U
⋅
(
∇
p
)
\frac{\partial p}{\partial t} +\bold{U}\cdot(\nabla p)
∂t∂p+U⋅(∇p)：
∂ ρ h ∂ t + ∇ ⋅ ( ρ h U ) = U ⋅ ( ∇ p ) + ∂ p ∂ t − ∇ ⋅ q ˙ ′ ′ + ρ g ⋅ U + ρ q ˙ ′ ′ ′ (6) \frac{\partial \rho h}{\partial t}+\nabla\cdot(\rho h\bold{U}) = \bold{U}\cdot(\nabla p) +\frac{\partial p}{\partial t} -\nabla\cdot\bold{\dot{q}''} +\rho\bold{g}\cdot\bold{U} +\rho \dot{q}''' \tag{6} ∂t∂ρh+∇⋅(ρhU)=U⋅(∇p)+∂t∂p−∇⋅q˙′′+ρg⋅U+ρq˙′′′(6)
对比公式(5)(6)，对于不可压流体或者压缩性较小的流体，(5)中的 − p ( ∇ ⋅ U ) -p(\nabla\cdot\bold{U}) −p(∇⋅U)会很小，这在数值中具有收敛性的优势，因此对于 ψ \psi ψ很小的流动状况，可以 使用内能形式的能量方程增强计算的收敛性 ，但收敛后两者的结果没有区别。
另外不同的求解器也会选择忽略(4)(5)中的一些做功项情况，在buoyantSimpleFoam的能量方程中忽略了粘性力的做功项:
fvScalarMatrix EEqn
(
fvm::div(phi, he)
+ (
he.name() == "e"
? fvc::div(phi, volScalarField("Ekp", 0.5*magSqr(U) + p/rho))
: fvc::div(phi, volScalarField("K", 0.5*magSqr(U)))
)
- fvm::laplacian(turbulence->alphaEff(), he)
==
rho*(U&g)
+ radiation->Sh(thermo, he)
+ fvOptions(rho, he)
);
而在rhoSimpleFoam的能量方程中同时忽略了粘性力和重力的做功:
fvScalarMatrix EEqn
(
fvm::div(phi, he)
+ (
he.name() == "e"
? fvc::div(phi, volScalarField("Ekp", 0.5*magSqr(U) + p/rho))
: fvc::div(phi, volScalarField("K", 0.5*magSqr(U)))
)
- fvm::laplacian(turbulence->alphaEff(), he)
==
fvOptions(rho, he)
);
在rhoCentralFoam的能量方程主要考虑公式(4)(5)中的粘性力做功:
volScalarField rhoE
(
IOobject
(
"rhoE",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
rho*(e + 0.5*magSqr(U))
);
surfaceScalarField sigmaDotU
(
"sigmaDotU",
(
fvc::interpolate(muEff)*mesh.magSf()*fvc::snGrad(U)
+ fvc::dotInterpolate(mesh.Sf(), tauMC)
)
& (a_pos*U_pos + a_neg*U_neg)
);
solve
(
fvm::ddt(rhoE)
+ fvc::div(phiEp)
- fvc::div(sigmaDotU)
);
扫一扫
3905
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
