@luofq-sysu 第一个问题:是的,我取的是流向速度时均值。
第二个问题,是的,是不同流向位置。我当时做的是一个平板流动,平板上流向的网格113个。150个不是一个必须的值,看你算例的情况跟需求。我的是二维流动,因此没有展向的问题。
fangyuanaza
帖子
-
-
@luofq-sysu 实现了,我意识到thickness 只能完成模拟后计算得到,因此,我是按照公式用python计算的。基本想法是在我感兴趣的区域监测了150个垂直于壁面的线,监测密度与速度。根据公式 得到每个位置的thickness,然后连接。具体做法是
1: 用python写一个sampleDict(不想手动写150个位置的监测信息)
2: 用python读取每个位置的信息,根据公式计算thickness. -
@李东岳 此问题已解决,做了很多测试后发现,不是残差不好没收敛的问题,而是后处理的问题。本算例计算域在利用对称边界条件后仅为原算例的1/4.而要监测的为duct中心的流动状况。也就是说监测点在计算域边界,因为计算机误差,无法精确取到这个点,计算机会进行插值,导致了这么多spikes。将监测点往计算域内部靠近能解决这个问题。
关于残差,尽管满足y+<1,但是有些监测的物理量的残差不乐观,这个可以加密网格解决。我在截面加密一倍网格后可以解决。但是个人认为没有必要。因为残差不好的物理量都是本身数值很小的,比如压力,Uy, Uz。本身这些量接近于零,因此数值的一点点波动,相对残差会比较大。这些并不影响最终结果。
-
@田畔的风 感谢您的回复~
不是纵横比的问题,调整纵横比之后:
而是因为计算域边界是曲线,当调用tricontourf时,在一个矩形域内进行了插值,导致边界变形,成了直线。同理,用contourf时,首先需要对mesh进行插值,也会有同样的问题。
这个问题已经解决了,对于计算域的两个block, 分开处理。先reshapre到对应的block网格,然后用contourf来画。可以得到物理量的云图。
-
@田畔的风 这个流场图真的很漂亮。请教一下老师,我的网格是O型+H型,Fortran代码算的,能导入paraView。想请教一下,这种情况如何导入到python画图呢?如果直接用contourf 或者tricontourf都会导致网格变形。比如在paraview中为:
在python中直接
ax1.tricontourf(x,y, u,linewidth=0.25, levels=levels,cmap='coolwarm',vmin=0, vmax=1),显示出来为:
造成的原因可能是因为您说的Matplotlib只能处理三角形的网格。想请问在不清楚网格具体信息的情况下,能不能从paraview导出成三角网格?尝试过您代码中的Slice,但由于本身是2D,不能切片。请问有什么建议么? -
@李东岳 看了您的算例设置以后发现添加源项的方式是一样的,您的算例是准直接模拟,所以会有细小的涡。我的就是steady的RANS,所以流场比较平均
-
@李东岳 谢谢李老师的回复!一般圆管临界雷诺数是2000-2600, 这个雷诺数是5600. 所以认为应该是能发展成湍流的如果不用周期边界条件,发现>70h 才发展为完全湍流。从tke的分布判断的。除了这些波动,整体趋势实际上与DNS数据非常接近。速度流场看起来正常。模拟的四分之一个计算域,可以看到转角的流动速度更小:
-
请教各位老师,最近模拟一个矩形管道流动,结果与DNS对比。刚开始采用一般性边界条件,入口匀速流动,出口定压。但是因为流动雷诺数很低,需要很长的区间发展成完全的湍流,但是残差情况较好,如下:
因此后来改成周期性边界条件cyclic。在控制方程中添加动量源项来驱动流动,通过fvOptions添加momentumSource,发现这样残差变差:
体现在模拟结果上为,速度曲线很多波动:
请教是什么原因呢? -
@李东岳 是的,我当时不清楚该怎么做来用python处理整个流场数据,最后是用paraview导入tecplot数据,然后导出成csv等方便读取数据的格式,用python进行处理。因为tecplot数据结构是分块的,而csv等格式是坐标点对应物理量的形式,方便读取。
第二个可能对大家有用的地方在于:有时候tecplot数据导入paraview会一直卡顿,无法完全显示,需要删掉一些符号:
• X, Y and Z array headers MUST NOT have other characters like [mm]. It must be simply “x”, “y” and “z”.Paraview also accepts this: “Z, mm” and this “x/c” • This parameter is not allowed by Paraview (it must be deleted from the file): “ZONETYPE=…” • “STRANDID” and “SOLUTION TIME” are currently unsupported by Paraview but they can be left in the file (it gives only a warning message) • Comments are allowed by using “#” (like Python)我当时是删除了ZONETYPE, 才能导入paraview的。
-
@chszkc 感谢提供这个资源,试一试~
-
@hy1112006 我有点记不清了,现在是可以打开的。可能是路径中的dyfluid-7 中的-不易识别?你可以试一试换一个路径编译,或者改成dyfluid_7试一试
-
请教各位老师,在做几组流场图对比的时候,不用软件导出来的图,而是导出数据,用python进行绘制,改如何做呢?比如tecplot形式生成的数据。
-
是,问题是如何让新装的python和open foam在一个docker 的container里面,按官网教程装好open foam后,要在生成的container里面加上python,用代码添加,这部分不知道怎么处理,请问有经验分享么?感觉docker是好,但是不太用户友好
-
@李东岳 感谢老师的回答,抱歉才看到回信,最终解决办法是强制读取边界信息,检测距离边界最近的网格点数据并给出。黄色震荡部分是因为插值了内部网格点和边界值的结果。
-
请教各位老师,有没有尝试过在Mac 系统上 M1 芯片上跑过python与openfoam. 用python命令自动运行openfoam. 我目前是在电脑上用docker成功安装了openfoam7,能运行算例,也能编译新的模型,接着在电脑上安装了anaconda2 和anaconda3,因为有时候需要调试旧的代码。在终端测试2和3都能用,但是进入open foam的环境下是没法运行anaconda的,根本原因应该是这三个软件不在一个container里面。那么如何在现有的安装了openfoam的container下安装python2&3,能让python脚本运行open foam算例呢?
-
请教各位老师,我设置了一个hump算例,需要监测其表面压力与摩擦力分布。其中一个方法是在paraview中实现,仅选中hump底面,用plotoverintersection 可以给定一个平面从而在相交曲线上plot实现。但是我需要这是一个自动化的过程,希望在sampleDict中实现,因此,我取出intersection的坐标,直接剪测这些点。出现的问题是这两种方法结果不一致,后者曲线震荡,怀疑是给的点不精确,给出的监测结果是边界与内部网格插值结果,这种情况怎么处理呢?请各位老师指教。下图红线是前者方法画的摩擦力系数图,黄色线是用后者方法画的:
-
@fangyuanaza 解决了这个问题,代码分享如下,以免大家有需要:
const volScalarFields& p = this->db().objectRegistry::lookupObject<volScalarField>("p")
经过Info输出发现,确实是每一步的压力场
-
@李东岳 学生也include 了header volFields.H, 应该就是两步,一是#include “volFields.H”, 第二是:
volScalarField p ( IOobject ( "p", this->runTime_.timeName(), this->mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), this->mesh_ );
控制方程也求解了每步的pressure,不清楚问题出在哪里,请老师指点一下
-
@李东岳 有的,计算结果中
Time = 1 GAMG: Solving for p, Initial residual = 1, Final residual = 0.00833422.
从残差控制可见求解了控制方程的压力
-
请教各位老师~ 想在湍流模型中读取压力场,进行运算后置入湍流输运方程中,采取createField.H中的形式,如下:
volScalarField p ( IOobject ( "p", this->runTime_.timeName(), this->mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), this->mesh_ );
可以编译通过,但是运行时候报错:
cannot find file".../../case/1/p".目前是steady case, 1100步保存一次结果。 尝试过用READ_IF_MODIFIED. NO WRITE,都不行。受编程指南启发,尝试加上
dimensionedScalar ( "p", dimensionSet (1,-1,-2,0,0,0,0), lookupObject<volScalarField>("p") )
报错,没有lookupObject. 请问该如何解决这个问题?p应该是public变量,而且确实在不同时间步都储存了。
-
@李东岳 谢谢老师回复~ 离散形式应该是
本质上我是想请教如果对某个方向的量进行操作,当对某个量的网格点进行计算是是取parameter[celli], 若x方向不进行操作,有没有类似于cell(i,j)这种表达,可以仅对j的量进行操作? -
请教各位老师,学生最近在模拟平板流动,想请教有没有在openfoam内部编程得到mmomentum thickness 而不是提取整个场的数据再在其他地方后处理的资料?
-
请教各位老师~
****目标:***希望给雷诺应力增加高阶项来修改湍流模型,高阶项的形式是固定的(标量对称张量)
***想法:***因为雷诺应力本身就是对称张量,其中的标量是随时变化的,机器学习生成的。如果用include **.H的方式来修改湍流模型,则每次需要编译,比较耗时,(因为增加项会不断变化)。所以想通过在constant文件夹中以.H形式存储标量表达式,从而在turbulenceProperties中读取文件,以系数的形式修改表达式,由于系数是实时读取的,因此只需要编译一次。.H文件中就是一个标量表达式,命名为zetai.H
****问题:****如果仅有zetai.H,内容如下(I1 I2是已经在湍流模型中计算好了的体标量场):
并在湍流模型中读取:
如何能将zetai中的内容转换为体标量场进行计算? -
@李东岳 谢谢李老师的回复~ 教师节快乐~ 是一个标准的channel。但是雷诺数非常高,Re=80 million, 没有8000万网格,是NASA官网上的一个标准算例,独特的地方就是雷诺数非常大
-
设置了一个高雷诺数(Re=80million based on channel's height)的管道流动算例,在加密网格过程中,发现残差曲线出现了一些小幅度振荡,用的求解器是simpleFoam,
设置的simple算法收敛要求针对各个量是1e-4,所以1600迭代步就收敛了,但是相比于粗网格或者中等网格的结果,发现残差振荡更厉害,粗网格,中等网格与细网格残差对比如下:
而当使用可压缩求解器时,残差平缓很多,只是e的残差下降很慢。残差如下:
所以想请教大家,稳态不可压算例残差达到要求的标准是什么?除了各个量达到残差要求,有没有对曲线是振荡还是平缓有限制? -
@李东岳 我重新定义了k_, 将其强制变成const,但是还是报错,错误如下:
error: no matching function for call to ‘ddt(const alphaField&, const rhoField&, const volScalarField&)’ + fvc::ddt(alpha, rho, k_) ~~~~~~~~^~~~~~~~~~~~~~~~
可见k_前面已经有const,但是还是不行,可见应该还是要修改alpha与rho的类型?但是我重新赋值定义alpha, rho为volScalarField, 会报初始化错误,不能直接这样改,修改方式如下:
const volScalarField alphat = this->alpha_; 再在fvc中调用fvc::ddt(alphat, rhot, k_) -
@李东岳 谢谢李老师的回复~ 但是学生不太理解,volScalarField是继承于GeometricField的类型,k_是volScalarField类型,应该不冲突把?而且为什么fvm中也是这样定义的,也是传入这三个变量就可以用呢?
-
后续的故事一直在这里自问自答 把里面的量重新定义成带边界的量之后就可以用fvc了(除kSource()与fvOpetions(alpha, rho, k_)以外),现在还有一项没有搞定,请教下各位CFDer, 就是fvc::ddt(alpha, rho, k_)这一项,报错是
error: no matching function for call to ‘ddt(const alphaField&, const rhoField&, Foam::volScalarField&)’ + fvc::ddt(alpha, rho, k_) ~~~~~~~~^~~~~~~~~~~~~~~~
问题是fvc与fvm中此项定义是一样的,都是:
ddt ( const volScalarField& alpha, const volScalarField& rho, const GeometricField<Type, fvPatchField, volMesh>& vf ) { return fv::ddtScheme<Type>::New ( vf.mesh(), vf.mesh().ddtScheme ( "ddt(" + alpha.name() + ',' + rho.name() + ',' + vf.name() + ')' ) ).ref().fvmDdt(alpha, rho, vf); }
虽然说alphg的类型是const alphaField&, rho的类型是rhoField&, 它们俩是geometricOneField 不是 volScalarField,但是为什么fvm可以进行运算但是fvc不可以呢?
-
@fangyuanaza 我知道什么原因了,以tke方程第一项为例,在OpenFOAM V240版本中,
fvc::SuSp(2.0/3.0*rho_*divU, k_)
这里面的量都是volScalarField, 定义rho_与divU的时候也是这样,但是在新版本V7中fvc::SuSp(2.0/3.0*alpha()*rho()*divU, k_)
中,这些量都是内部场的量,没有包含边界数据,改成fvc::SuSp(2.0/3.0*alpha*rho*divU, k_)
,并且divU的定义也从volScalarField::Internel 改成volScalarField, 就可以编译成功。但是这样太麻烦了吧,比如epsilonByk的传入变量也都是volScalarField::Internel形式的,要全部重新定义一边,这样感觉改起来太笨了,请教大家有什么好的改进方式么?从volScalarField::Internal的量,除了重新定义的方法,怎么能改成volScalarField,让显式离散fvc重新能用上,求tke方程残差呢? -
@xpqiu 感谢老师这么详细的指导!我这个算例是二维平板流动的,第一层y+<1,在粘性底层,所以nut为零,应该用nu。所以其实我根本不需要编程计算uTau,因为它可以由壁面监测的壁面摩擦力除以密度然后开根号得到。因为对于u+, y+曲线图来说,uTau在壁面固定地方是固定的。我也对比了这两个值,在边界上编程用nu代替nut,壁面摩擦速度数值是2.5657. 而用监测的壁面摩擦力及密度计算的壁面摩擦速度值是2.5734. 非常接近。再次感谢老师的指导,也感觉到看书的时候认为理解的概念,实际操作发现其实根本没懂,还是要多实践啊。谢谢老师~
-
@xpqiu 请教老师,为什么不写成这样呢:
forAll(uTau_, celli){ uTau_[celli] = sqrt(nut[celli] * mag(tgradU()[celli].yx())); } uadd_ = u/(uTau_ + uTausmall); yadd_ = y_*uTau_/((this->mu()/this->rho_));
此处uTau定义为一个标量场。
是学生对这个公式理解不对吗? -
@xpqiu 理解了,因此我应该取uTau的边界网格进行遍历
forAll(uTau.boundaryFieldRef(),patchi)
或者
forAll(uTau.boundaryField(),patchi)
可以运行成功!但是这里对这种编程思路不太理解:
由李老师提供的公式来看,这里
uTau应该是一个体标量场,这个值应该在近壁面是较大的,因为速度梯度较大,也就是边界层内部值比较大。那么为什么仅给网格边界赋值呢?而且,壁面上nut为零。如果只对边界进行赋值的话,uTau在壁面也是零。我输出uTau对比过了,uTau模拟出来结果在远场上边界才有值,其余地方为零,这是因为壁面nut为零,而场内部初始话赋予的零为初始值。所以学生认为应该直接对整个网格场赋值。不知道理解对不对,请老师指点一下。
-
@xpqiu 谢谢老师~按您这样修改以后,确实没有sqrt的错误了,但是还有错误没有解决:
#0 Foam::error::printStack(Foam::Ostream&) at ??:? #1 Foam::sigSegv::sigHandler(int) at ??:? #2 ? in "/lib/x86_64-linux-gnu/libc.so.6" #3 Foam::RASModels::kOmegaSSTyuPlus<Foam::EddyDiffusivity<Foam::ThermalDiffusivity<Foam::CompressibleTurbulenceModel<Foam::fluidThermo> > > >::correct() at ??:? #4 ? in "/home/dyfluid/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/rhoSimpleFoam" #5 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6" #6 ? in "/home/dyfluid/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/rhoSimpleFoam" Segmentation fault (core dumped)
#3应该是最关键的,但是学生不太理解的是编译的时候仅仅是在基本的kOmegaSST湍流模型修改的BasicTurbulencemodel,这里为什么会出现fluidThermo名称空间下的错误?
-
@李东岳 请教下李老师~对于这里计算uTau,我用dimensionedScalar 是可以运算的,但是用这块代码计算uTau会报错,可能理解不太对。对于forAll中的×××,认为是uTau, 在这之前,我先定义了uTau是一个体标量场,初始值给的0,您代码中的u.boundary中的u认为是U.component(0)(), 这样代码能编译成功,但是运行几步会报错,请老师指教一下,万分感激(Ps 我在湍流模型里面实现这几个量,离壁面距离已经给出了,就是y_):
volScalarField u(U.component(0)()); volScalarField uTau(0*sqrt(k_)); //dimensionedScalar uTau("uTau", dimVelocity, 1.003); dimensionedScalar uTausmall( "0", dimensionSet(0, 1, -1, 0, 0, 0, 0), 1e-20 ); forAll(uTau,patchi) { uTau.boundaryFieldRef()[patchi] = sqrt(nut.boundaryField()[patchi]*u.boundaryField()[patchi].snGrad()); } uadd_ = u/(uTau + uTausmall); yadd_ = y_*uTau/((this->mu()/this->rho_));
运行报错为:
#0 Foam::error::printStack(Foam::Ostream&) at ??:? #1 Foam::sigFpe::sigHandler(int) at ??:? #2 ? in "/lib/x86_64-linux-gnu/libc.so.6" #3 Foam::sqrt(Foam::Field<double>&, Foam::UList<double> const&) at ??:? #4 Foam::sqrt(Foam::tmp<Foam::Field<double> > const&) at ??:? #5 Foam::RASModels::kOmegaSSTyuPlus<Foam::EddyDiffusivity<Foam::ThermalDiffusivity<Foam::CompressibleTurbulenceModel<Foam::fluidThermo> > > >::correct() at ??:? #6 ? in "/home/dyfluid/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/rhoSimpleFoam" #7 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6" #8 ? in "/home/dyfluid/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/rhoSimpleFoam" Floating point exception (core dumped)
说是这个里面的sqrt传入变量不对?
-
请教一下各位CFDer, 很多算例的设置是无量纲化的,那么有没有热物性参数是怎么无量纲化的介绍呢?或者比如有理想状态下空气热物性参数的无量纲化结果?因为不同的无量纲化过程可能结果是有差异的,比如说压力,有的用rhoU^2来无量纲化,有的用0.5rho*U^2来无量纲化,导致压力相差一倍。
-
@chengan-wang 抱歉,才看到消息。我建议跟着18.04的步骤尝试一下?我对比了不同版本Ubuntu上OF2.4.0的安装,是差不多的, 不同的是有些的用的gcc版本,有些用gcc 5 有的用4.7等,因为用新的Ubuntu装以前的OF,主要问题就是gcc不兼容,所以试试跟着18.04装一下?可以有问题再讨论。只是个小建议哈,毕竟我也没这么安装过,建议尝试的话之前也保存一下原来的系统。
-
@李东岳 这是我在V 2.4.0上的做法,思路是一样的,就是方程变显示,左右相减,只是V7中湍动能方程更复杂一点,以下是V2.4.0编译成功的写法:
-
@李东岳 李老师,不好意思,是我没有说清楚。我之前没有用internal,直接定义velScalarField residual = ..., 计算完之后我再update boundary的值。我在V 2.4.0也是这么做的,这里加Internal我也是测试看是不是加了能成功,但是没有,后来也没去掉了,上传图片上保留了。刚刚我也去掉测试了,还是同上面一样的错误。
-
@fangyuanaza 后面保留fvc也会出现跟前面一样的错误 no matching function for call .....
还测试过仿照李老师之前的建议改成k_(), 也是不行 -
@李东岳 李老师好~ 这样的话((2.0/3.0)*alpha()*rho()*divU, k_()) + (alpha()*rho()*epsilonByk(F1, F23), k_) 不报错啦,但是后面的ddt(alphat, rho, k_) 去掉fvc::,就不能用ddt来进行对时间的求导了,后面的div(alphaRhoPhi, k_)也会报错:
-
@李东岳 谢谢李老师一早就回复学生! 改成k_()还是报一样的错误:
-
请教各位老师~学生想看看湍动能方程左右两边的差值,于是将fvm改为fvc,就不是隐式方程了,而是直接显式求差值。这种方法在openfoamV 2.4 0上成功了,但是在openfoamV7上会报错,代码如下:
报错如下,显示SuSp中(2.0/3.0)*alpha()*rho()*divU与 k_类型不一致?我开始以为是SuSp的问题,因为这个代表根据rho的正负决定用隐式还是显式,但是改成Sp, 还是同样的报错,请各位老师指点下迷津~
-
@bestucan 谢谢回复~ 主要是2%的条件放在原来的网格上又能算,原来的网格跟现在的网格只是相当于平移了一小段距离而已。
-
@李东岳 是的~编译成功了! 给李老师打call!
-
请教各种大神,想通过X方向速度的大小限制一个对称张量nonlinearStress的作用范围,满足小于39.86或者大于41.8,就将这个张量设置成零。如下:
报错说不应该写成U.component(0)[celli], 我用mag(U.component(0))或者直接用U.component(0)都不行。那应该怎么写呢?我理解的U.component(0)已经是标量了。可以写成k_[celli],为什么不能写成U.component(0)[celli]?烦请各位老师指点一下 -
@xpqiu 您是对的,是我犯低级错误了,我定义的小量small_val在这个要修改的式子后面,放到前面就可以了。谢谢
-
@xpqiu 其他行代码也调试完了,就剩这个地方了,是我读取信息有问题么?
aij_ ( IOobject ( "aij", runTime_.timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ ),
我是这样读取的,我用Info是可以打印出来aij各个分量的值的
-
@xpqiu
问题就出现在这里,我把这个注释掉就能用 -
@xpqiu 我还以为这里是在说类型不一致的问题,我是一行行排查错误的,把其他的改动都注释掉了,现在排查定位到了这一行。
后处理得到momentum thickness
后处理得到momentum thickness
考虑动量源项后残差变差
python进行OpenFOAM流场后处理
python进行OpenFOAM流场后处理
考虑动量源项后残差变差
考虑动量源项后残差变差
考虑动量源项后残差变差
python进行OpenFOAM流场后处理
python进行OpenFOAM流场后处理
paraFOAM时的dlopen error
python进行OpenFOAM流场后处理
如何在docker上同时运行python和open foam
如何监测边界与平面相交线的参数
如何在docker上同时运行python和open foam
如何监测边界与平面相交线的参数
如何在湍流模型中读取压力场
如何在湍流模型中读取压力场
如何在湍流模型中读取压力场
如何在湍流模型中读取压力场
后处理得到momentum thickness
后处理得到momentum thickness
以读取系数的方式修改湍流模型
如何判断残差曲线是否达到要求?
如何判断残差曲线是否达到要求?
新旧版本编程差异问题
新旧版本编程差异问题
新旧版本编程差异问题
新旧版本编程差异问题
yPlus在openfoam代码里面的实现
yPlus在openfoam代码里面的实现
yPlus在openfoam代码里面的实现
yPlus在openfoam代码里面的实现
yPlus在openfoam代码里面的实现
无量纲算例中热物性参数设置
分享Ubuntu上同时安装openfoam2.4.0与7
分享Ubuntu上同时安装openfoam2.4.0与7
新旧版本编程差异问题
新旧版本编程差异问题
新旧版本编程差异问题
新旧版本编程差异问题
新旧版本编程差异问题
新旧版本编程差异问题
遇到个诡异的问题
速度分量写成条件判断的问题
速度分量写成条件判断的问题
场计算类型的问题
场计算类型的问题
场计算类型的问题
场计算类型的问题