YEqn.H
Go to the documentation of this file.
1 tmp<fv::convectionScheme<scalar>> mvConvection
2 (
4  (
5  mesh,
6  fields,
7  phi,
8  mesh.divScheme("div(phi,Yi_h)")
9  )
10 );
11 
12 {
13  combustion->correct();
14  dQ = combustion->dQ();
15  label inertIndex = -1;
16  volScalarField Yt(0.0*Y[0]);
17 
18  forAll(Y, i)
19  {
20  if (Y[i].name() != inertSpecie)
21  {
22  volScalarField& Yi = Y[i];
23 
24  fvScalarMatrix YEqn
25  (
26  mvConvection->fvmDiv(phi, Yi)
27  - fvm::laplacian(turbulence->muEff(), Yi)
28  ==
29  parcels.SYi(i, Yi)
30  + combustion->R(Yi)
31  + fvOptions(rho, Yi)
32  );
33 
34  YEqn.relax();
35 
36  fvOptions.constrain(YEqn);
37 
38  YEqn.solve(mesh.solver("Yi"));
39 
40  fvOptions.correct(Yi);
41 
42  Yi.max(0.0);
43  Yt += Yi;
44  }
45  else
46  {
47  inertIndex = i;
48  }
49  }
50 
51  Y[inertIndex] = scalar(1) - Yt;
52  Y[inertIndex].max(0.0);
53 }
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
surfaceScalarField & phi
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Info<< "Creating combustion model\n"<< endl;autoPtr< combustionModels::psiCombustionModel > combustion(combustionModels::psiCombustionModel::New(mesh))
dQ
Definition: YEqn.H:14
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
fv::options & fvOptions
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
Y[inertIndex]
Definition: YEqn.H:51
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label inertIndex
Definition: YEqn.H:15
forAll(Y, i)
Definition: YEqn.H:18
tmp< fv::convectionScheme< scalar > > mvConvection(fv::convectionScheme< scalar >::New(mesh, fields, phi, mesh.divScheme("div(phi,Yi_h)")))
const word inertSpecie(thermo.lookup("inertSpecie"))
volScalarField Yt(0.0 *Y[0])