CFD中文网

    CFD中文网

    • 登录
    • 搜索
    • 最新

    LES直流槽道边界层模拟,如何得到正则化速度u+以及正则化坐标y+?

    OpenFOAM
    槽道湍流 大涡模拟
    5
    45
    428
    正在加载更多帖子
    • 从旧到新
    • 从新到旧
    • 最多赞同
    回复
    • 在新帖中回复
    登录后回复
    此主题已被删除。只有拥有主题管理权限的用户可以查看。
    • C
      coolhhh @学流体的小明 最后由 编辑

      @学流体的小明 收到,非常感谢

      1 条回复 最后回复 回复 引用
      • 李东岳
        李东岳 管理员 @学流体的小明 最后由 李东岳 编辑

        Komen那个是的

        但是我看EugeneDeVilliers博士论文那个用的6万网格计算的u+y+还可以(图5.5)

        我觉得这个是很重要的内容,如果openfoam官方的算例,出现了网格分辨率不够导致的错误,这个bug必须应该处理

        CFD高性能服务器 http://dyfluid.com/servers.html

        C 1 条回复 最后回复 回复 引用
        • H
          hongjiewang @学流体的小明 最后由 编辑

          @学流体的小明 先前我做过这个工作,用于大涡模拟入口的前置算例,我当时画u+--y+使用的你说的方法二得到的,效果还行

          1 条回复 最后回复 回复 引用
          • H
            hongjiewang @学流体的小明 最后由 李东岳 编辑

            @学流体的小明

            对于此,我想请问,Um即平均速度如何给定,雷诺数的计算,假如我是在半宽通道中跑的,即下表面为wall,上表面为symmetry,此时,雷诺数公式里面应该把2h换成水力半径吗

            学流体的小明 1 条回复 最后回复 回复 引用
            • C
              coolhhh @李东岳 最后由 编辑

              @李东岳 在 LES直流槽道边界层模拟,如何得到正则化速度u+以及正则化坐标y+? 中说:

              Komen那个是的

              但是我看EugeneDeVilliers博士论文那个用的6万网格计算的u+y+还可以(图5.5)

              我觉得这个是很重要的内容,如果openfoam官方的算例,出现了网格分辨率不够导致的错误,这个bug必须应该处理

              李老师,EugeneDeVilliers博士论文图5.5,图名写normalised by the DNS shear velocity,是否有可能正则化用的$U_{\tau }$不一样?
              1.jpg

              学流体的小明 李东岳 2 条回复 最后回复 回复 引用
              • 学流体的小明
                学流体的小明 @hongjiewang 最后由 学流体的小明 编辑

                @hongjiewang
                Um是通过在动量方程中直接添加源项实现的,在v2012版本中是:

                /*--------------------------------*- C++ -*----------------------------------*\
                | =========                 |                                                 |
                | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
                |  \\    /   O peration     | Version:  v2012                                 |
                |   \\  /    A nd           | Website:  www.openfoam.com                      |
                |    \\/     M anipulation  |                                                 |
                \*---------------------------------------------------------------------------*/
                FoamFile
                {
                    version     2.0;
                    format      ascii;
                    class       dictionary;
                    location    "constant";
                    object      fvOptions;
                }
                // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
                
                momentumSource
                {
                    type            meanVelocityForce;
                
                    selectionMode   all;
                
                    fields          (U);
                    Ubar            (0.4 0 0);
                }
                
                
                // ************************************************************************* //
                

                程序会自动调整压力梯度以达到Ubar的预设值。
                之后我是打算用给定的动量方程源项,相当于给定$\partial p / \partial x$,通过vectorSemiImplicitSource实现。也是在fvOptions中调用。4800是随便给的,和上面0.4m/s的平均速度没有联系哈。$\partial p / \partial x$和$\tau_w$的关系见方法四。
                有一些参考:
                http://xiaopingqiu.github.io/2016/03/20/fvOptions2/
                https://caefn.com/openfoam/fvoptions-semiimplicitsource

                momentumSource
                {
                   type vectorSemiImplicitSource;
                   active on;
                   
                
                   vectorSemiImplicitSourceCoeffs
                   {
                      selectionMode all;
                      volumeMode        specific; //absolute
                      injectionRateSuSp
                      {
                         U           ( (4800 0 0) 0); //partial p / partial x
                      }
                   }
                }
                
                李东岳 1 条回复 最后回复 回复 引用
                • 学流体的小明
                  学流体的小明 @coolhhh 最后由 编辑

                  那EugeneDeVilliers博士论文的结果也是很好的,你看我方法一算出来的$u_{\tau}$,也相当于DNS的结果吧,画出来之后粘性底层完全对不上。
                  而且就算真的有DNS结果,我画出来的曲线的形状,也无法和理论解对上。就关键是形状它就对不上。

                  1 条回复 最后回复 回复 引用
                  • 学流体的小明
                    学流体的小明 最后由 编辑

                    @hongjiewang
                    槽道的雷诺数公式都是用两倍半高来,不管实际模拟了多少的。
                    摩擦雷诺数是用的一倍半高。

                    其中h为槽道半高,uτ为摩擦速度,ν为运动粘度。根据Pope在Turbulent Flows一书中提到的公式

                    这个公式里面的Re是用的两倍半高,Re_tau是一倍半高。

                    1 条回复 最后回复 回复 引用
                    • 李东岳
                      李东岳 管理员 @coolhhh 最后由 李东岳 编辑

                      @coolhhh

                      normalised by the DNS shear velocity

                      我觉得这句话提供了很多信息,他们应该是用方法1计算的$u_\tau$。

                      如果这样的话,应该是个小bug,就是网格分辨率不够。

                      很难相信会留这么一个小bug在里面

                      CFD高性能服务器 http://dyfluid.com/servers.html

                      1 条回复 最后回复 回复 引用
                      • Referenced by  李东岳 李东岳 
                      • 李东岳
                        李东岳 管理员 @学流体的小明 最后由 编辑

                        @学流体的小明 在 LES直流槽道边界层模拟,如何得到正则化速度u+以及正则化坐标y+? 中说:

                        通过vectorSemiImplicitSource实现

                        你可以试一下这个,我直觉感觉这个可能会导致流速持续的增加下不来了

                        CFD高性能服务器 http://dyfluid.com/servers.html

                        学流体的小明 1 条回复 最后回复 回复 引用
                        • 学流体的小明
                          学流体的小明 @李东岳 最后由 编辑

                          那不会的,我最开始就是用vectorSemiImplicitSource的,另外一个算例的正则化的结果和

                          方法二(方法三和四也是这样子的图)摩擦速度归一化结果是

                          类似。您说的流速持续增加很可能是压力梯度的具体数值没设置好,设置得大了。我之前设置的小了,发现平均速度、以及采样的某个点的速度总是在下降。就是驱动力不够,有耗散,速度就降下来了。


                          槽道的平均速度是会趋于稳定的。
                          $h=0.005m$,$U=4.0m/s$

                          bd4ea1d3-b703-4e0a-adce-823944d6b06f-image.png

                          1 条回复 最后回复 回复 引用
                          • 李东岳
                            李东岳 管理员 最后由 李东岳 编辑

                            在动量方程添加一个固定的源项,直觉理解会一直往上增加速度。你这么一说,那就是要设置的刚刚好,这个设置的压力梯度与本身的压力梯度相抗衡。好像是需要这样。不能大、不能小。有道理啊!:146:

                            CFD高性能服务器 http://dyfluid.com/servers.html

                            学流体的小明 1 条回复 最后回复 回复 引用
                            • C
                              coolhhh 最后由 编辑

                              @李东岳 我用OpenFOAM-8也算了channel395,用tutorial算例计算,采用 @学流体的小明 提供的MATLAB程序中method2算的$U_\tau$:

                              1. tutorial算例用的WALE模型,0文件夹其实不需设置k和nuTilda,只需要设置p、U和nut。
                              2. tutorial算例是给定了一个湍流初始场。一开始没注意删掉了,没有初始场到1000s计算不出湍流。初始湍流场的不同对结果是否影响会比较大?
                              3. 一般认为开始计算的一段时间的湍流还没充分发展,计算平均场是都会剔除。tutorial算例提供的统计是从0时刻就开始统计。不过我尝试从1000s才开始统计,共计算5000s,结果与默认算例结果一致。
                                同时还提取了EugeneDeVilliers博士论文图5.5中6万网格的结果,其实跟我们的计算结果一致的,都是有偏差!推测他也是用method2算的。他的论文没有给出参考线,导致看起来结果很好。 图中Spalding曲线是根据《无痛苦NS方程笔记》中提到的公式画出对比。
                                1.jpg

                              (1)tutorial算例默认设置结果:从0秒开始统计

                              u+ - y+.png
                              utao-t.png


                              (2)从1000秒开始统计。1000s的结果其实是1000s这个时刻的瞬时结果,但对整个底面平均后,与后面时刻的结果基本一致。
                              u+ - y+.png
                              utao-t.png
                              4. tutorial算例默认设置,method2计算结果都是在粘性子层结果较好,log区结果不吻合。尝试$U_\tau$按照method1取值,即$U_\tau=0.0079$,结果如下。说明了目前算例结果,改变$U_\tau$,只能满足粘性子层或者log区之一,是无法同时满足这两个区。想要同时满足,可能需要尝试加密网格之类方式
                              u+ - y+.png

                              C 2 条回复 最后回复 回复 引用
                              • C
                                coolhhh @coolhhh 最后由 编辑

                                @学流体的小明 @李东岳 不好意思,刚再仔细看了下,我的getdata数据好像有点问题,我再核实一下

                                1 条回复 最后回复 回复 引用
                                • C
                                  coolhhh @coolhhh 最后由 编辑

                                  @学流体的小明 @李东岳
                                  很抱歉撒上面的getdata数据不小心搞错了,更正如下

                                  1. 先说结论,EugeneDeVilliers博士论文图5.5中6万网格的结果,是按照method1计算的。
                                    按照method1,$U_{\tau}=0.0079$,结果如下图所示。EugeneDeVilliers结果也是在粘性子层偏差较大,和method1计算结果基本一致
                                    u+ - y+.png

                                  2. tutorial算例默认设置结果:从0秒开始统计,结果更正如下:
                                    u+ - y+.png

                                  3. tutorial算例从1000秒开始统计,结果更正如下
                                    u+ - y+.png

                                  学流体的小明 1 条回复 最后回复 回复 引用
                                  • 学流体的小明
                                    学流体的小明 @李东岳 最后由 学流体的小明 编辑

                                    槽道湍流是无法一直流动下去的,它里面总会有耗散,所以需要施加驱动。
                                    施加驱动的两种方式本质上都是给一个压力梯度,稳定之后,压力梯度就平衡了壁面摩擦力。
                                    其实是这么来的:

                                    1959afe8-18ce-45e1-8206-0f7c65a9f830-9c197790f8b2d16fe0b0eafc7d86828.jpg

                                    也不是不能大、不能小,设置得大了,槽道流速增加,会达到新的平衡的。

                                    不能大、不能小。

                                    L 1 条回复 最后回复 回复 引用
                                    • 学流体的小明
                                      学流体的小明 @coolhhh 最后由 编辑

                                      做的真好,谢谢你的探索和分享
                                      所以EugeneDeVilliers结果也不是那么好,只是他没给出参考线😂

                                      1 条回复 最后回复 回复 引用
                                      • 李东岳
                                        李东岳 管理员 最后由 李东岳 编辑

                                        嗯,经过你们2个的验证,应该就是自带tutorials的问题了。

                                        另外压力梯度那个,我最近在方程里面加了科氏力,也是一个规规矩矩的槽道,只不过尺度非常大(大气层)。目前算出来的结果就是速度持续的在旋转稳定不下来,与类似压力梯度的效果类似,我深入研究下咋回事,看看你那个图能不能解释一下

                                        CFD高性能服务器 http://dyfluid.com/servers.html

                                        1 条回复 最后回复 回复 引用
                                        • 李东岳
                                          李东岳 管理员 最后由 李东岳 编辑

                                          http://dyfluid.com/boundaryFoam.html 我顺着这个里面的方程2,在y方向做积分,继续推了一下:

                                          $$
                                          \int\left(-\frac{\partial \tau_{xy}}{\partial y}\right)\rd y=-\int\left(\frac{\partial p}{ \partial x}\right)\rd y
                                          $$

                                          $$
                                          \tau_{xy}^h -\tau_{w}=-\tau_{w}=\frac{\partial p}{ \partial x}h
                                          $$

                                          所以充分发展的管道,如果壁面剪切力知道的话,应该这样给一个压力梯度。

                                          也不是不能大、不能小,设置得大了,槽道流速增加,会达到新的平衡的。

                                          对,应该是。:146: :146: :146: 还得持续学习啊!

                                          CFD高性能服务器 http://dyfluid.com/servers.html

                                          1 条回复 最后回复 回复 引用
                                          • 李东岳
                                            李东岳 管理员 最后由 编辑

                                            我最近在方程里面加了科氏力,也是一个规规矩矩的槽道,只不过尺度非常大(大气层)。目前算出来的结果就是速度持续的在旋转稳定不下来

                                            我发现添加压力梯度跟平均速度,在我这里还不太一样。添加速度的话,如果是x方向速度,因为科室力会产生旋转,x方向速度一直不变,但是y会增加,导致速度增加。压力梯度却没有。因为跟本帖关系也不大,我先研究一下,然后新开一个帖子讨论。

                                            CFD高性能服务器 http://dyfluid.com/servers.html

                                            李东岳 1 条回复 最后回复 回复 引用
                                            • Referenced by  李东岳 李东岳 
                                            • 李东岳
                                              李东岳 管理员 @李东岳 最后由 编辑

                                              @李东岳 https://www.cfd-china.com/topic/6259/添加动量方程源项导致速度持续增加

                                              CFD高性能服务器 http://dyfluid.com/servers.html

                                              1 条回复 最后回复 回复 引用
                                              • Referenced by  学流体的小明 学流体的小明 
                                              • L
                                                luofq-sysu @学流体的小明 最后由 编辑

                                                @学流体的小明 感谢此贴,感谢几位的研究,解答了我的很多疑问。
                                                并请教:槽道流LES计算的入口inflow速度场影响大吗?
                                                因为我近期也开展了多个雷诺数Re180~5200的槽道流计算,但计算得到的壁面摩擦速度偏小(大约只有目标值的一半)。
                                                对比了一下你的channel1000案例文件,发现一个主要的差别是初始场0/U,看起来你是采用了一种生成脉动速度的入口条件吗?

                                                学流体的小明 1 条回复 最后回复 回复 引用
                                                • 学流体的小明
                                                  学流体的小明 @luofq-sysu 最后由 编辑

                                                  @luofq-sysu

                                                  并请教:槽道流LES计算的入口inflow速度场影响大吗?

                                                  我的计算域没有入口出口之分,沿流动方向的两个边界都是周期性的边界条件,展向也是周期性的边界条件。所以无法回答你的问题。

                                                  看起来你是采用了一种生成脉动速度的入口条件吗?

                                                  是的,用的就是EugeneDeVilliers在他博士论文中提到的初始化槽道流流场的方法。

                                                  L 1 条回复 最后回复 回复 引用
                                                  • L
                                                    luofq-sysu @学流体的小明 最后由 编辑

                                                    谢谢。@学流体的小明

                                                    我的计算域没有入口出口之分,沿流动方向的两个边界都是周期性的边界条件,展向也是周期性的边界条件。所以无法回答你的问题。

                                                    我的算例也是周期边界的槽道流,其实应该是internalField初始场的影响。
                                                    前面我采用的是uniform的初始速度场,计算湍流雷诺数远达不到目标值。最近按你算例的codeStream代码加入了初始脉动速度,计算雷诺数比较接近目标值了。
                                                    以Re_tau=1000为例,我现在计算得到的雷诺数数值是910

                                                    是的,用的就是EugeneDeVilliers在他博士论文中提到的初始化槽道流流场的方法。

                                                    EugeneDeVilliers大佬的论文提到了这种湍流需要加入初始扰动才能发展,这个初始化方法也有一个 perturbU的开源代码。不过codeStream的实现更方便:146:

                                                    1 条回复 最后回复 回复 引用
                                                    • 李东岳
                                                      李东岳 管理员 最后由 编辑

                                                      最近按你算例的codeStream代码加入了初始脉动速度

                                                      这个codeStream写的扰动,可以保证$\nabla\cdot\bfU=0$么

                                                      CFD高性能服务器 http://dyfluid.com/servers.html

                                                      L 1 条回复 最后回复 回复 引用
                                                      • L
                                                        luofq-sysu @李东岳 最后由 编辑

                                                        这个codeStream写的扰动,可以保证$\nabla\cdot\bfU=0$么

                                                        我只看到codeStream跟perturbU的代码是一致的,具体的公式在EugeneDeVilliers博士论文5.1.2章也能找到

                                                        1 条回复 最后回复 回复 引用
                                                        • First post
                                                          Last post