在输出时间步的同时,输出颗粒所受drag force文件
-
RT, 好像现有版本的OP没有相应的utility可以实现这一功能吧。是不是只能在求解器的creatFields.H文件里里添加了?
老实说本人的C++功底太差了。。所以想麻烦各位帮我举个例子怎么修改?比如我使用的是ErgunWenYuDrag,想在每个输出的时间步中添加这个颗粒的ErgunWenYuDragForce数据文件,请问怎么实现?
我试着添加如下,但是报错一大堆。
volScalarField ErgunWenYuDrag ( IOobject ( IOobject::groupName("ErgunWenYuDragForce", continuousPhaseName), runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh, );
-
说实话, 我也不是太清楚, 但是你可以参考一下src/lagrangian/basic 中的cloudIO.H和IOPosition, 我觉得.
-
@qingdong_wang @东岳 谢谢一楼和东岳老师的回复。还想请问东岳老师,这个曳力项,应该写在哪个文件里输出?
我发现求解器(比如DPMFoam)里面的creatFields.H文件,输出的都是流体相的场文件,比如Uc, p, phic之类的,并没有颗粒相的场文件。我也查找了一楼所说的src/lagrangian/basic目录下的文件,也没发现控制输出颗粒场的文件。
所以还想麻烦您进一步指明,这个曳力项应该写在哪里?
-
@zhangxc0223 老铁实现了吗?如何不通过createFields添加新变量来输出变量?或者怎么添加不需要read的变量来输出?我试了no_read还是提示没有对应文件
-
如果你要欧拉场,可以看
cloudVolSuSp
怎么输出的,如果你要拉格朗日场,去KinematicParcelIO.C
文件中,添加你的dragForce:参考下面都是要输出的信息
template<class ParcelType> template<class CloudType> void Foam::KinematicParcel<ParcelType>::readFields(CloudType& c) { bool write = c.size(); ParcelType::readFields(c); IOField<label> active ( c.fieldIOobject("active", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, active); IOField<label> typeId ( c.fieldIOobject("typeId", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, typeId); IOField<scalar> nParticle ( c.fieldIOobject("nParticle", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, nParticle); IOField<scalar> d ( c.fieldIOobject("d", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, d); IOField<scalar> dTarget ( c.fieldIOobject("dTarget", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, dTarget); IOField<vector> U ( c.fieldIOobject("U", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, U); IOField<scalar> rho ( c.fieldIOobject("rho", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, rho); IOField<scalar> age ( c.fieldIOobject("age", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, age); IOField<scalar> tTurb ( c.fieldIOobject("tTurb", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, tTurb); IOField<vector> UTurb ( c.fieldIOobject("UTurb", IOobject::MUST_READ), write ); c.checkFieldIOobject(c, UTurb); label i = 0; forAllIter(typename CloudType, c, iter) { KinematicParcel<ParcelType>& p = iter(); p.active_ = active[i]; p.typeId_ = typeId[i]; p.nParticle_ = nParticle[i]; p.d_ = d[i]; p.dTarget_ = dTarget[i]; p.U_ = U[i]; p.rho_ = rho[i]; p.age_ = age[i]; p.tTurb_ = tTurb[i]; p.UTurb_ = UTurb[i]; i++; } }
-
这个题我会
第一步,在ParticleForce类中添加虚函数virtual forceSuSp calcDragForce ( const typename CloudType::parcelType& p, const typename CloudType::parcelType::trackingData& td, const scalar dt, const scalar mass, const scalar Re, const scalar muc ) const { return forceSuSp(Zero); }
第二步,在ErgunWenYuDragForce类添加此函数,并返回
... { return this->calcCoupled ( p, td, dt, mass, Re, muc );
第三步,在ParticleForceList类中添加函数
forceSuSp calcDragForce ( const typename CloudType::parcelType& p, const typename CloudType::parcelType::trackingData& td, const scalar dt, const scalar mass, const scalar Re, const scalar muc ) const { forceSuSp value(Zero); if (calcCoupled_) { forAll(*this, i) { value += this->operator[](i).calcDragForce(p, td, dt, mass, Re, muc); } } return value; }
当然也可以这么写
... { forceSuSp value(Zero); if (calcCoupled_) { forAll(*this, i) { if(isA<ErgunWenYuDragForce>(this->operator[](i))) value = this->operator[](i).calcDragForce(p, td, dt, mass, Re, muc); } } return value; }
这么写的话,前面两步就可以不用做了
第四步,在KinematicParcel类的calcVelocity函数中添加template<class ParcelType> template<class TrackCloudType> const Foam::vector Foam::KinematicParcel<ParcelType>::calcVelocity ( TrackCloudType& cloud, trackingData& td, const scalar dt, const scalar Re, const scalar mu, const scalar mass, const vector& Su, vector& dUTrans, vector& drag, //添加返回参数 scalar& Spu ) const { ... const forceSuSp FDrag = forces.calcDragForce(p, ttd, dt, mass, Re, mu); ... drag = -(FDrag.Sp()*td.Uc() + FDrag.Su()); ... }
在calc函数中添加
template<class ParcelType> template<class TrackCloudType> void Foam::KinematicParcel<ParcelType>::calc ( TrackCloudType& cloud, trackingData& td, const scalar dt ) { ... vector drag = Zero; // Calculate new particle velocity this->U_ = calcVelocity(cloud, td, dt, Re, td.muc(), mass0, Su, dUTrans, darg, Spu); if (cloud.solution().coupled()) { // Update momentum transfer cloud.UTrans()[this->cell()] += np0*dUTrans; cloud.dragForce()[this->cell()] += np0*drag; // Update momentum transfer coefficient cloud.UCoeff()[this->cell()] += np0*Spu; }
第五步,在KinematicCloud类中添加变量
autoPtr<volVectorField::Internal> drag_;
构造中添加
drag_ ( new volVectorField::Internal ( IOobject ( this->name() + ":drag", this->db().time().timeName(), this->db(), IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE ), mesh_, dimensionedVector(dimForce, Zero) ) ),
以及函数
template<class CloudType> inline Foam::DimensionedField<Foam::vector, Foam::volMesh>& Foam::KinematicCloud<CloudType>::drag() { return *drag_; } template<class CloudType> inline const Foam::DimensionedField<Foam::vector, Foam::volMesh>& Foam::KinematicCloud<CloudType>::drag() const { return *drag_; }
第六步,主程序添加
auto cloudDrag(kinematicCloud.drag()); volVectorField cloudVolDrag ( IOobject ( "cloudVolDrag", runTime.timeName(), mesh ), mesh, dimensionedVector(cloudDrag.dimensions()/dimVolume, Zero), zeroGradientFvPatchVectorField::typeName ); cloudVolDrag.primitiveFieldRef() = -cloudDrag/mesh.V(); cloudVolSUSu.correctBoundaryConditions();
大概是这样,没编译过,新鲜的。里面的量纲和曳力应该还需要调整。
...
其实大部分情况下,一般都只会用到曳力,所以可以直接用
kinematicCloud.SU(Uc)
就行了,那这样的话,前面所有的都不用做了 -
非常感谢大佬们的对话,受益匪浅啊