Skip to content
  • 最新
  • 版块
  • 东岳流体
  • 随机看[请狂点我]
皮肤
  • Light
  • Cerulean
  • Cosmo
  • Flatly
  • Journal
  • Litera
  • Lumen
  • Lux
  • Materia
  • Minty
  • Morph
  • Pulse
  • Sandstone
  • Simplex
  • Sketchy
  • Spacelab
  • United
  • Yeti
  • Zephyr
  • Dark
  • Cyborg
  • Darkly
  • Quartz
  • Slate
  • Solar
  • Superhero
  • Vapor

  • 默认(不使用皮肤)
  • 不使用皮肤
折叠
CFD中文网

CFD中文网

  1. CFD中文网
  2. OpenFOAM
  3. Stokes流中的圆球绕流的阻力系数,怎么设置才能算准。

Stokes流中的圆球绕流的阻力系数,怎么设置才能算准。

已定时 已固定 已锁定 已移动 OpenFOAM
56 帖子 4 发布者 64.8k 浏览
  • 从旧到新
  • 从新到旧
  • 最多赞同
回复
  • 在新帖中回复
登录后回复
此主题已被删除。只有拥有主题管理权限的用户可以查看。
  • xiezhuoyuX 离线
    xiezhuoyuX 离线
    xiezhuoyu
    写于 最后由 李东岳 编辑
    #41

    @Prometheus10

    77c032de-a705-453e-899c-0e00994b6570-image.png

    请教,symmetry边界,直观上理解是没有界面法向的速度变化吧,没有速度变化,为什么有切应力?

    Prometheus10P 1 条回复 最后回复
  • Prometheus10P 离线
    Prometheus10P 离线
    Prometheus10
    在 中回复了 xiezhuoyu 最后由 李东岳 编辑
    #42

    @xiezhuoyu 您好 @xiezhuoyu ,先尝试回答您的提问:
    这个图带来了误解(这个图其实是partial-Slip的图,我放这个图主要是想请大家帮忙厘清下u·0和气泡“滑移界面”以及partial-Slip之间的关系)。

    那么,回到问题上来:
    u·η=0(也就是symmetry边界的数学表达,COMSOL里可以看到对称边界的数学方程就是u·η=0)的作用是什么? 我想象对于一个气泡界面处的流体而言,施加这个u·η=0的条件,并不是约束沿气泡界面法向的速度变化,而是要求流体沿壁面“贴体”进行运动,也就是说靠近壁面处的流体,不能沿法向运动,只能沿着切向运动。
    但是!切向有没有运动,需要来看切应力的大小,也就是看有没有对切应力做约束。展开来说:如果此时在u·η=0的条件上继续添加切应力为0的约束条件,那么就形成free-slip边界条件;反之,如果添加的是切应力无限大的约束条件,那么就形成noSlip边界条件(u=0)。

    1 条回复 最后回复
  • Prometheus10P 离线
    Prometheus10P 离线
    Prometheus10
    在 中回复了 李东岳 最后由 编辑
    #43

    @李东岳 感谢东岳老师,还有Pengdi老师,我有一些想明白的地方,同时也还有一些疑惑的地方,同时还有我最新的验证结果。我缕一缕,继续补充上来!:xinxin:

    1 条回复 最后回复
  • 李东岳李 在线
    李东岳李 在线
    李东岳 管理员
    写于 最后由 编辑
    #44

    妥了,实锤了,就是slip与symmetry的问题,估计slip没有做到完全的无剪切。我这几天详细看一下

    针对我上面说的,这个不对,slip可以做到无剪切。


    在测试算例,我尝试把球表面的值输出出来,第一个是slip,后面是symmetry:

    U & n = -4.2351647e-22
    max(U & n) = 4.3368087e-19
    average(U & n) = -2.494512e-20
    average(mag(U)) = 0.0045558926
    p face value = 0.0002584113
    p cell value = 0.0002584113
    average(p face) = average(p) [0 2 -2 0 0 0 0] -1.1678857e-07
    
    U & n = 4.2351647e-22
    max(U & n) = 4.3368087e-19
    average(U & n) = -8.4957405e-21
    average(mag(U)) = 0.0030605553
    p face value = 0.00026689782
    p cell value = 0.00026689782
    average(p face) = average(p) [0 2 -2 0 0 0 0] -5.620298e-09
    

    简单嫩看出:

    1. 边界处slip跟symmetry的$\bfU_f\cdot\bfn_f\approx 0$,所以二者都可以把边界处法向速度限制为0。
    2. 压力严格的$\nabla p\cdot\bfn_f=0$,但是可以看出压力的值差距还是不小。看average(p face)这个值
    3. 因为slip是衍生来自symmetry,因此切向方向二者处理是一样的,都是零梯度。

    所以我觉得slip与symmetry,二者都可以满足数学上的法向速度为0以及切向零梯度(在这里我不太清楚sci里面说的free slip是不是这个意思)。但是在求解的时候边界条件导致算法求解流程导致结果不一致,这体现在average(p face)以及average(mag(U))都有差异。

    在求解的时候,边界条件导致算法求解流程导致结果不一致

    这方面需要更详细看一下。

    http://dyfluid.com/index.html
    需要帮助debug算例的看这个 https://cfd-china.com/topic/8018

    P 1 条回复 最后回复
  • P 离线
    P 离线
    pengdi
    在 中回复了 李东岳 最后由 编辑
    #45

    @李东岳 补充一下我今天做的测试。我用雷诺数5000的二维管道流开展了测试,将左下角其中一个边设置为symmetry或者slip边界开展了对比。计算域如图:
    a9a4dd2b-26ad-47bb-8c1f-f5b7bfe5cff1-image.png
    计算得到沿着这条边的速度和压力分布基本一致,不过检测到了微弱的法向流速,计算结果如图:
    a77c2644-014c-40d3-90de-125071ee0b69-image.png
    8b272cd4-8aae-41d0-8987-12c8b5de2257-image.png 。
    东岳老师,我也感觉应该是算法的问题,边界条件本身应该没啥问题。
    备注:我测试的是平直边界,弧形边界可能区别比较大。

    P Prometheus10P 2 条回复 最后回复
  • P 离线
    P 离线
    pengdi
    在 中回复了 pengdi 最后由 编辑
    #46

    @pengdi 我再对上幅图进一步补充一下,上一幅图数据通过paraview中的Plot on Intersection Curve截取的,这种方法截取的是边界上一层网格的数据。在边界处,使用slip边界计算出来的是没有法向速度的,只是边界上一层网格有。使用symmetry边界,我还不知道如何检测边界面的速度。

    1 条回复 最后回复
  • 李东岳李 在线
    李东岳李 在线
    李东岳 管理员
    写于 最后由 编辑
    #47

    使用symmetry边界,我还不知道如何检测边界面的速度

    这个需要在求解器层面输出出来原始数据,用Info(麻烦一点)

    http://dyfluid.com/index.html
    需要帮助debug算例的看这个 https://cfd-china.com/topic/8018

    P 1 条回复 最后回复
  • P 离线
    P 离线
    pengdi
    在 中回复了 李东岳 最后由 编辑
    #48

    @李东岳 好的,东岳老师。明天我对比一下不同边界类型的圆柱绕流和方柱绕流计算,感觉曲面边界的影响可能会很大。

    1 条回复 最后回复
  • Prometheus10P 离线
    Prometheus10P 离线
    Prometheus10
    在 中回复了 pengdi 最后由 编辑
    #49

    @pengdi Pengdi老师,@李东岳 东岳老师
    我来补充一下低雷诺数条件下(Re=0.1时)对圆球型壁面设置为symmetry或者slip边界的切向速度的结果:
    依然是Wedge楔形网格的圆球绕流文件,然后将仿真结果在Paraview中用 plotOverLine提取:

    首先是symmetry边界条件,第一个网格体心的速度大概是主流速度的1/2;
    

    53ea6b1b-1607-4765-8649-5944793a4c49-image.png

    然后是slip边界条件,其第一个网格体心的速度就是主流速度;
    

    6c0e1e3a-d6b3-452c-a6c2-6fdc3084323c-image.png

    沿用一下这个示意图,symmetry边界也就是u·n=0确实很好描述了气泡界面处的速度分布,和Multiphase flow handbook第25页给出的公式也很吻合。(我这里给出的结果,准确来说,应该要说成是紧挨界面处,而不是界面上的切向速度,因为我只会用Paraview显示网格体心速度,还没能掌握如何用Info速度界面上的切向速度):zoule:
    Michaelides, E., Crowe, C.T., Schwarzkopf, J.D., 2016. Multiphase flow handbook. CRC Press.

    af9c4c27-dca5-4e83-8810-3935ac9df36a-image.png

    参考文献:
    Liu, B., Manica, R., Xu, Z., Liu, Q., 2020. The boundary condition at the air–liquid interface and its effect on film drainage between colliding bubbles. Curr. Opin. Colloid Interface Sci. 50, 101374.
    
    P 1 条回复 最后回复
  • P 离线
    P 离线
    pengdi
    在 中回复了 Prometheus10 最后由 编辑
    #50

    @Prometheus10 您好,关于圆柱绕流的slip和symmetry边界计算我也计算对比过,我计算的雷诺数等于200,差别确实很大。我怀疑可能正方形绕流差别就没那么大了,今天比较忙,还没验证,我看明天能不能验证一下,稍后补充上来。

    1 条回复 最后回复
  • Prometheus10P 离线
    Prometheus10P 离线
    Prometheus10
    写于 最后由 编辑
    #51

    @李东岳 东岳老师好,我现在想根据运行时计算得到的Cd值的大小反过来适配速度进口的速度大小,直至阻力和浮力平衡,通过看您其他帖子,好像需要在0/U文件夹下使用codedFixedValue边界来完成,想进一步咨询您:

    在codedFixedValue如何调用当前时刻的Cd?
    
        inlet
        {
            type            codedFixedValue;
            value           $internalField;
            name            dragForceBalanceInlet; //name of new BC type
            code
            #{
                const fvPatch& inletPatch = this->patch();
                
                vectorField& vf = *this; 
                forAll(vf, i)
                { 
                    scalar z = inletPatch.Cf()[i].y();
                    scalar rhoL = 998.13;
                    scalar rhoG = 1.225;
                    scalar myG = 9.81;
                    scalar myD = 0.0000505;
                    vf[i].x() = sqrt(4*(rhoL-rhoG)*myG*myD/(3*Cd*rhoL));//想在这里调用Cd
                    vf[i].y() = 0.0; 
                    vf[i].z() = 0.0; 
                }
            #};
           
            codeOptions
            #{
                  -I$(LIB_SRC)/finiteVolume/lnInclude \
                  -I$(LIB_SRC)/meshTools/lnInclude
             #};
      
            codeInclude
             #{
                    #include "fvCFD.H"
                    #include <cmath>
                   #include <iostream>
            #};
        }
    

    上述代码是通过浮力和阻力的力平衡,得到阻力系数和速度的等式。
    ea0fbbed-4b03-4419-9554-8fb343f250ac-image.png
    Tomiyama, A., Kataoka, I., Zun, I., Sakaguchi, T., 1998. Drag Coefficients of Single Bubbles under Normal and Micro Gravity Conditions. JSME International Journal Series B 41, 472–479.

    1 条回复 最后回复
  • 李东岳李 在线
    李东岳李 在线
    李东岳 管理员
    写于 最后由 编辑
    #52

    Cd这个在forceCoeff.C里面计算出来的,这个我看没有公开的接口,也不知道哪个是先声明的。另外一个方法是是你自己计算一下Cd。

    $$
    Cd=\frac{\left( \sum p_f\bfS_f + \sum \bfS_f\cdot \tau_f \right)\cdot\mathbf{direction}_{drag}}{(\sum|\bfS_f|) \frac{1}{2}\rho|\bfU|^2}
    $$

    里面的$\tau$定义在eddyViscosity里面的sigma(),被称之为R。不过我不确定能否在codedFixedValue里面去引用(不确定这个前后关系是哪个先声明的)

    http://dyfluid.com/index.html
    需要帮助debug算例的看这个 https://cfd-china.com/topic/8018

    Prometheus10P 1 条回复 最后回复
  • Prometheus10P 离线
    Prometheus10P 离线
    Prometheus10
    在 中回复了 李东岳 最后由 编辑
    #53

    @李东岳 好的,谢谢东岳老师。:xinxin:

    1 条回复 最后回复
  • 李东岳李 李东岳 被引用 于这个主题
  • Prometheus10P 离线
    Prometheus10P 离线
    Prometheus10
    写于 最后由 Prometheus10 编辑
    #54

    @pengdi @李东岳
    东岳老师, PengDi老师,您们好

    我目前仍然有一点没有想明白,想再请两位老师指点:

    首先是Slip边界条件的定义:
    【在界面处,法向速度为0,同时切向速度的梯度为0.】

    沿用知乎博主的解析神秘色彩OpenFOAM 边界条件系列解析—Slip边界
    99e75478-1ff3-4522-87cd-5ac546b52cf6-image.png
    COMSOL的帮助文档也是一致的内容
    Theory for the Wall Boundary Condition

    再有是东岳老师在#44楼所提到的:

    因为slip是衍生来自symmetry,因此切向方向二者处理是一样的,都是零梯度。
    

    那么,slip和symmetry按道理都应该同时符合这几个公式
    $\mathbf{U}_{b\perp}=(\mathbf{U}_b\cdot\mathbf{n})\mathbf{n}=0$
    $\mathbf{U}_b{|}=\mathbf{U}_c{|}$
    $\mathbf{U}_c{|}=\mathbf{U}_c-\left(\mathbf{U}_c\cdot\mathbf{n}\right)\mathbf{n}$

    但是,当Slip和Symmetry应用到弯曲边界的时候:为什么速度有这么大的差异,并由此得出的阻力系数也大不相同。
    雷诺数Re=1时,阻力系数:
    CD(slip)=13.49,
    CD(symmetry)=17.6

    c76cb706-941e-450b-b24c-0fefb60c1ed6-image.png
    9386b544-c256-430a-8381-fffa6379bf91-image.png

    从以上两个图结果上来看:
    Symmetry的速度明显是和主流速度之间渐变过渡的(左图蓝线,大概从U/2逐渐发展到U),也就意味着流体在界面处是受到切应力的;
    而Slip的速度是和主流速度一致(左图蓝线,Slip界面处就是和主流速度一样的U),因此我认为:

    symmetry条件应用到弯曲边界的时候,并不会对弯曲边界的切应力做限制,因此弯曲边界的symmetry条件只是以下表达式:
    $\mathbf{U}_{b\perp}=(\mathbf{U}_b\cdot\mathbf{n})\mathbf{n}=0$

    此时Symmetry已经不再等同于Slip条件,其真实意义已经为:

    Symmetry:【在界面处,法向速度为0,同时切向速度的梯度可以不为0.】
    而Slip依然是:【在界面处,法向速度为0,同时切向速度的梯度为0.】

    不知道这样理解是否正确?

    P 1 条回复 最后回复
  • Prometheus10P Prometheus10 被引用 于这个主题
  • 李东岳李 在线
    李东岳李 在线
    李东岳 管理员
    写于 最后由 编辑
    #55

    从你的结论来看,symmetry确实存在切应力,这与我之前的理解是相反的。但是基金会那本书里面也写了,symmetry是所有的场都是symmetry,是一种更严厉的限制。但是目前来看像是相反的。目前我在备课还不能细看。你这个问题(symmetry与slip边界),还有论坛里面有个时间步长越小结果越差的,成了老大难了。。感觉是非常底层的bug。

    http://dyfluid.com/index.html
    需要帮助debug算例的看这个 https://cfd-china.com/topic/8018

    1 条回复 最后回复
  • P 离线
    P 离线
    pengdi
    在 中回复了 Prometheus10 最后由 编辑
    #56

    @Prometheus10 对 如果这样理解的话 为什么symmetry计算结果和理论解更相近呢

    1 条回复 最后回复

  • 登录

  • 登录或注册以进行搜索。
  • 第一个帖子
    最后一个帖子
0
  • 最新
  • 版块
  • 东岳流体
  • 随机看[请狂点我]