alphaEqn.H
Go to the documentation of this file.
1 {
2  #include "alphaScheme.H"
3 
4  Pair<tmp<volScalarField>> vDotAlphal = mixture.vDotAlphal();
5  const volScalarField& vDotcAlphal = vDotAlphal[0]();
6  const volScalarField& vDotvAlphal = vDotAlphal[1]();
7  const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal);
8 
9  tmp<surfaceScalarField> talphaPhi;
10 
11  if (MULESCorr)
12  {
13  fvScalarMatrix alpha1Eqn
14  (
15  fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
16  + fv::gaussConvectionScheme<scalar>
17  (
18  mesh,
19  phi,
20  upwind<scalar>(mesh, phi)
21  ).fvmDiv(phi, alpha1)
22  - fvm::Sp(divU, alpha1)
23  ==
25  + vDotcAlphal
26  );
27 
28  alpha1Eqn.solve();
29 
30  Info<< "Phase-1 volume fraction = "
31  << alpha1.weightedAverage(mesh.Vsc()).value()
32  << " Min(" << alpha1.name() << ") = " << min(alpha1).value()
33  << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
34  << endl;
35 
36  talphaPhi = alpha1Eqn.flux();
37  }
38 
39  volScalarField alpha10("alpha10", alpha1);
40 
41  for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
42  {
43  tmp<surfaceScalarField> talphaPhiCorr
44  (
45  fvc::flux
46  (
47  phi,
48  alpha1,
49  compressionScheme.rewind()
50  )
51  );
52 
53  if (MULESCorr)
54  {
55  talphaPhiCorr.ref() -= talphaPhi();
56 
57  volScalarField alpha100("alpha100", alpha10);
58  alpha10 = alpha1;
59 
61  (
62  geometricOneField(),
63  alpha1,
64  talphaPhi(),
65  talphaPhiCorr.ref(),
67  (
68  divU*(alpha10 - alpha100)
69  - vDotvmcAlphal*alpha10
70  )(),
71  oneField(),
72  zeroField()
73  );
74 
75  // Under-relax the correction for all but the 1st corrector
76  if (aCorr == 0)
77  {
78  talphaPhi.ref() += talphaPhiCorr();
79  }
80  else
81  {
82  alpha1 = 0.5*alpha1 + 0.5*alpha10;
83  talphaPhi.ref() += 0.5*talphaPhiCorr();
84  }
85  }
86  else
87  {
89  (
90  geometricOneField(),
91  alpha1,
92  phi,
93  talphaPhiCorr.ref(),
95  (divU*alpha1 + vDotcAlphal)(),
96  oneField(),
97  zeroField()
98  );
99 
100  talphaPhi = talphaPhiCorr;
101  }
102 
103  alpha2 = 1.0 - alpha1;
104  }
105 
107 
108  Info<< "Liquid phase volume fraction = "
109  << alpha1.weightedAverage(mesh.V()).value()
110  << " Min(" << alpha1.name() << ") = " << min(alpha1).value()
111  << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
112  << endl;
113 }
ITstream compressionScheme(compressionSchemes.found(alphaScheme) ? mesh.divScheme(divAlphaName) :ITstream(divAlphaName, tokenList { word("Gauss"), word("interfaceCompression"), alphaScheme, alphaControls.lookup< scalar >("cAlpha") }))
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
volScalarField alpha10("alpha10", alpha1)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
tmp< surfaceScalarField > talphaPhi
Definition: alphaEqn.H:9
alpha2
Definition: alphaEqn.H:115
phi
Definition: pEqn.H:104
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dynamicFvMesh & mesh
rhoPhi
Definition: alphaEqn.H:106
volScalarField & alpha1(mixture.alpha1())
const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal)
const volScalarField & vDotcAlphal
Definition: alphaEqn.H:5
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-8/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 > &)
messageStream Info
const label nAlphaCorr(alphaControls.lookup< label >("nAlphaCorr"))
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
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)
const bool MULESCorr(alphaControls.lookupOrDefault< Switch >("MULESCorr", false))
zeroField Sp
Definition: alphaSuSp.H:2