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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Application
27  PDRFoam
28 
29 Description
30  Solver for compressible premixed/partially-premixed combustion with
31  turbulence modelling.
32 
33  Combusting RANS code using the b-Xi two-equation model.
34  Xi may be obtained by either the solution of the Xi transport
35  equation or from an algebraic expression. Both approaches are
36  based on Gulder's flame speed correlation which has been shown
37  to be appropriate by comparison with the results from the
38  spectral model.
39 
40  Strain effects are incorporated directly into the Xi equation
41  but not in the algebraic approximation. Further work need to be
42  done on this issue, particularly regarding the enhanced removal rate
43  caused by flame compression. Analysis using results of the spectral
44  model will be required.
45 
46  For cases involving very lean Propane flames or other flames which are
47  very strain-sensitive, a transport equation for the laminar flame
48  speed is present. This equation is derived using heuristic arguments
49  involving the strain time scale and the strain-rate at extinction.
50  the transport velocity is the same as that for the Xi equation.
51 
52  For large flames e.g. explosions additional modelling for the flame
53  wrinkling due to surface instabilities may be applied.
54 
55  PDR (porosity/distributed resistance) modelling is included to handle
56  regions containing blockages which cannot be resolved by the mesh.
57 
58 \*---------------------------------------------------------------------------*/
59 
60 #include "fvCFD.H"
61 #include "dynamicFvMesh.H"
62 #include "psiuReactionThermo.H"
64 #include "laminarFlameSpeed.H"
65 #include "XiModel.H"
66 #include "PDRDragModel.H"
67 #include "ignition.H"
68 #include "bound.H"
69 #include "dynamicRefineFvMesh.H"
70 #include "pimpleControl.H"
71 
72 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
73 
74 int main(int argc, char *argv[])
75 {
76  argList::addNote
77  (
78  "Solver for compressible premixed/partially-premixed combustion with"
79  " turbulence modelling."
80  );
81 
82  #include "setRootCaseLists.H"
83  #include "createTime.H"
84  #include "createDynamicFvMesh.H"
85 
86  pimpleControl pimple(mesh);
87 
88  #include "readCombustionProperties.H"
89  #include "readGravitationalAcceleration.H"
90  #include "createFields.H"
91  #include "initContinuityErrs.H"
92  #include "createTimeControls.H"
93  #include "compressibleCourantNo.H"
94  #include "setInitialDeltaT.H"
95 
96  turbulence->validate();
97  scalar StCoNum = 0.0;
98 
99  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
100 
101  Info<< "\nStarting time loop\n" << endl;
102 
103  bool hasChanged = false;
104 
105  while (runTime.run())
106  {
107  #include "readTimeControls.H"
108  #include "compressibleCourantNo.H"
109  #include "setDeltaT.H"
110 
111  // Indicators for refinement.
112  // Note: before ++runTime only for post-processing reasons.
113  tmp<volScalarField> tmagGradP = mag(fvc::grad(p));
114  volScalarField normalisedGradP
115  (
116  "normalisedGradP",
117  tmagGradP()/max(tmagGradP())
118  );
119  normalisedGradP.writeOpt() = IOobject::AUTO_WRITE;
120  tmagGradP.clear();
121 
122  ++runTime;
123 
124  Info<< "\n\nTime = " << runTime.timeName() << endl;
125 
126  {
127  // Make the fluxes absolute
129 
130  // Test : disable refinement for some cells
131  bitSet& protectedCell =
132  refCast<dynamicRefineFvMesh>(mesh).protectedCell();
133 
134  if (protectedCell.empty())
135  {
136  protectedCell.setSize(mesh.nCells());
137  protectedCell = false;
138  }
139 
140  forAll(betav, celli)
141  {
142  if (betav[celli] < 0.99)
143  {
144  protectedCell.set(celli);
145  }
146  }
147 
148  // Flux estimate for introduced faces.
149  volVectorField rhoU("rhoU", rho*U);
150 
151  // Do any mesh changes
152  bool meshChanged = mesh.update();
153 
154 
155  if (meshChanged)
156  {
157  hasChanged = true;
158  }
159 
160  if (runTime.write() && hasChanged)
161  {
162  betav.write();
163  Lobs.write();
164  CT.write();
165  drag->writeFields();
166  flameWrinkling->writeFields();
167  hasChanged = false;
168  }
169 
170  // Make the fluxes relative to the mesh motion
172  }
173 
174 
175  #include "rhoEqn.H"
176 
177  // --- Pressure-velocity PIMPLE corrector loop
178  while (pimple.loop())
179  {
180  #include "UEqn.H"
181 
182 
183  // --- Pressure corrector loop
184  while (pimple.correct())
185  {
186  #include "bEqn.H"
187  #include "ftEqn.H"
188  #include "huEqn.H"
189  #include "hEqn.H"
190 
191  if (!ign.ignited())
192  {
193  hu == h;
194  }
195 
196  #include "pEqn.H"
197  }
198 
199  if (pimple.turbCorr())
200  {
201  turbulence->correct();
202  }
203  }
204 
205  runTime.write();
206 
207  runTime.printExecutionTime(Info);
208  }
209 
210  Info<< "End\n";
211 
212  return 0;
213 }
214 
215 
216 // ************************************************************************* //
runTime
engineTime & runTime
Definition: createEngineTime.H:13
drag
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:165
psiuReactionThermo.H
p
volScalarField & p
Definition: createFieldRefs.H:8
turbulence
Info<< "Reading field U\n"<< 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\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
Foam::fac::grad
tmp< GeometricField< typename outerProduct< vector, Type >::type, faPatchField, areaMesh >> grad(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Definition: facGrad.C:56
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
laminarFlameSpeed.H
rho
rho
Definition: readInitialConditions.H:88
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
pimpleControl.H
Foam::fvc::makeRelative
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:77
pimple
pimpleControl & pimple
Definition: setRegionFluidFields.H:56
setRootCaseLists.H
Foam::constant::universal::h
const dimensionedScalar h
Planck constant.
Definition: setRegionSolidFields.H:33
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
bound.H
Bound the given scalar field if it has gone unbounded.
ignition.H
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
dynamicRefineFvMesh.H
betav
const volScalarField & betav
Definition: setRegionSolidFields.H:35
flameWrinkling
autoPtr< XiModel > flameWrinkling
Create the flame-wrinkling model.
Definition: createFields.H:175
Foam::volVectorField
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:62
XiModel.H
U
U
Definition: pEqn.H:72
createTimeControls.H
Read the control parameters used by setDeltaT.
readTimeControls.H
Read the control parameters used by setDeltaT.
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
createTime.H
dynamicFvMesh.H
fvCFD.H
PDRDragModel.H
Foam::fvc::makeAbsolute
void makeAbsolute(surfaceScalarField &phi, const volVectorField &U)
Make the given flux absolute.
Definition: fvcMeshPhi.C:116
createDynamicFvMesh.H
turbulentFluidThermoModel.H