PDRFoam.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  PDRFoam
26 
27 Description
28  Solver for compressible premixed/partially-premixed combustion with
29  turbulence modelling.
30 
31  Combusting RANS code using the b-Xi two-equation model.
32  Xi may be obtained by either the solution of the Xi transport
33  equation or from an algebraic exression. Both approaches are
34  based on Gulder's flame speed correlation which has been shown
35  to be appropriate by comparison with the results from the
36  spectral model.
37 
38  Strain effects are incorporated directly into the Xi equation
39  but not in the algebraic approximation. Further work need to be
40  done on this issue, particularly regarding the enhanced removal rate
41  caused by flame compression. Analysis using results of the spectral
42  model will be required.
43 
44  For cases involving very lean Propane flames or other flames which are
45  very strain-sensitive, a transport equation for the laminar flame
46  speed is present. This equation is derived using heuristic arguments
47  involving the strain time scale and the strain-rate at extinction.
48  the transport velocity is the same as that for the Xi equation.
49 
50  For large flames e.g. explosions additional modelling for the flame
51  wrinkling due to surface instabilities may be applied.
52 
53  PDR (porosity/distributed resistance) modelling is included to handle
54  regions containing blockages which cannot be resolved by the mesh.
55 
56  The fields used by this solver are:
57 
58  betav: Volume porosity
59  Lobs: Average diameter of obstacle in cell (m)
60  Aw: Obstacle surface area per unit volume (1/m)
61  CR: Drag tensor (1/m)
62  CT: Turbulence generation parameter (1/m)
63  Nv: Number of obstacles in cell per unit volume (m^-2)
64  nsv: Tensor whose diagonal indicates the number to subtract from
65  Nv to get the number of obstacles crossing the flow in each
66  direction.
67 
68 \*---------------------------------------------------------------------------*/
69 
70 #include "fvCFD.H"
71 #include "psiuReactionThermo.H"
73 #include "laminarFlameSpeed.H"
74 #include "XiModel.H"
75 #include "PDRDragModel.H"
76 #include "ignition.H"
77 #include "Switch.H"
78 #include "bound.H"
79 #include "pimpleControl.H"
80 #include "fvOptions.H"
81 
82 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
83 
84 int main(int argc, char *argv[])
85 {
86  #include "postProcess.H"
87 
88  #include "setRootCaseLists.H"
89  #include "createTime.H"
90  #include "createMesh.H"
91  #include "createControl.H"
92  #include "readCombustionProperties.H"
93  #include "readGravitationalAcceleration.H"
94  #include "createFields.H"
95  #include "createFieldRefs.H"
96  #include "initContinuityErrs.H"
97  #include "createTimeControls.H"
98  #include "compressibleCourantNo.H"
99  #include "setInitialDeltaT.H"
100 
101  turbulence->validate();
102  scalar StCoNum = 0.0;
103 
104  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
105 
106  Info<< "\nStarting time loop\n" << endl;
107 
108  while (runTime.run())
109  {
110  #include "readTimeControls.H"
111  #include "compressibleCourantNo.H"
112  #include "setDeltaT.H"
113 
114  runTime++;
115  Info<< "\n\nTime = " << runTime.timeName() << endl;
116 
117  #include "rhoEqn.H"
118 
119  // --- Pressure-velocity PIMPLE corrector loop
120  while (pimple.loop())
121  {
122  #include "UEqn.H"
123 
124  // --- Pressure corrector loop
125  while (pimple.correct())
126  {
127  #include "bEqn.H"
128  #include "ftEqn.H"
129  #include "EauEqn.H"
130  #include "EaEqn.H"
131 
132  if (!ign.ignited())
133  {
134  thermo.heu() == thermo.he();
135  }
136 
137  #include "pEqn.H"
138  }
139 
140  if (pimple.turbCorr())
141  {
142  turbulence->correct();
143  }
144  }
145 
146  runTime.write();
147 
148  Info<< "\nExecutionTime = "
149  << runTime.elapsedCpuTime()
150  << " s\n" << endl;
151  }
152 
153  Info<< "\n end\n";
154 
155  return 0;
156 }
157 
158 
159 // ************************************************************************* //
pimpleNoLoopControl & pimple
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Read the control parameters used by setDeltaT.
rhoReactionThermo & thermo
Definition: createFields.H:28
Bound the given scalar field if it has gone unbounded.
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
messageStream Info
Execute application functionObjects to post-process existing results.
Read the control parameters used by setDeltaT.