一个关于 Axisymmetric Subsonic Jet 的问题

King_RVM

有个关于rhoPimpleFoam求解亚音速射流的问题求助CFDer。我刚入手OpenFoam，最近在做一个亚音速射流的问题（出口流速大概在165m/s），采用轴对称方法计算（旋转角度1deg，采用wedge边界条件），需要给定总压力入口边界条件（totalPressure）。参考实验来自NASA官网：https://turbmodels.larc.nasa.gov/jetsubsonic_val.html
算例的网格和边界如下：
批注 2020-06-27 103534.jpg
我刚开始采用稳态的求解可压缩问题的rhoSimpleFoam求解器，射流的速度云图结果还不错，如下图：

现在我想采用瞬态的rhoPimpleFoam求解器，网格和边界条件设置均相同，但是得到的速度云图看起来是错误的，没有出现明显的射流界面，时间步1e-6s计算稳定，在0.85s时的结果如下图：

想请问CFDer有没有做过类似计算，发生这个错误的可能原因有些什么呢？我想多尝试尝试改一下，十分感谢！！！

李东岳

瞬态用的是RANS么

King_RVM

@东岳是的李老师，采用RANS。还有就是网格是采取的NASA官网提供的比较粗糙的网格，所以时间步稍微大一些（1e-6s）。细网格对时间步（目前1e-9以下）要求太高了，算的比较慢，目前算了0.01s，发现云图和粗网格的类似，所以现在认为原因不太可能是网格影响的。

李东岳

你用rhoPimpleFoam同时在controlDict里面计算一下meanU之类的数据看看呢？

King_RVM

@东岳谢谢李老师，我监测了物理时间1.06s的Umean，如下图，这个云图也不正确，和rhoSimpleFoam结果相差太大。
批注 2020-06-27 154350.jpg

King_RVM

@东岳我现在主要怀疑会不会是边界条件造成的，但是边界条件和rhoSimpleFoam是一样的，所以挺疑惑的。

李东岳

网格多少万？

King_RVM

@东岳网格数是4520，目前用的是NASA提供的比较粗糙的网格，如下图：批注 2020-06-27 173033.jpg
不过刚刚Umean的结果是在已有结果的基础上，加上Umean后跑了0.05s的结果，这样的结果可能是不正确的，因为没有统计之前时间的速度。我刚刚重新从0开始统计Umean的结果，晚些时候我把结果发到论坛上给您看一下。

附：关于Umean的监测我采用的是知乎上的一个方法：https://zhuanlan.zhihu.com/p/52667661
在controlDict里加入的代码如下：

functions
{
    fieldAverage1
    {
        type            fieldAverage;
        functionObjectLibs ( "libfieldFunctionObjects.so" );
        enabled         true;
        outputControl   outputTime;

        fields
        (
            U
            {
                mean        on;
                prime2Mean  on;
                base        time;
            }
        );
    }
}

在中文网上也看到了关于这个问题的讨论：https://www.cfd-china.com/topic/1453/想在程序中调用umean这个值-应该如何实现
我还想问一下李老师，这个Umean是单个时间步内的所有迭代的平均，还是计算从0s开始到当前时间的一个总的速度平均？如果是后者的话，采用瞬态方法时，如果刚开始计算得到的速度值变化较大的话，对Umean的影响不会很大吗？

King_RVM

@东岳李老师，这是我重新从0s开始计算得到的Umean的云图，当前时间步为0.1s，如下图：

李东岳

Umean是从0开始的值，一般要去掉开始一段时间步，从后面开始平均

网格不多，可以的话，把你的rhoSimpleFoam+rhoPimpleFoam都传上来我给你看看

King_RVM

@东岳太感谢了，李老师，，文件稍微有点大，我就用坚果云了，这是坚果云链接：
https://www.jianguoyun.com/p/DUyH968Qr-LHBxjQ5aYD
https://www.jianguoyun.com/p/DRvCxh0Qr-LHBxjV5aYD

李东岳

这两个网格的boundary名字为什么不一样？为什么不用同样的网格？

King_RVM

@东岳很不好意思李老师，这个算例是我和另一个同学一块讨论着来做的，之前上传的是我俩用相同网格，分别用rhoSimpleFoam和rhoPimpleFoam来计算的结果，所以边界命名不一样。
晚上我重新用自己的命名方式，用两种求解器重新跑了一遍，网格图如下：

rhoSimpleFoam（10w步迭代结果）和rhoPimpleFoam（0.08s时刻结果）结果分别如下图2和图3所示：

这是我重新跑了一遍的算例文件，存在了坚果云中。再次感谢！
https://www.jianguoyun.com/p/DdnpeEMQ1IPKCBiFqqcD
https://www.jianguoyun.com/p/DdA5ZSwQ1IPKCBiGqqcD

cccrrryyy

@King_RVM rhoSimpleFoam和参考结果相比怎么样？我一般用rhoPimpleFoam和unsteady RANS都给很大的时间步让它快速到稳态解，最终的这个稳态解和rhoSimpleFoam类似。你这个不知道从方法来说有没有问题，比如说每一个时间步收敛了么？unsteady RANS每个时间步需要不少的迭代步，你用的是多少？

你可以用flowRatePatch测试一下入口和出口的流量是否还守恒，rhoSimpleFoam应该没有问题。rhoPimpleFoam我的经验是每一个时间步都是守恒的，如果不守恒那说明有问题。

King_RVM

@cccrrryyy 谢谢您的回复。rhoSimpleFoam和实验结果对比还可以，我刚处理了一下这个粗糙网格的结果，图1是沿轴向的速度U分布对比，图2和图3是数值与试验的速度云图对比（图2是U的对比，图3是Umean的对比），可以看出Umean的云图更为光滑一些，图2速度U的云图能看到出口处速度不光滑有波动，现在我们怀疑是出口处壁面边界引起的非物理的泻涡造成的，但是图3Umean的出口处速度云图就很光滑：

对于您的第二个问题，因为我刚开始学OpenFoam，还没找到内迭代步在哪儿里设置，这是我的controlDict文件，里面没有内迭代步的设置，这个我再找找资料看看怎么设置。

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     rhoPimpleFoam;

startFrom       latestTime;

startTime       0;

stopAt          endTime;

endTime         2.0;

deltaT          1e-6;

writeControl    adjustableRunTime;

writeInterval   0.005;

purgeWrite      0;

writeFormat     ascii;

writePrecision   8;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;

functions
{
    #includeFunc MachNo
    #includeFunc residuals
    #includeFunc probes
	#includeFunc R
}

functions
{
    fieldAverage1
    {
        type            fieldAverage;
        functionObjectLibs ( "libfieldFunctionObjects.so" );
        enabled         true;
        outputControl   outputTime;

        fields
        (
            U
            {
                mean        on;
                prime2Mean  on;
                base        time;
            }
        );
    }
}

// ************************************************************************* //

第三个问题我也尝试一下，然后再来讨论。

King_RVM

@cccrrryyy 你好，关于第二个建议，我再fvSolution文件中设置了最大迭代步maxIter为100，计算结果没有什么变化，fvSolution中设置的代码如下：

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    p
    {
        solver          GAMG;
        smoother        GaussSeidel;
        tolerance       1e-6;
        relTol          0.01;
		maxIter 		100;
    }

    pFinal
    {
        $p;
        relTol          0;
    }

    "(rho|U|k|omega|e)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
		maxIter 		100;
    }

    "(rho|U|k|omega|e)Final"
    {
        $U;
        relTol          0;
    }
}

关于第三个建议，我监测了3个入口和1个出口的flowRatePatch数据，两者的总和基本上保持为0，流量也是守恒的。

李东岳

进口omega为什么是0.11？这是计算出来的吧？

李东岳

https://www.jianguoyun.com/p/DUG7u1MQ9s3ZBhi8g6cD

你处理一下这个数据，出一个沿轴向的速度U分布对比，看看是否有改善？

cccrrryyy

@King_RVM 不是这个，是PIMPLE里面的nOuterCorrectors，本来设的是多少？

King_RVM

@东岳是的李老师，这个是计算出来的，公式如下：

CodeCogsEqn (1).gif
其中k取的是0.001，这个值比较小是因为当入口湍动能值较大时，得到的速度值偏大，改小后结果有所改善（但是相比大的k值，小的k值残差振荡比较严重）。turbulent viscosity 取得是0.001，这个是参考的NASA官网提供的值。

King_RVM

@东岳李老师，我做了一个rhoSimpleFoam 10w步的结果与您提供的13000的结果的速度云图对比，因为没监测Ux（轴向速度），目前只能提供一个瞬时的速度Umagnitude的对比结果，如下图：
图片1.jpg
13000这个结果有所改善，但同样在大概2倍的出口半径处仍然存在小的波动。李老师可否提供postProcessing中的监测点数据，这样我可以画个沿轴向的速度曲线来对比。下图是之前我算rhoSimpleFoam的速度曲线，也可以看出在射流出口处速度波动较大，和云图表现出的结果是一致（我们猜测会不会是射流喷嘴处壁面的影响）。
图片2.jpg

King_RVM

@cccrrryyy 不好意思，我刚开始学还没弄明白这些值。pimple里这个的设置如下：

SIMPLE
{
    residualControl
    {
        p               1e-4;
        U               1e-4;
        "(k|omega|e)"   1e-4;
    }

    nNonOrthogonalCorrectors 0;
    pMinFactor      0.1;
    pMaxFactor      2;
}

PIMPLE
{
    nCorrectors              2;
    nNonOrthogonalCorrectors 1;
    nOuterCorrectors         1;
    pMinFactor      0.1;
    pMaxFactor      2;
}

relaxationFactors
{
    equations
    {
        ".*"     1;
    }
}

nOuterCorrectors默认的是1，这个算例是我根据OpenFoam自带的tutorials-rhoPimpleFoam 中的aerofoilNACA0012的算例修改的，因此这个值就没变。

李东岳

你可以直接把13000的数据放到你的算例下面，后处理数据再出图
https://www.jianguoyun.com/p/DTTBLxIQ9s3ZBhjxi6cD 或者你下载我配置后的rhoSimpleFoam的算例文件出个结果图我看看这几天我更新在rhoSimpleFoam解析里面
rhoPimpleFoam那个我感觉是湍流模型的问题，算的太慢了不好调试，你试试调整下湍流边界条件，看看有无改善

李东岳

捕获.JPG

这个是rhoPimpleFoam算的，结果只出到0.02秒，不过看起来也行了。我更新在rhoPimpleFoam解析，然后你去下载一下玩玩

等你的plot图啊

King_RVM

@东岳太感谢了李老师，我用您的配置算例，然后监测了singleGraph一条线上的速度值，U的15000步迭代后的结果如图1和2所示，射流喷嘴出口处仍存在小的波动，但整体计算结果相比NASA的数值结果更靠近实验数据（NASA的实验数据是在不同半径处沿360度周向取了20个监测点得到的速度均值），这个在速度云图2中也可以看出来：

经过时间光滑后Umean的结果如图3和4所示，出口处速度波动会消失，只是远离喷嘴处结果和实验偏离的有些大：

King_RVM

@东岳这个rhoPimpleFoam的结果也太好了，是修改了湍流模型吗？真期待您的计算结果。
还想问一下李老师，rhoSimpleFoam里面您修改了outlet的边界条件，并且增加了Ux，Uy，和Uz的边界值和初始值，主要是outlet边界条件对结果有重要影响吗？

King_RVM

这个是计算结果文件：https://www.jianguoyun.com/p/DQHmutkQ1IPKCBjH46cD

李东岳

https://www.jianguoyun.com/p/Dbf3umcQ9s3ZBhjhk6cD 这个是瞬态的文件，你对比一下吧，我改了一些内容，记不太清哪些有影响了，欢迎你的理解与分享！

感谢你的plot

King_RVM

@东岳好的，我学习一下，也非常感谢东岳大哥的指导与交流。

King_RVM

@东岳经过对比验证，我们发现fvSolution中绝对残差tolerance对结果的影响非常大，默认值是1e-6，调整成1e-9后结果有很大改善，改变这个值可以提高计算精度。

李东岳

@King_RVM 感谢分享！！

李东岳

@king_rvm 在一个关于 Axisymmetric Subsonic Jet 的问题中说：

我们发现fvSolution中绝对残差tolerance对结果的影响非常大，默认值是1e-6，调整成1e-9后结果有很大改善，改变这个值可以提高计算精度。

找了半天可算找到了。没想到1e-6变成1e-9的差距这么大。2年过去了。这个是我给你们改的么，我自己都忘了..

King_RVM

@李东岳是的李老师，之后就是用您调的结果算的

King_RVM

@李东岳对，我们也没想到影响会这么大，有时候这类原因都想不到