alphaEqn.H
Go to the documentation of this file.
1 {
2  word alphaScheme("div(phi,alpha)");
3  word alpharScheme("div(phirb,alpha)");
4 
5  surfaceScalarField phir("phir", phic*mixture.nHatf());
6 
7  Pair<tmp<volScalarField>> vDotAlphal = mixture.vDotAlphal();
8  const volScalarField& vDotcAlphal = vDotAlphal[0]();
9  const volScalarField& vDotvAlphal = vDotAlphal[1]();
10  const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal);
11 
12  tmp<surfaceScalarField> talphaPhi;
13 
14  if (MULESCorr)
15  {
16  fvScalarMatrix alpha1Eqn
17  (
18  fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
19  + fv::gaussConvectionScheme<scalar>
20  (
21  mesh,
22  phi,
23  upwind<scalar>(mesh, phi)
24  ).fvmDiv(phi, alpha1)
25  - fvm::Sp(divU, alpha1)
26  ==
28  + vDotcAlphal
29  );
30 
31  alpha1Eqn.solve();
32 
33  Info<< "Phase-1 volume fraction = "
34  << alpha1.weightedAverage(mesh.Vsc()).value()
35  << " Min(" << alpha1.name() << ") = " << min(alpha1).value()
36  << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
37  << endl;
38 
39  talphaPhi = alpha1Eqn.flux();
40  }
41 
42  volScalarField alpha10("alpha10", alpha1);
43 
44  for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
45  {
46  tmp<surfaceScalarField> talphaPhiCorr
47  (
48  fvc::flux
49  (
50  phi,
51  alpha1,
52  alphaScheme
53  )
54  + fvc::flux
55  (
57  alpha1,
59  )
60  );
61 
62  if (MULESCorr)
63  {
64  talphaPhiCorr.ref() -= talphaPhi();
65 
66  volScalarField alpha100("alpha100", alpha10);
67  alpha10 = alpha1;
68 
70  (
71  geometricOneField(),
72  alpha1,
73  talphaPhi(),
74  talphaPhiCorr.ref(),
76  (
77  divU*(alpha10 - alpha100)
78  - vDotvmcAlphal*alpha10
79  )(),
80  oneField(),
81  zeroField()
82  );
83 
84  // Under-relax the correction for all but the 1st corrector
85  if (aCorr == 0)
86  {
87  talphaPhi.ref() += talphaPhiCorr();
88  }
89  else
90  {
91  alpha1 = 0.5*alpha1 + 0.5*alpha10;
92  talphaPhi.ref() += 0.5*talphaPhiCorr();
93  }
94  }
95  else
96  {
98  (
99  geometricOneField(),
100  alpha1,
101  phi,
102  talphaPhiCorr.ref(),
104  (divU*alpha1 + vDotcAlphal)(),
105  oneField(),
106  zeroField()
107  );
108 
109  talphaPhi = talphaPhiCorr;
110  }
111 
112  alpha2 = 1.0 - alpha1;
113  }
114 
116 
117  Info<< "Liquid phase volume fraction = "
118  << alpha1.weightedAverage(mesh.V()).value()
119  << " Min(" << alpha1.name() << ") = " << min(alpha1).value()
120  << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
121  << endl;
122 }
tmp< surfaceScalarField > talphaPhi
Definition: alphaEqn.H:12
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
surfaceScalarField & phi
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
word alpharScheme("div(phirb,alpha)")
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
const volScalarField & vDotcAlphal
Definition: alphaEqn.H:8
label nAlphaCorr(readLabel(alphaControls.lookup("nAlphaCorr")))
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
rhoPhi
Definition: alphaEqn.H:115
volScalarField alpha10("alpha10", alpha1)
const volScalarField & alpha1
alpha2
Definition: alphaEqn.H:115
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal)
surfaceScalarField phic(mixture.cAlpha() *mag(phi/mesh.magSf()))
surfaceScalarField phir(IOobject("phir", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), mixture.cAlpha() *mag(phi/mesh.magSf()) *mixture.nHatf())
messageStream Info
const dimensionedScalar & rho2
Definition: createFields.H:44
zeroField divU
Definition: alphaSuSp.H:3
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
phaseChangeTwoPhaseMixture & mixture
Definition: createFields.H:38
const dimensionedScalar & rho1
Definition: createFields.H:43
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
bool MULESCorr(alphaControls.lookupOrDefault< Switch >("MULESCorr", false))
zeroField Sp
Definition: alphaSuSp.H:2