Navigation

    CFD中文网

    CFD中文网

    • Login
    • Search
    • 最新

    关于correctPhi.H这个函数

    OpenFOAM
    6
    16
    10896
    Loading More Posts
    • Oldest to Newest
    • Newest to Oldest
    • Most Votes
    Reply
    • Reply as topic
    Log in to reply
    This topic has been deleted. Only users with topic management privileges can see it.
    • 金石为开
      金石为开 last edited by

      对于多相流VOF方法要先运行一下correctPhi.H这个文件,我见功能是修改通量,Required during start-up, restart, mesh-motion etc.在开始的地方运行,这个怎么理解?采用什么公式方法,请教大家~~~

      1 Reply Last reply Reply Quote
      • 李东岳
        李东岳 管理员 last edited by

        correctPhi.H:

        Flux correction functions to ensure continuity.
        
        Required during start-up, restart, mesh-motion etc. when non-conservative fluxes may adversely affect the prediction-part of the solution algorithm (the part before the first pressure solution which would ensure continuity). This is particularly important for VoF and other multi-phase solver in which non-conservative fluxes cause unboundedness of the phase-fraction.
        

        interFoam有一个动网格求解器叫interDyMFoam,两个求解器可以互相调用。从描述来看:correctPhi强制连续性方程。因此通过下述代码:

        while (pimple.correctNonOrthogonal())
            {
                // Solve for pcorr such that the divergence of the corrected flux
                // matches the divU provided (from previous iteration, time-step...)
                fvScalarMatrix pcorrEqn
                (
                    fvm::laplacian(rAUf, pcorr) == fvc::div(phi) - divU
                );
        
                pcorrEqn.setReference(0, 0);
                pcorrEqn.solve();
        
                if (pimple.finalNonOrthogonalIter())
                {
                    phi -= pcorrEqn.flux();
                }
            }
        

        求解由连续性方程构造的泊松方程后,重组守恒的phi场:phi -= pcorrEqn.flux();。

        这是一个数值上的考虑,可有可无。

        线上CFD课程 7月1日报名截止 http://dyfluid.com/class.html
        CFD高性能服务器 http://dyfluid.com/servers.html

        金石为开 zhanghan 2 Replies Last reply Reply Quote
        • 金石为开
          金石为开 @李东岳 last edited by

          @李东岳 东岳老师,我看过你写的pimpleDyFoam里面的讲解,对于动网格,变动之后有体积守恒定律,网格动了之后要进行通量修正。代码中的divU其实就是0场(公式里面对1做时间导数为0),不知道理解的对不对?
          还有我认为这是求解之前强制的做了一个压力求解和通量修正,我不清楚这是出于什么的数值上的考虑或者是物理上的考虑?发散?稳定?精度还是其他谢谢老师。另外您的fluent理论指南好了没有呀?很期待看看有什么学习参考的,谢谢老师

          1 Reply Last reply Reply Quote
          • 李东岳
            李东岳 管理员 last edited by 李东岳

            @金石为开 fvm::laplacian(rAUf, pcorr) == fvc::div(phi) - divU divU在没初始化的时候就是0。即GeometricZeroField

            还有我认为这是求解之前强制的做了一个压力求解和通量修正,我不清楚这是出于什么的数值上的考虑

            我觉得是数值上的。因为correctPhi有的时候可以关闭。并不是必须要开。不过从教材上来看,物理上也必须符合空间守恒。所以。我也不是很明确。

            另外,pimpleDyMFoam有个方程写错了,已经修改,参考Peric书方程12.13。不过代码还对不上,从方程12.13来看代码应该为:

            fvm::laplacian(rAU, pcorr) == fvc::div(mesh.phi()) - divU
            

            感兴趣你可以找个简单测试一下?要把correctPhi设置为on。

            线上CFD课程 7月1日报名截止 http://dyfluid.com/class.html
            CFD高性能服务器 http://dyfluid.com/servers.html

            金石为开 1 Reply Last reply Reply Quote
            • 李东岳
              李东岳 管理员 last edited by

              哦对了还有fluent理论指南。最近琐事太多了:sad: 。都已经翻译完了。就是格式太难统一了。因为是各自翻译的。用的是word。处理格式比较费时。

              线上CFD课程 7月1日报名截止 http://dyfluid.com/class.html
              CFD高性能服务器 http://dyfluid.com/servers.html

              1 Reply Last reply Reply Quote
              • zhanghan
                zhanghan @李东岳 last edited by

                @李东岳 李老师,您有研究过把ansys的计算结果导到openfoam中,或者是其他软件中的数据导到openfoam中吗??望您给一些指点!

                1 Reply Last reply Reply Quote
                • 金石为开
                  金石为开 @李东岳 last edited by

                  @李东岳 coreectPhi.C出错了??phi变成mesh.phi()??如果这里错了,那,,所有用到correctPhi的地方都错了呀?。。是不是呢?。。。

                  H 1 Reply Last reply Reply Quote
                  • H
                    Haining LUO @金石为开 last edited by

                    @金石为开 我觉得这里是因为网格动到了新的位置,而当前的phi相对于当前的网格就不再divergence free,网格更新了而phi场也需要更新,在当前网格下对absolute phi进行修正,而不是对应的方程12.13。假设了一个pcorrEqn场的拉普拉斯等于当前phi的不为零散度,然后用phi -= pcorrEqn.flux();修正。

                    1 Reply Last reply Reply Quote
                    • 程
                      程迪 last edited by 程迪

                      理论上是为了满足压力方程的相容性条件,实际上调用了adjustPhi.H。

                      参考

                      • CS205b/CME306 Lecture 16
                      • Jasak的回复

                      adjustPhi和pRef针对的是压力方程为全Neuman边界的情况,这会有两个问题:

                      • 方程系数矩阵M为奇异,解有无穷多,相差一个常数,所以需要pRef
                      • 方程相容性问题,需要adjustPhi

                      pRef

                      当方程为:
                      $$
                      \Delta p = f \text{ on }\Omega
                      $$
                      $$
                      \nabla p \cdot \mathbf n = g \text{ on }\partial \Omega
                      $$
                      首先,如果$p_1$是一个解,那么$p_2=p_1+const$也是方程的解,方程解不唯一,所以离散方程会有奇异性,无法求解。
                      解决方法是固定一个点的压力为给定值,就是设置pRef。

                      adjustPhi

                      p变量前面是直接加了微分的,所以需要对于全Neumann条件的情形还需要满足相容性条件:

                      $$
                      \int_{\Omega}{f dV} =\int_{\Omega}{\Delta p dV} = \int_{ \partial\Omega}{\nabla p \cdot\mathbf n dS}= \int_{ \partial\Omega}{g dS}
                      $$

                      这就是AdjustPhi去强制满足这个条件。
                      带入f和g的表达式

                      $$
                      f= \nabla\cdot\mathbf u^*
                      $$
                      $$
                      g=0
                      $$

                      可得:
                      $$
                      \int_{\Omega}{\nabla \cdot \mathbf u^* dV} =\int_{\partial \Omega}{\mathbf u^*dS}= 0
                      $$

                      你再看看adjustPhi出现的位置,是在解压力方程之前,对phiHbyA进行的调整。其实也就是这里的$\mathbf u^*$项。必须使其满足以上积分的相容性条件。

                      github: chengdi123000
                      网站:chengdi123000.github.io
                      本作品采用知识共享署名-非商业性使用-相同方式共享 4.0 国际许可协议进行许可。

                      H 我 2 Replies Last reply Reply Quote
                      • H
                        Haining LUO @程迪 last edited by

                        @程迪 多谢分享,不过上面的帖子说的不是correctPhi.H吗?
                        你这里提到的是adjustPhi(phiHbyA, U, p),我看icoFoam里面是在求解Poisson方程前,在pimpleoam和pimpleDyMFoam都是在pEqn.H里面,你提到的pRef也是,所以是为了求解Poisson方程作的准备。我觉得我们说的不是一个东西,这里是ajustPhi的帖子adjustPhi的作用是检查边界条件?

                        程 1 Reply Last reply Reply Quote
                        • 程
                          程迪 @Haining LUO last edited by

                          @Haining-LUO
                          correctPhi调用的就是adjustPhi。不过adjustPhi通常针对p,correctPhi是针对的pcorr。道理都一样的。保证微分方程相容性。

                          github: chengdi123000
                          网站:chengdi123000.github.io
                          本作品采用知识共享署名-非商业性使用-相同方式共享 4.0 国际许可协议进行许可。

                          H 1 Reply Last reply Reply Quote
                          • H
                            Haining LUO @程迪 last edited by

                            @程迪 明白你的意思了,在pimpleFoam和icoFoam里面都是在解p方程前调用了adjustPhi以保证相容性(通过修正流量保证连续方程的方式:其中涉及了两种边界条件fixedValue和inletOutletFvPatchField,修改的是inletOutletFvPatchField上的phi)。

                            但仔细看pimpleDyMFoam里面的correctPhi.H

                            if (mesh.changing()) // 如果是动网格,肯定会执行的部分
                            {
                                forAll(U.boundaryField(), patchI)
                                {
                                    if (U.boundaryField()[patchI].fixesValue())
                                    {
                                        U.boundaryField()[patchI].initEvaluate();
                                    }
                                }
                            
                                forAll(U.boundaryField(), patchI)
                                {
                                    if (U.boundaryField()[patchI].fixesValue())
                                    {
                                        U.boundaryField()[patchI].evaluate();
                            
                                        phi.boundaryField()[patchI] =
                                            U.boundaryField()[patchI]
                                          & mesh.Sf().boundaryField()[patchI];
                                    }
                                }
                            }
                            
                            {  // 这里我想说妈蛋……这个大括号是什么鬼
                                volScalarField pcorr
                                (
                                    IOobject
                                    (
                                        "pcorr",
                                        runTime.timeName(),
                                        mesh,
                                        IOobject::NO_READ,
                                        IOobject::NO_WRITE
                                    ),
                                    mesh,
                                    dimensionedScalar("pcorr", p.dimensions(), 0.0),
                                    pcorrTypes
                                );
                            
                                dimensionedScalar rAUf("rAUf", dimTime, 1.0);  // 构建了一个单位且均匀的有时间量纲的场
                            
                                while (pimple.correctNonOrthogonal()) // 这里其实还没有进入真正的pimpleFoam循环(while (pimle.loop())),但corrNonOrtho_已经存在。但是这里没有adjustPhi !! 为什么 ?? 这里虽然求解的不是物理场p,但是pcorr也得有个边界条件吧? (于是我在createPcorrTypes.H里面找到pcorrTypes几乎是zeroGradient,如果p设定的条件中有fixedValue,对应的pcorr也是fixedValue,所以如果没有fixedValue的p边界这里也是需要adjustPhi的!)
                                {
                                    fvScalarMatrix pcorrEqn
                                    (
                                        fvm::laplacian(rAUf, pcorr) == fvc::div(phi)  
                                    );
                            
                                    pcorrEqn.setReference(pRefCell, pRefValue);
                                    pcorrEqn.solve();
                            
                                    if (pimple.finalNonOrthogonalIter())
                                    {
                                        phi -= pcorrEqn.flux();
                                    }
                                }
                            
                                #include "continuityErrs.H"
                            }
                            

                            在这一段代码之后才是pimple.loop(),所以我觉得这是因为

                                    mesh.update();
                            
                                    // Calculate absolute flux from the mapped surface velocity
                                    phi = mesh.Sf() & Uf;
                            

                            用新网格算了新phi,而新的phi不是无散度的,因此用以上来修正之后再进入N-S方程求解。不过@李东岳 指出的代码和我这里OpenFOAM/2.3.1-foss-2016a不大一样,我觉得这里求解的并不是Peric书方程12.13,你觉得呢?

                            程 1 Reply Last reply Reply Quote
                            • 程
                              程迪 @Haining LUO last edited by

                              @Haining-LUO

                              Peric 12.13是空间守恒律吧。网格的东西,和这玩意儿有啥关系?

                              据我所知,volScalarField存的是体心的物理量,而不是体心物理量×体积,网格变形相当于物理量分布变化了,新phi和老phi值是一样,分布不同而已。

                              github: chengdi123000
                              网站:chengdi123000.github.io
                              本作品采用知识共享署名-非商业性使用-相同方式共享 4.0 国际许可协议进行许可。

                              H 1 Reply Last reply Reply Quote
                              • H
                                Haining LUO @程迪 last edited by Haining LUO

                                @程迪

                                phi是surfaceScalarField,在phi = mesh.Sf() & Uf;之前,phi还是定义在老网格的位置,通过Uf=fvc::interpolate(U)得到Uf后(我觉得是相对新网格的插值),这样就定义了一个新phi:phi = mesh.Sf() & Uf;,于是新的phi定义在新的网格面上,不再满足连续方程(旧phi是满足的),所以在进入pimple循环前通过pcorr来修正。

                                1 Reply Last reply Reply Quote
                                • 李东岳
                                  李东岳 管理员 last edited by

                                  { // 这里我想说妈蛋……这个大括号是什么鬼

                                  这个是OpenFOAM的代码风格。

                                  其他内容,我需要周末研究下才能跟你们讨论 :cheeky:

                                  线上CFD课程 7月1日报名截止 http://dyfluid.com/class.html
                                  CFD高性能服务器 http://dyfluid.com/servers.html

                                  1 Reply Last reply Reply Quote
                                  • 我
                                    我是河滩 @程迪 last edited by

                                    @程迪 压力都是零梯度边界条件,压力泊松方程的系数矩阵是奇异的,adjustphi让矩阵不奇异了吗?奇异矩阵怎么解的?

                                    动边界

                                    1 Reply Last reply Reply Quote
                                    • First post
                                      Last post

                                    CFD中文网 | 东岳流体 | 京ICP备15017992号-2
                                    论坛登录问题反馈可联系 li.dy@dyfluid.com