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-2023 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 expression. 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 "argList.H"
71 #include "timeSelector.H"
76 #include "laminarFlameSpeed.H"
77 #include "XiModel.H"
78 #include "PDRDragModel.H"
79 #include "ignition.H"
80 #include "pimpleControl.H"
81 #include "pressureReference.H"
82 #include "findRefCell.H"
83 #include "constrainPressure.H"
84 #include "constrainHbyA.H"
85 #include "adjustPhi.H"
87 #include "fvModels.H"
88 #include "fvConstraints.H"
89 
90 #include "fvcDdt.H"
91 #include "fvcGrad.H"
92 #include "fvcFlux.H"
93 #include "fvcReconstruct.H"
94 #include "fvcMeshPhi.H"
95 
96 #include "fvmDdt.H"
97 #include "fvmDiv.H"
98 #include "fvmLaplacian.H"
99 
100 using namespace Foam;
101 
102 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
103 
104 int main(int argc, char *argv[])
105 {
106  #include "postProcess.H"
107 
108  #include "setRootCase.H"
109  #include "createTime.H"
110  #include "createMesh.H"
111  #include "createControl.H"
112  #include "readCombustionProperties.H"
113  #include "readGravitationalAcceleration.H"
114  #include "createFields.H"
115  #include "createFieldRefs.H"
116  #include "initContinuityErrs.H"
117  #include "createTimeControls.H"
118  #include "compressibleCourantNo.H"
119  #include "setInitialDeltaT.H"
120 
121  turbulence->validate();
122  scalar StCoNum = 0.0;
123 
124  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
125 
126  Info<< "\nStarting time loop\n" << endl;
127 
128  while (pimple.run(runTime))
129  {
130  #include "readTimeControls.H"
131  #include "compressibleCourantNo.H"
132  #include "setDeltaT.H"
133 
134  runTime++;
135  Info<< "\n\nTime = " << runTime.name() << endl;
136 
137  #include "rhoEqn.H"
138 
139  // --- Pressure-velocity PIMPLE corrector loop
140  while (pimple.loop())
141  {
142  fvModels.correct();
143 
144  if (pimple.predictTransport())
145  {
146  turbulence->predict();
147  thermophysicalTransport.predict();
148  }
149 
150  #include "UEqn.H"
151 
152  // --- Pressure corrector loop
153  while (pimple.correct())
154  {
155  #include "bEqn.H"
156  #include "ftEqn.H"
157  #include "EauEqn.H"
158  #include "EaEqn.H"
159 
160  if (!ign.ignited())
161  {
162  thermo.heu() == thermo.he();
163  }
164 
165  #include "pEqn.H"
166  }
167 
168  if (pimple.correctTransport())
169  {
170  turbulence->correct();
171  thermophysicalTransport.correct();
172  }
173  }
174 
175  runTime.write();
176 
177  Info<< "\nExecutionTime = "
178  << runTime.elapsedCpuTime()
179  << " s\n" << endl;
180  }
181 
182  Info<< "\n end\n";
183 
184  return 0;
185 }
186 
187 
188 // ************************************************************************* //
int main(int argc, char *argv[])
Definition: PDRFoam.C:101
StCoNum
Definition: StCourantNo.H:36
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
Solve the continuity for density.
Reset the timestep to maintain a constant maximum courant Number. Reduction of time-step is immediate...
virtual const volScalarField & he() const =0
Enthalpy/Internal energy [J/kg].
Finite volume models.
Definition: fvModels.H:65
virtual void correct()
Correct the fvModels.
Definition: fvModels.C:380
bool run(Time &time)
Time run loop.
bool correctTransport() const
Flag to indicate whether to correct the transport models.
bool correct()
Piso loop within outer loop.
bool predictTransport() const
Flag to indicate whether to predict the transport models.
Calculates and outputs the mean and maximum Courant Numbers.
pimpleControl pimple(mesh)
Read the control parameters used by setDeltaT.
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic,...
Calculate the first temporal derivative.
Calculate the face-flux of the given field.
Calculate the gradient of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Reconstruct volField from a face flux field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
Declare and initialise the cumulative continuity error.
Info<< "Creating thermophysical transport model\n"<< endl;turbulenceThermophysicalTransportModels::unityLewisEddyDiffusivity< RASThermophysicalTransportModel< ThermophysicalTransportModel< compressibleMomentumTransportModel, fluidThermo > >> thermophysicalTransport(turbulence(), thermo, true)
Namespace for OpenFOAM.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
Execute application functionObjects to post-process existing results.
Read the control parameters used by setDeltaT.
Set the initial timestep corresponding to the timestep adjustment algorithm in setDeltaT but only if ...
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))
fluidMulticomponentThermo & thermo
Definition: createFields.H:31