原有的湍流模型加上非线性项雷诺应力的问题
-
各位有经验的老师们,我想问一下加入修正雷诺应力项的修正湍流模型的问题。目前,我已经尝试了模仿ShihQuadraticKE湍流模型写法,重新写了一个新的湍流模型,具体修改的地方如下:
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // void kEpsilonTBNN::correctNut() { correctNonlinearStress(fvc::grad(U_)); } void kEpsilonTBNN::correctNonlinearStress(const volTensorField& gradU) { timeScale_=k_/epsilon_; // Linear (nut) nut_ = g1_*k_*timeScale_; nut_.correctBoundaryConditions(); // nolinear(tau_NL) volSymmTensorField S(timeScale_*symm(gradU)); volTensorField W(timeScale_*skew(gradU)); nonlinearStress_ = 2*k_ *( g2_ * twoSymm(S&W) + g3_ * dev(innerSqr(S)) + g4_ * dev(symm(W&W)) + g5_ * twoSymm(W&innerSqr(S)) + g6_ * dev(twoSymm(W&W&S)) + g7_ * twoSymm(W&S&W&W) + g8_ * twoSymm(S&W&innerSqr(S)) + g9_ * dev(twoSymm(innerSqr(S)&W&W)) + g10_ * twoSymm(W&innerSqr(S)&W&W) ); }
这是我主要的修正项,即Aij
其中,g1到g10是我通过机器学习出来的数据,
g(i)跟T(i)都是有关Sij跟Ri的公式。目前,按照这种方法计算是发散的,其主要的原因是g(i)跟T(i)是随每次迭代求解出的U,k,epsilon发生变化的,但是只通过读取0文件下的g(i)场只是一个不变的标量,从而导致求解发散,因此需要将学习出来的神经网络结构耦合到OpenFOAM的环境中。因此我打算修改simpleFoam,linearViscousStrees和新加的reynoldsNet神经网络结构。
-
目前,是重新编写了simpleFoam,可以是在计算迭代中同时生成Sij跟Rij数据,为后面的linearViscousStrees提供相应数据,同时修改了动量方程的有关雷诺应力项。
tmp<fvVectorMatrix> tUEqn ( fvm::div(phi, U) + MRF.DDt(U) + turbulence->divDevRhoReff(U, S, R) == fvOptions(U)
我参照了之前的有关文献,同时修改了incompressible文件下的,incompressibleTurbulenceModel.H。
//- Return the source term for the momentum equation using NN surrogate virtual tmp<fvVectorMatrix> divDevReff ( volVectorField& U, volTensorField& S, volTensorField& R ) { }
incompressible/IncompressibleTurbulenceModel文件夹下的,IncompressibleTurbulenceModel.C和.H文件,添加了以下内容
template<class TransportModel> Foam::tmp<Foam::fvVectorMatrix> Foam::IncompressibleTurbulenceModel<TransportModel>::divDevReff ( volVectorField& U, volTensorField& S, volTensorField& R ) const { return divDevRhoReff(U, S, R); }
//- Return the source term for the momentum equation virtual tmp<fvVectorMatrix> divDevReff ( volVectorField& U, volTensorField& S, volTensorField& R ) const;
并修改了linearViscousStrees.C和.H文件,添加了以下内容
template<class BasicTurbulenceModel> Foam::tmp<Foam::fvVectorMatrix> Foam::linearViscousStress<BasicTurbulenceModel>::divDevRhoReff ( volVectorField& U, volTensorField& S, volTensorField& R ) const { NN结构:a_dd } return ( - fvm::laplacian(this->alpha_*this->rho_*this->nu(), U) - fvc::div((this->alpha_*this->rho_*this->nu())*dev2(T(fvc::grad(U))) - this->a_dd) );
a_dd就是我要加入的非线性项,里面的内容主要是链接了reynoldsNet神经网络。
-
@李东岳 多谢李老师对我目前所做的项目有个极高的评价,由于我已经是研三了,这是我课题的最后一部分,“联合培养”项目开展时间有点晚了,多谢李老师的好意。对于我目前所做的东西,感觉李老师理解稍有偏差。g(i),也就是g1到g10是一个有关gardU,k,epsilon的函数,目前无法写出其显式表达,因此打算用神经网络结构来表示,所以按之前步骤只放到0文件夹下面就是一个不变量,迭代一次就应该更新了,每一次OF求解迭代出来的g(i)应该是不一样的。之前我是想跟python耦合的,即求出一步就用学习完成的神经网络更新一次g(i),但是对于高网络数量的模型计算就十分缓慢,因此我才打算将在python的神经网络结构转换为c++能识别的神经网络结构,并重新编写linearViscousStrees和IncompressibleTurbulenceModel来进行OF耦合。目前,修改的代码写的差不多了,就是不知linearViscousStrees和IncompressibleTurbulenceModel的有效编译步骤该如何进行,还需请教一下李老师。
-
@SHUKK
你可以参考一下这个项目:https://github.com/furstj/myTurbulenceModels里面有非线性湍流模型的例子,比如这个 EARSM 模型,https://github.com/furstj/myTurbulenceModels/tree/main/turbulenceModels/RAS/EARSM
你可以在这个 EARSM 模型的基础上来修改,基本上只需要其中的 correctNonlinearStress 函数就可以了。
-
@xpqiu 非常感谢,您的建议。我之前也是修改kepsilon湍流模型的,跟EARSM 模型一样,加了correctNonlinearStress函数。之前算发散的时候,我从g1到g10一个一个添加进去计算,发现只带入g1时候才能计算,并且得出来的效果与kOmegaSST的结果相似的,但是只要带入g2就开始发散了。后来,看之前做的过程时候发现,g(i),也就是g1到g10是一个有关gardU,k,epsilon的函数,目前无法写出其显式表达,应该用外接神经网络结构来计算更新。那为什么只带g1能算呢,因为g1*T1是等于OF的原有线性项假设,相当于没有改变其假设,只是修改了Cmu的值而已。
EARSM 模型的非线性项中,我看它的非线性项写法与我的类似,但是它的bate(i)是有明确的显式公式的,是能随着求解器迭代来计算更新的,而我的是隐式公式,目前没有明确公式。volScalarField beta1 = -N*(2.0*sqr(N) - 7.0*IIW) / Q; volScalarField beta3 = -12.0 * IV / (N * Q); volScalarField beta4 = -2.0 * (sqr(N) - 2.0*IIW) / Q; volScalarField beta6 = -6.0 * N / Q; volScalarField beta9 = 6.0 / Q; volScalarField Cmu = - 0.5 * (beta1 + IIW * beta6); this->nut_ = Cmu * this->k_ * tau; this->nut_.correctBoundaryConditions(); this->nonlinearStress_ = this->k_ * symm( beta3 * ( (W & W) - (1.0/3.0) * IIW * I ) + beta4 * ( (S & W) - (W & S) ) + beta6 * ( (S & W & W) + (W & W & S) - IIW * S - (2.0/3.0) * IV * I) + beta9 * ( (W & S & W & W) - (W & W & S & W) )
现在我要做的内容就让我的g(i)跟它bate(i)一样跟随迭代,为了达到我这个效果,不能单修改湍流模型,还要外加上我的c++神经网络来更新g(i)。目前就是编译的一些内容和步骤有点疑惑,需要请教。
-
@李东岳 在 原有的湍流模型加上非线性项雷诺应力的问题 中说:
每个步长首先通过pytorch获取g(i),g(i)是一些列的标量
把g(i)的数值传输到OpenFOAM的非线性湍流模型里
每个时间步封闭,速度压力耦合求解
进入到下一个时间步,回归到第一步所以这是机器学习+CFD的一种方法。但是因为每一步都要训练出来一个g(i),所以速度慢。
一个有关g(i)隐式的公式,若能被OF识别并有效输出,那我就不用每次迭代计算调用pytorch来进行计算
你们打算找出来g(i)的一个公式,写成$g(i)=f(\bfU)$类似的函数,然后流程就变成了:
- 最开始训练出$g(i)=f(\bfU)$函数表达式
- 进行CFD求解(在这一步里面,不需要花时间训练找参数)
我这样理解对么?
-
@李东岳 李老师,我遇到一个wmake的问题。我把我训练完成的神经网络写成一个reynoldsNet库,打算在linearViscousStrees中调用,但是显示了无法识别的问题这是什么导致的?
wmakeLnIncludeAll: running wmakeLnInclude on dependent libraries: unknown option: '-I/home/user3/OpenFOAM/user3-9/reynoldsNet/lnInclude' Usage: wmakeLnInclude [OPTION] [dir] options: -update | -u update -silent | -s use 'silent' mode (do not echo command) -help | -h print the usage Link all the source files in the <dir> into <dir>/lnInclude Note The '-u' option forces an update when the lnInclude directory already exists and changes the default linking from 'ln -s' to 'ln -sf'.
-
@李东岳 上面的步骤是在编译MomentumTransportModels文件夹时候的。同时,我在simpleNNFoam中修改了我的需要连接的linearViscousStrees中的divDevTau(U, S, R),S跟R都已经在 createFields.H中创建了。
// Momentum predictor MRF.correctBoundaryVelocity(U); tmp<fvVectorMatrix> tUEqn ( fvm::div(phi, U) + MRF.DDt(U) + turbulence->divDevTau(U, S, R) // Data driven R-S + viscous component == fvModels.source(U) ); fvVectorMatrix& UEqn = tUEqn.ref(); UEqn.relax(); fvConstraints.constrain(UEqn); if (simple.momentumPredictor()) { solve(UEqn == -fvc::grad(p)); fvConstraints.constrain(U); }
但是我在wmake我的simpleNNFoam中出现了,数据类型不对的问题。
In file included from /home/user3/libtorch/include/c10/util/Exception.h:7, from /home/user3/libtorch/include/ATen/core/Generator.h:11, from /home/user3/libtorch/include/ATen/CPUGeneratorImpl.h:3, from /home/user3/libtorch/include/ATen/Context.h:3, from /home/user3/libtorch/include/ATen/ATen.h:7, from /home/user3/libtorch/include/torch/csrc/api/include/torch/types.h:3, from /home/user3/libtorch/include/torch/script.h:3, from /home/user3/OpenFOAM/user3-9/reynoldsNet/lnInclude/reynoldsNet.H:42, from /home/user3/OpenFOAM/user3-9/src/MomentumTransportModels/momentumTransportModels/lnInclude/linearViscousStress.H:38, from /home/user3/OpenFOAM/user3-9/src/MomentumTransportModels/momentumTransportModels/lnInclude/Stokes.H:39, from /home/user3/OpenFOAM/user3-9/src/MomentumTransportModels/momentumTransportModels/lnInclude/laminarModel.C:27, from /home/user3/OpenFOAM/user3-9/src/MomentumTransportModels/momentumTransportModels/lnInclude/laminarModel.H:193, from /home/user3/OpenFOAM/user3-9/src/MomentumTransportModels/incompressible/lnInclude/kinematicMomentumTransportModel.H:47, from simpleNNFoam.C:66: /home/user3/libtorch/include/c10/util/variant.h:1204:11: error: no matching function for call to 'Foam::data::data(const c10::detail_::impl<c10::SmallVector<c10::SymInt, 5>, at::Tensor>&)' 1204 | data(lib::forward<V>(v)), | ^~~~~~~~~~~~~~~~~~~~~~~~
我查看了 Foam::data::data(const objectRegistry& obr),要的是objectRegistry类,但是在我的reynoldsNet库中已经将libtorch输出的torch::Tensor转换成std::vector<float>然后带入输出,没有将是torch的数据直接带入OF中。这是有数据没按要求来转换吗?
-