EauEqn.H
Go to the documentation of this file.
1if (ign.ignited())
2{
3 volScalarField& heau = thermo.heu();
4
5 fvScalarMatrix heauEqn
6 (
7 fvm::ddt(rho, heau) + mvConvection->fvmDiv(phi, heau)
8 + (fvc::ddt(rho, K) + fvc::div(phi, K))*rho/thermo.rhou()
9 + (
10 heau.name() == "eau"
11 ? fvc::div
12 (
13 fvc::absolute(phi/fvc::interpolate(rho), U),
14 p,
15 "div(phiv,p)"
16 )*rho/thermo.rhou()
17 : -dpdt*rho/thermo.rhou()
18 )
19 - fvm::laplacian(turbulence->alphaEff(), heau)
20
21 // These terms cannot be used in partially-premixed combustion due to
22 // the resultant inconsistency between ft and heau transport.
23 // A possible solution would be to solve for ftu as well as ft.
24 //- fvm::div(muEff*fvc::grad(b)/(b + 0.001), heau)
25 //+ fvm::Sp(fvc::div(muEff*fvc::grad(b)/(b + 0.001)), heau)
26
27 ==
28 fvOptions(rho, heau)
29 );
30
31 fvOptions.constrain(heauEqn);
32
33 heauEqn.solve();
34
35 fvOptions.correct(heau);
36}
CGAL::Exact_predicates_exact_constructions_kernel K
tmp< fv::convectionScheme< scalar > > mvConvection(fv::convectionScheme< scalar >::New(mesh, fields, phi, mesh.divScheme("div(phi,Yi_h)")))
fv::options & fvOptions
surfaceScalarField & phi
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:72
volScalarField & p
volScalarField & dpdt
compressible::turbulenceModel & turbulence