q-DNS计算槽道流遇到了一些问题,求大神们指点
-
我想用q-DNS模拟一个Re_tau=180的槽道湍流,网格与这篇链接文本文献中基本相同。
我更改了一下icoFoam的计算公式,代码包:0_1546247060398_pGradIcoFoam.zip
具体改法如下:
- createFields.H中,添加一个pGrad场,量纲与压力梯度相同,大小直接给出:
volVectorField pGrad ( IOobject ( "pGrad", runTime.timeName(), mesh, IOobject::NO_READ, //MUST_READ, IOobject::AUTO_WRITE ), mesh, dimensionedVector("pGrad",dimensionSet(0, 1, -2, 0, 0),vector(-1,0,0)) );
- 动量预测方程中把pGrad项加入到方程右边,与压力梯度项并列
{ solve(UEqn == -fvc::grad(p)-pGrad); }
-
压力泊松方程中也加入pGrad项带来的影响
fvScalarMatrix pEqn ( fvm::laplacian(rAU, p) + fvc::div(rAU\*pGrad) == fvc::div(phiHbyA) );
-
还有最后一步的速度修正:
U = HbyA - rAU\*fvc::grad(p)- rAU\*pGrad;
- 另外,还有按照苏老师博客里的方法,把icoFoam改为可以变时间步长的求解器
- 加入计算并输出平均速度的代码:
scalar magUbarAve = 0.0; scalar V = 0.0; const scalarField& cv = mesh.V(); vector flowDir_ = vector(1,0,0); for(label celli=0; celli < mesh.C().size() ; celli++) { scalar volCell = cv[celli]; V += volCell; magUbarAve += (flowDir_ & U[celli])*volCell; } reduce(magUbarAve, sumOp<scalar>()); reduce(V, sumOp<scalar>()); magUbarAve /= V; Info << "Ubar = " << magUbarAve << endl;
计算算例:
0_1546249122420_channel180.zip
计算初始场用perturbU工具进行打乱,然后计算一些时间步获得完全发展的湍流场
我现在遇到的问题是,计算一段时间后,发现靠近边界层的速度算得还比较准,但是槽道中部远离边界层的地方出现奇怪的脉动
放大可看到这个脉动大概以三个cells为一个周期。
最奇怪的是平均速度型竟然上下不对称已做时间足够长的平均
我的改动中应该没有会引起这种不对称的地方?
壁面律图像和雷诺应力图像也同样,壁面处算得还好,中间算不好
似乎可能是中间网格不够密的原因?但是我使用的网格应该和那篇参考文献中相同啊,他也用的是openFoam做的q-DNS,但他得到的结果就很好,所以我现在不知道错在什么地方。
我修改后的求解器、算例的基本设定都已上传,希望各位老师能帮忙指点一下,不胜感激! -
@东岳 谢谢东岳老师
这个meanVelocityForce的源项我用过,但是和我想要的不太一样,meanVelocityForce也是通过添加一个压力梯度来控制平均速度的,但它给添加的压力梯度项是随时间步不断计算调整的,我想添加一个简单的常数的压力梯度。不过我确实可以通过修改meanVelocityForce来达到这个目的,今后我会试一下这个。
另外,我那个平均速度的代码就是从meanVelocityForce中copy过来的,稍微改了一下下让他能在求解器中运行,现在来看运行的还可以。而且这个平均速度只是为了在log中输出一个当前流场的平均速度,方便我检测求解状况,不会影响迭代求解的,所以还是不太清楚结果中的问题出在了哪里。ORZ -
@xiexiaoyang 有道理
-
@XieXiaoyang 您好,我想请教一下,像这种 DNS / LES 槽道流的结果,我想画速度剖面,进行正则化时候的摩擦速度${u_{\tau}}$和粘性长度${\delta_{\nu}=\nu/u_{\tau}}$是怎样计算的?
我当前的做法是:使用postChannel(就像channel395那个算例一样)将空间时间都做平均之后,可以得到${u(y)}$这样一个分布,然后计算壁面上的速度曲线的斜率,计算壁面摩擦应力,计算摩擦速度:
\begin{equation}
{u_\tau } = \sqrt {\frac{{{\tau _w}}}{\rho }} {\rm{ = }}\sqrt {\nu \frac{{\partial u}}{{\partial y}}}
\end{equation}
但算出来的摩擦速度和预设的相差很大。所以最后的曲线是这样的
算了很久,对数率那里也无法和理论解对上。计算结果的整个曲线也与理论解的形状不符合。
感觉我使用的方法不太对,想请问一下您的做法。 -
-
@SSSSK
两种方法估计摩擦应力:
第一种是经验公式:
$${{\mathop {\rm Re} \nolimits} _\tau } \approx 0.09 {{\mathop {\rm Re} \nolimits} ^ {0.88}}$$
$${\mathop{\rm Re}\nolimits} = {U_m}\left( {2h} \right)/\nu $$
$h$为槽道半高,$U_m$为槽道平均速度。
第二种是完全发展槽道流中的计算:
$$ \frac{{\partial p}}{{\partial x}} \times h = {\tau _w} $$