PDRFoamAutoRefine.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2015 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 \*---------------------------------------------------------------------------*/
57 
58 #include "fvCFD.H"
59 #include "dynamicFvMesh.H"
60 #include "psiuReactionThermo.H"
62 #include "laminarFlameSpeed.H"
63 #include "XiModel.H"
64 #include "PDRDragModel.H"
65 #include "ignition.H"
66 #include "Switch.H"
67 #include "bound.H"
68 #include "dynamicRefineFvMesh.H"
69 #include "pimpleControl.H"
70 
71 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
72 
73 int main(int argc, char *argv[])
74 {
75  #include "setRootCase.H"
76 
77  #include "createTime.H"
78  #include "createDynamicFvMesh.H"
79 
80  pimpleControl pimple(mesh);
81 
82  #include "readCombustionProperties.H"
83  #include "readGravitationalAcceleration.H"
84  #include "createFields.H"
85  #include "initContinuityErrs.H"
86  #include "createTimeControls.H"
87  #include "compressibleCourantNo.H"
88  #include "setInitialDeltaT.H"
89 
90  scalar StCoNum = 0.0;
91 
92  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
93 
94  Info<< "\nStarting time loop\n" << endl;
95 
96  bool hasChanged = false;
97 
98  while (runTime.run())
99  {
100  #include "createTimeControls.H"
101  #include "compressibleCourantNo.H"
102  #include "setDeltaT.H"
103 
104  // Indicators for refinement. Note: before runTime++
105  // only for postprocessing reasons.
106  tmp<volScalarField> tmagGradP = mag(fvc::grad(p));
107  volScalarField normalisedGradP
108  (
109  "normalisedGradP",
110  tmagGradP()/max(tmagGradP())
111  );
112  normalisedGradP.writeOpt() = IOobject::AUTO_WRITE;
113  tmagGradP.clear();
114 
115  runTime++;
116 
117  Info<< "\n\nTime = " << runTime.timeName() << endl;
118 
119  {
120  // Make the fluxes absolute
122 
123  // Test : disable refinement for some cells
124  PackedBoolList& protectedCell =
125  refCast<dynamicRefineFvMesh>(mesh).protectedCell();
126 
127  if (protectedCell.empty())
128  {
129  protectedCell.setSize(mesh.nCells());
130  protectedCell = 0;
131  }
132 
133  forAll(betav, cellI)
134  {
135  if (betav[cellI] < 0.99)
136  {
137  protectedCell[cellI] = 1;
138  }
139  }
140 
141  // Flux estimate for introduced faces.
142  volVectorField rhoU("rhoU", rho*U);
143 
144  // Do any mesh changes
145  bool meshChanged = mesh.update();
146 
147 
148  if (meshChanged)
149  {
150  hasChanged = true;
151  }
152 
153  if (runTime.write() && hasChanged)
154  {
155  betav.write();
156  Lobs.write();
157  CT.write();
158  drag->writeFields();
159  flameWrinkling->writeFields();
160  hasChanged = false;
161  }
162 
163  // Make the fluxes relative to the mesh motion
165  }
166 
167 
168  #include "rhoEqn.H"
169 
170  // --- Pressure-velocity PIMPLE corrector loop
171  while (pimple.loop())
172  {
173  #include "UEqn.H"
174 
175 
176  // --- Pressure corrector loop
177  while (pimple.correct())
178  {
179  #include "bEqn.H"
180  #include "ftEqn.H"
181  #include "huEqn.H"
182  #include "hEqn.H"
183 
184  if (!ign.ignited())
185  {
186  hu == h;
187  }
188 
189  #include "pEqn.H"
190  }
191 
192  if (pimple.turbCorr())
193  {
194  turbulence->correct();
195  }
196  }
197 
198  runTime.write();
199 
200  Info<< "\nExecutionTime = "
201  << runTime.elapsedCpuTime()
202  << " s\n" << endl;
203  }
204 
205  Info<< "\n end\n";
206 
207  return 0;
208 }
209 
210 
211 // ************************************************************************* //
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag
Definition: createFields.H:183
dimensioned< scalar > mag(const dimensioned< Type > &)
void makeAbsolute(surfaceScalarField &phi, const volVectorField &U)
Make the given flux absolute.
Definition: fvcMeshPhi.C:113
Bound the given scalar field if it has gone unbounded.
const volScalarField & betav
messageStream Info
autoPtr< compressible::turbulenceModel > turbulence
Definition: createFields.H:23
dynamicFvMesh & mesh
const dictionary & pimple
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
volScalarField & p
Definition: createFields.H:51
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:74
surfaceScalarField & phi
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
Read the control parameters used by setDeltaT.
int main(int argc, char *argv[])
Definition: postCalc.C:54
autoPtr< XiModel > flameWrinkling
Create the flame-wrinkling model.
Definition: createFields.H:193
U
Definition: pEqn.H:82
const dimensionedScalar h
Planck constant.
Definition: createFields.H:6