关于interFoam 中的alphaEqn里面的两个系数cAlpha和icAlpha
-
版本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)); }
查看了一下源文件, 这两个系数都是用户可以指定的,小弟不太明白,这两个系数的作用的什么?有知道的朋友告诉一下,谢谢!
-
@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; }
-
在论文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 α. -
@supersoldier thanks
-
@supersoldier 非常有用的回复,直接帮我解决了很多疑问。谢谢!!
-
@supersoldier 你好,可以再详细解释一下不?