CFD中文网

    CFD中文网

    • 登录
    • 搜索
    • 最新

    关于interFoam 中的alphaEqn里面的两个系数cAlpha和icAlpha

    OpenFOAM
    5
    8
    5453
    正在加载更多帖子
    • 从旧到新
    • 从新到旧
    • 最多赞同
    回复
    • 在新帖中回复
    登录后回复
    此主题已被删除。只有拥有主题管理权限的用户可以查看。
    • Y
      yuan_neu 最后由 CFD中文网 编辑

      版本3.0.x中的interFoam里面的alphaEqn.H line46-54

      // Standard face-flux compression coefficient
      surfaceScalarField phic(mixture.cAlpha()*mag(phi/mesh.magSf()));
      
      // Add the optional isotropic compression contribution
      if (icAlpha > 0)
      {
          phic *= (1.0 - icAlpha);
          phic += (mixture.cAlpha()*icAlpha)*fvc::interpolate(mag(U));
      }
      

      查看了一下源文件, 这两个系数都是用户可以指定的,小弟不太明白,这两个系数的作用的什么?有知道的朋友告诉一下,谢谢!

      1 条回复 最后回复 回复 引用
      • C
        CFD中文网 最后由 编辑

        麻烦把所有alphaeqn.H里面的内容发一下:confused:

        CFD中国标准用户测试帐号
        目前由徐笑笑登录

        1 条回复 最后回复 回复 引用
        • Y
          yuan_neu 最后由 编辑

          @cfd-china alphaeqn.H里面的内容如下

          {
              word alphaScheme("div(phi,alpha)");
              word alpharScheme("div(phirb,alpha)");
          
              tmp<fv::ddtScheme<scalar> > ddtAlpha
              (
                  fv::ddtScheme<scalar>::New
                  (
                      mesh,
                      mesh.ddtScheme("ddt(alpha)")
                  )
              );
          
              // Set the off-centering coefficient according to ddt scheme
              scalar ocCoeff = 0;
              if
              (
                  isType<fv::EulerDdtScheme<scalar> >(ddtAlpha())
               || isType<fv::localEulerDdtScheme<scalar> >(ddtAlpha())
              )
              {
                  ocCoeff = 0;
              }
              else if (isType<fv::CrankNicolsonDdtScheme<scalar> >(ddtAlpha()))
              {
                  if (nAlphaSubCycles > 1)
                  {
                      FatalErrorIn(args.executable())
                          << "Sub-cycling is not supported "
                             "with the CrankNicolson ddt scheme"
                          << exit(FatalError);
                  }
          
                  ocCoeff =
                      refCast<fv::CrankNicolsonDdtScheme<scalar> >(ddtAlpha()).ocCoeff();
              }
              else
              {
                  FatalErrorIn(args.executable())
                      << "Only Euler and CrankNicolson ddt schemes are supported"
                      << exit(FatalError);
              }
          
              scalar cnCoeff = 1.0/(1.0 + ocCoeff);
          
              // Standard face-flux compression coefficient
              surfaceScalarField phic(mixture.cAlpha()*mag(phi/mesh.magSf()));
          
              // Add the optional isotropic compression contribution
              if (icAlpha > 0)
              {
                  phic *= (1.0 - icAlpha);
                  phic += (mixture.cAlpha()*icAlpha)*fvc::interpolate(mag(U));
              }
          
              // Do not compress interface at non-coupled boundary faces
              // (inlets, outlets etc.)
              forAll(phic.boundaryField(), patchi)
              {
                  fvsPatchScalarField& phicp = phic.boundaryField()[patchi];
          
                  if (!phicp.coupled())
                  {
                      phicp == 0;
                  }
              }
          
              tmp<surfaceScalarField> phiCN(phi);
          
              // Calculate the Crank-Nicolson off-centred volumetric flux
              if (ocCoeff > 0)
              {
                  phiCN = cnCoeff*phi + (1.0 - cnCoeff)*phi.oldTime();
              }
          
              if (MULESCorr)
              {
                  fvScalarMatrix alpha1Eqn
                  (
                      (
                          LTS
                        ? fv::localEulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
                        : fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
                      )
                    + fv::gaussConvectionScheme<scalar>
                      (
                          mesh,
                          phiCN,
                          upwind<scalar>(mesh, phiCN)
                      ).fvmDiv(phiCN, alpha1)
                  );
          
                  alpha1Eqn.solve();
          
                  Info<< "Phase-1 volume fraction = "
                      << alpha1.weightedAverage(mesh.Vsc()).value()
                      << "  Min(" << alpha1.name() << ") = " << min(alpha1).value()
                      << "  Max(" << alpha1.name() << ") = " << max(alpha1).value()
                      << endl;
          
                  tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
                  alphaPhi = talphaPhiUD();//based on updated alpha from alpha1Eqn.solve() to correct alphaPhi
          
                  if (alphaApplyPrevCorr && talphaPhiCorr0.valid())
                  {
                      Info<< "Applying the previous iteration compression flux" << endl;
                      MULES::correct(alpha1, alphaPhi, talphaPhiCorr0(), 1, 0);
          
                      alphaPhi += talphaPhiCorr0();
                  }
          
                  // Cache the upwind-flux
                  talphaPhiCorr0 = talphaPhiUD;
          
                  alpha2 = 1.0 - alpha1;
          
                  mixture.correct();
              }
          
          
              for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
              {
                  surfaceScalarField phir(phic*mixture.nHatf());
          
                  tmp<surfaceScalarField> talphaPhiUn
                  (
                      fvc::flux
                      (
                          phi,
                          alpha1,
                          alphaScheme
                      )
                    + fvc::flux
                      (
                         -fvc::flux(-phir, alpha2, alpharScheme),
                          alpha1,
                          alpharScheme
                      )
                  );
          
                  // Calculate the Crank-Nicolson off-centred alpha flux
                  if (ocCoeff > 0)
                  {
                      talphaPhiUn =
                          cnCoeff*talphaPhiUn + (1.0 - cnCoeff)*alphaPhi.oldTime();
                  }
          
                  if (MULESCorr)
                  {
                      tmp<surfaceScalarField> talphaPhiCorr(talphaPhiUn() - alphaPhi);
                      volScalarField alpha10("alpha10", alpha1);
          
                      MULES::correct(alpha1, talphaPhiUn(), talphaPhiCorr(), 1, 0);
          
                      // Under-relax the correction for all but the 1st corrector
                      if (aCorr == 0)
                      {
                          alphaPhi += talphaPhiCorr();
                      }
                      else
                      {
                          alpha1 = 0.5*alpha1 + 0.5*alpha10;
                          alphaPhi += 0.5*talphaPhiCorr();
                      }
                  }
                  else
                  {
                      alphaPhi = talphaPhiUn;
          
                      MULES::explicitSolve(alpha1, phiCN, alphaPhi, 1, 0);
                  }
          
                  alpha2 = 1.0 - alpha1;
          
                  mixture.correct();
              }
          
              if (alphaApplyPrevCorr && MULESCorr)
              {
                  talphaPhiCorr0 = alphaPhi - talphaPhiCorr0;
              }
          
              if
              (
                  word(mesh.ddtScheme("ddt(rho,U)"))
               == fv::EulerDdtScheme<vector>::typeName
              )
              {
                  rhoPhi = alphaPhi*(rho1 - rho2) + phiCN*rho2;
              }
              else
              {
                  if (ocCoeff > 0)
                  {
                      // Calculate the end-of-time-step alpha flux
                      alphaPhi = (alphaPhi - (1.0 - cnCoeff)*alphaPhi.oldTime())/cnCoeff;
                  }
          
                  // Calculate the end-of-time-step mass flux
                  rhoPhi = alphaPhi*(rho1 - rho2) + phi*rho2;
              }
          
              Info<< "Phase-1 volume fraction = "
                  << alpha1.weightedAverage(mesh.Vsc()).value()
                  << "  Min(" << alpha1.name() << ") = " << min(alpha1).value()
                  << "  Max(" << alpha1.name() << ") = " << max(alpha1).value()
                  << endl;
          }
          
          1 条回复 最后回复 回复 引用
          • supersoldier
            supersoldier 最后由 李东岳 编辑

            在论文DYNAMIC MESHING AROUND FLUID-FLUID INTERFACES WITH APPLICATIONS TO DROPLET TRACKING IN CONTRACTION GEOMETRIES第58页最上部有写cAlpha的作用。

            The cAlpha keyword specified in the file <case>/system/fvSolution is
            a factor that controls the compression of the interface, where 0 corresponds to no
            compression. After solving the volume fluid equation, the density and viscosity are
            modified using the new values of α.

            Y yhdthu 3 条回复 最后回复 回复 引用
            • Y
              yuan_neu @supersoldier 最后由 编辑

              @supersoldier thanks

              1 条回复 最后回复 回复 引用
              • Y
                yuan_neu @supersoldier 最后由 编辑

                @supersoldier 非常有用的回复,直接帮我解决了很多疑问。谢谢!!

                C 1 条回复 最后回复 回复 引用
                • yhdthu
                  yhdthu 讲师 @supersoldier 最后由 编辑

                  @supersoldier 你好,可以再详细解释一下不?

                  长风破浪会有时,直挂云帆济沧海

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

                    @yuan_neu 请问icAlpha的作用是什么?

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