alphaEqn.H
Go to the documentation of this file.
1 {
2  word alphaScheme("div(phi,alpha)");
3  word alpharScheme("div(phirb,alpha)");
4 
5  if (MULESCorr)
6  {
7  fvScalarMatrix alpha1Eqn
8  (
9  fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
10  + fv::gaussConvectionScheme<scalar>
11  (
12  mesh,
13  phi,
14  upwind<scalar>(mesh, phi)
15  ).fvmDiv(phi, alpha1)
16  );
17 
18  solve(alpha1Eqn);
19 
20  Info<< "Phase-1 volume fraction = "
21  << alpha1.weightedAverage(mesh.Vsc()).value()
22  << " Min(" << alpha1.name() << ") = " << min(alpha1).value()
23  << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
24  << endl;
25 
26  tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
27  alphaPhi = talphaPhiUD();
28 
29  if (alphaApplyPrevCorr && talphaPhiCorr0.valid())
30  {
31  Info<< "Applying the previous iteration correction flux" << endl;
32 
34  (
35  geometricOneField(),
36  alpha1,
37  alphaPhi,
38  talphaPhiCorr0.ref(),
39  UniformField<scalar>(mixture.alphaMax()),
40  zeroField()
41  );
42 
43  alphaPhi += talphaPhiCorr0();
44  }
45 
46  // Cache the upwind-flux
47  talphaPhiCorr0 = talphaPhiUD;
48  }
49 
50  for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
51  {
52  tmp<surfaceScalarField> talphaPhiUn
53  (
54  fvc::flux
55  (
56  phi,
57  alpha1,
59  )
60  + fvc::flux
61  (
62  phir,
63  alpha1,
65  )
66  );
67 
68  if (MULESCorr)
69  {
70  tmp<surfaceScalarField> talphaPhiCorr(talphaPhiUn() - alphaPhi);
71  volScalarField alpha10("alpha10", alpha1);
72 
74  (
75  geometricOneField(),
76  alpha1,
77  talphaPhiUn(),
78  talphaPhiCorr.ref(),
79  UniformField<scalar>(mixture.alphaMax()),
80  zeroField()
81  );
82 
83  // Under-relax the correction for all but the 1st corrector
84  if (aCorr == 0)
85  {
86  alphaPhi += talphaPhiCorr();
87  }
88  else
89  {
90  alpha1 = 0.5*alpha1 + 0.5*alpha10;
91  alphaPhi += 0.5*talphaPhiCorr();
92  }
93  }
94  else
95  {
96  alphaPhi = talphaPhiUn;
97 
99  (
100  geometricOneField(),
101  alpha1,
102  phi,
103  alphaPhi,
104  UniformField<scalar>(mixture.alphaMax()),
105  zeroField()
106  );
107  }
108  }
109 
111  {
113  }
114 
115  alpha2 = 1.0 - alpha1;
116 
117  Info<< "Phase-1 volume fraction = "
118  << alpha1.weightedAverage(mesh.Vsc()).value()
119  << " Min(" << alpha1.name() << ") = " << min(alpha1).value()
120  << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
121  << endl;
122 }
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
volScalarField alpha10("alpha10", alpha1)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const bool alphaApplyPrevCorr(alphaControls.lookupOrDefault< Switch >("alphaApplyPrevCorr", false))
tmp< surfaceScalarField > talphaPhiCorr0
Definition: createFields.H:132
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
rhoEqn solve()
alpha2
Definition: alphaEqn.H:115
phi
Definition: pEqn.H:104
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dynamicFvMesh & mesh
volScalarField & alpha1(mixture.alpha1())
surfaceScalarField alphaPhi(phi.name()+alpha1.name(), fvc::flux(phi, alpha1, alphaScheme))
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 > &)
const word alphaScheme(mesh.divScheme(divAlphaName)[1].wordToken())
word alpharScheme("div(phirb,alpha)")
messageStream Info
const label nAlphaCorr(alphaControls.lookup< label >("nAlphaCorr"))
surfaceScalarField phir(IOobject("phir", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), mixture.cAlpha() *mag(phi/mesh.magSf()) *mixture.nHatf())
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
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))