UEqns.H
Go to the documentation of this file.
1 #include "MRFCorrectBCs.H"
2 
3 PtrList<fvVectorMatrix> UEqns(fluid.phases().size());
4 autoPtr<multiphaseSystem::dragCoeffFields> dragCoeffs(fluid.dragCoeffs());
5 
6 int phasei = 0;
7 forAllIter(PtrDictionary<phaseModel>, fluid.phases(), iter)
8 {
9  phaseModel& phase = iter();
10  const volScalarField& alpha = phase;
11  volVectorField& U = phase.U();
12 
13  volScalarField nuEff(turbulence->nut() + iter().nu());
14 
15  UEqns.set
16  (
17  phasei,
18  new fvVectorMatrix
19  (
20  fvm::ddt(alpha, U)
21  + fvm::div(phase.alphaPhi(), U)
22 
23  + (alpha/phase.rho())*fluid.Cvm(phase)*
24  (
25  fvm::ddt(U)
26  + fvm::div(phase.phi(), U)
27  - fvm::Sp(fvc::div(phase.phi()), U)
28  )
29 
30  - fvm::laplacian(alpha*nuEff, U)
31  - fvc::div
32  (
33  alpha*(nuEff*dev(T(fvc::grad(U))) /*- ((2.0/3.0)*I)*k*/),
34  "div(Rc)"
35  )
36  ==
37  //- fvm::Sp(fluid.dragCoeff(phase, dragCoeffs())/phase.rho(), U)
38  //- (alpha*phase.rho())*fluid.lift(phase)
39  //+
40  (alpha/phase.rho())*fluid.Svm(phase)
41  - fvm::Sp
42  (
43  slamDampCoeff
44  *max
45  (
46  mag(U()) - maxSlamVelocity,
48  )
49  /pow(mesh.V(), 1.0/3.0),
50  U
51  )
52  )
53  );
54  MRF.addAcceleration
55  (
56  alpha*(1 + (1/phase.rho())*fluid.Cvm(phase)),
57  UEqns[phasei]
58  );
59  //UEqns[phasei].relax();
60 
61  phasei++;
62 }
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
U
Definition: pEqn.H:83
multiphaseSystem & fluid
Definition: createFields.H:10
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
autoPtr< multiphaseSystem::dragCoeffFields > dragCoeffs(fluid.dragCoeffs())
tmp< GeometricField< Type, fvPatchField, volMesh > > Sp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:67
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
IOMRFZoneList & MRF
forAllIter(PtrDictionary< phaseModel >, fluid.phases(), iter)
Definition: UEqns.H:7
const volScalarField & T
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
PtrList< fvVectorMatrix > UEqns(fluid.phases().size())
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
dimensioned< scalar > mag(const dimensioned< Type > &)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const dimensionSet dimVelocity