YEqn.H
Go to the documentation of this file.
1tmp<fv::convectionScheme<scalar>> mvConvection(nullptr);
2
3if (Y.size())
4{
5 mvConvection = tmp<fv::convectionScheme<scalar>>
6 (
7 fv::convectionScheme<scalar>::New
8 (
9 mesh,
10 fields,
11 phi,
12 mesh.divScheme("div(phi,Yi_h)")
13 )
14 );
15}
16
17{
18 reaction.correct();
19 Qdot = reaction.Qdot();
20 volScalarField Yt
21 (
22 IOobject("Yt", runTime.timeName(), mesh),
23 mesh,
24 dimensionedScalar("Yt", dimless, 0)
25 );
26
28 {
29 if (i != inertIndex && composition.active(i))
30 {
31 volScalarField& Yi = Y[i];
32
33 fvScalarMatrix YiEqn
34 (
35 fvm::ddt(rho, Yi)
36 + mvConvection->fvmDiv(phi, Yi)
37 - fvm::laplacian(turbulence.muEff(), Yi)
38 ==
39 reaction.R(Yi)
40 + fvOptions(rho, Yi)
41 );
42
43 YiEqn.relax();
44
45 fvOptions.constrain(YiEqn);
46
47 YiEqn.solve(mesh.solver("Yi"));
48
49 fvOptions.correct(Yi);
50
51 Yi.max(0.0);
52 Yt += Yi;
53 }
54 }
55
56 if (Y.size())
57 {
58 Y[inertIndex] = scalar(1) - Yt;
59 Y[inertIndex].max(0.0);
60 }
61}
fv::options & fvOptions
surfaceScalarField & phi
basicSpecieMixture & composition
PtrList< volScalarField > & Y
Qdot
Definition: YEqn.H:14
volScalarField Yt(0.0 *Y[0])
tmp< fv::convectionScheme< scalar > > mvConvection(fv::convectionScheme< scalar >::New(mesh, fields, phi, mesh.divScheme("div(phi,Yi_h)")))
dynamicFvMesh & mesh
engineTime & runTime
label inertIndex
compressible::turbulenceModel & turbulence
CombustionModel< rhoReactionThermo > & reaction
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333