All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
YEEqn.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  radiation->correct();
13  combustion->correct();
14  volScalarField Yt(0.0*Y[0]);
15 
16  forAll(Y, i)
17  {
18  if (i != inertIndex && composition.active(i))
19  {
20  volScalarField& Yi = Y[i];
21 
22  fvScalarMatrix YiEqn
23  (
24  fvm::ddt(rho, Yi)
25  + mvConvection->fvmDiv(phi, Yi)
26  + thermophysicalTransport->divj(Yi)
27  ==
28  parcels.SYi(i, Yi)
29  + surfaceFilm.Srho(i)
30  + combustion->R(Yi)
31  + fvOptions(rho, Yi)
32  );
33 
34  YiEqn.relax();
35 
36  fvOptions.constrain(YiEqn);
37 
38  YiEqn.solve("Yi");
39 
40  fvOptions.correct(Yi);
41 
42  Yi.max(0.0);
43  Yt += Yi;
44  }
45  }
46 
47  Y[inertIndex] = scalar(1) - Yt;
48  Y[inertIndex].max(0.0);
49 
50  volScalarField& he = thermo.he();
51 
53  (
54  fvm::ddt(rho, he) + mvConvection->fvmDiv(phi, he)
55  + fvc::ddt(rho, K) + fvc::div(phi, K)
56  + (
57  he.name() == "e"
58  ? fvc::div
59  (
61  p,
62  "div(phiv,p)"
63  )
64  : -dpdt
65  )
66  + thermophysicalTransport->divq(he)
67  ==
68  combustion->Qdot()
69  + radiation->Sh(thermo, he)
70  + parcels.Sh(he)
71  + surfaceFilm.Sh()
72  + fvOptions(rho, he)
73  );
74 
75  EEqn.relax();
76 
77  fvOptions.constrain(EEqn);
78 
79  EEqn.solve();
80 
81  fvOptions.correct(he);
82 
83  thermo.correct();
84 
85  Info<< "min/max(T) = "
86  << min(T).value() << ", " << max(T).value() << endl;
87 }
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
fv::options & fvOptions
Y[inertIndex]
Definition: YEEqn.H:47
basicSpecieMixture & composition
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Y [inertIndex] max(0.0)
rhoReactionThermo & thermo
Definition: createFields.H:28
CGAL::Exact_predicates_exact_constructions_kernel K
scalar Qdot
Definition: solveChemistry.H:2
phi
Definition: pEqn.H:104
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:103
label inertIndex
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:57
volScalarField & dpdt
dynamicFvMesh & mesh
rhoReactionThermophysicalTransportModel & thermophysicalTransport
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
forAll(Y, i)
Definition: YEEqn.H:16
regionModels::surfaceFilmModel & surfaceFilm
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const volScalarField & T
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Info<< "Creating combustion model\"<< endl;autoPtr< CombustionModel< psiReactionThermo > > combustion(CombustionModel< psiReactionThermo >::New(thermo, turbulence()))
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
U
Definition: pEqn.H:72
tmp< fv::convectionScheme< scalar > > mvConvection(fv::convectionScheme< scalar >::New(mesh, fields, phi, mesh.divScheme("div(phi,Yi_h)")))
autoPtr< radiationModel > radiation(radiationModel::New(T))
messageStream Info
fvScalarMatrix EEqn(fvm::ddt(rho, he)+mvConvection->fvmDiv(phi, he)+fvc::ddt(rho, K)+fvc::div(phi, K)+(he.name()=="e" ? fvc::div(fvc::absolute(phi/fvc::interpolate(rho), U), p, "div(phiv,p)") :-dpdt)+thermophysicalTransport->divq(he)==combustion->Qdot()+radiation->Sh(thermo, he)+parcels.Sh(he)+surfaceFilm.Sh()+fvOptions(rho, he))
volScalarField & p
volScalarField Yt(0.0 *Y[0])