All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
denseParticleFoam.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) 2013-2021 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  denseParticleFoam
26 
27 Description
28  Transient solver for the coupled transport of particle clouds including the
29  effect of the volume fraction of particles on the continuous phase, with
30  optional mesh motion and mesh topology changes.
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #include "NamedEnum.H"
35 
36 namespace Foam
37 {
38  enum class cloudForceSplit
39  {
40  faceExplicitCellImplicit, // Implicit part of the cloud force added to
41  // the cell momentum equation. Explicit part
42  // to the face momentum equation. This is the
43  // least likely to create staggering patterns
44  // in the velocity field, but it can create
45  // unphysical perturbations in cell
46  // velocities even when particles and flow
47  // have the similar velocities.
48 
49  faceExplicitCellLagged, // Entire cloud force evaluated explicitly
50  // and added to the face momentum equation.
51  // Lagged correction (i.e.,
52  // fvm::Sp(cloudSU.diag(), Uc) -
53  // cloudSU.diag()*Uc) added to the cell
54  // momentum equation. This creates physical
55  // cell velocities when particles and flow
56  // have the same velocity, but can also
57  // result in staggering patterns in packed
58  // beds. Unsuitable for MPPIC.
59 
60  faceImplicit // Implicit and explicit parts of the force
61  // both added to the face momentum equation.
62  // Behaves somewhere between the other two.
63  };
64 
65  template<>
67  {
68  "faceExplicitCellImplicit",
69  "faceExplicitCellLagged",
70  "faceImplicit"
71  };
72 
73  const NamedEnum<cloudForceSplit, 3> cloudForceSplitNames;
74 }
75 
76 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77 
78 #include "fvCFD.H"
79 #include "dynamicFvMesh.H"
82 #include "pimpleControl.H"
83 #include "pressureReference.H"
84 #include "CorrectPhi.H"
85 #include "fvModels.H"
86 #include "fvConstraints.H"
87 #include "parcelCloudList.H"
88 
89 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
90 
91 int main(int argc, char *argv[])
92 {
93  #include "postProcess.H"
94 
95  #include "setRootCaseLists.H"
96  #include "createTime.H"
97  #include "createDynamicFvMesh.H"
98  #include "createDyMControls.H"
99  #include "createFields.H"
100  #include "createUcfIfPresent.H"
101  #include "initContinuityErrs.H"
102 
103  Info<< "\nStarting time loop\n" << endl;
104 
105  while (pimple.run(runTime))
106  {
107  #include "readDyMControls.H"
108  #include "CourantNo.H"
109  #include "setDeltaT.H"
110 
111  runTime++;
112 
113  Info<< "Time = " << runTime.timeName() << nl << endl;
114 
115  // Store the particle positions
116  clouds.storeGlobalPositions();
117 
118  mesh.update();
119 
120  if (mesh.changing())
121  {
122  if (correctPhi)
123  {
124  #include "correctPhic.H"
125  }
126 
127  if (checkMeshCourantNo)
128  {
129  #include "meshCourantNo.H"
130  }
131  }
132 
133  continuousPhaseTransport.correct();
134  muc = rhoc*continuousPhaseTransport.nu();
135 
136  clouds.evolve();
137 
138  // Update continuous phase volume fraction field
139  alphac = max(1.0 - clouds.theta(), alphacMin);
140  alphac.correctBoundaryConditions();
141  alphacf = fvc::interpolate(alphac);
142  alphaPhic = alphacf*phic;
143 
144  // Cloud forces
145  fvVectorMatrix cloudSU(clouds.SU(Uc));
146  volVectorField cloudSUu
147  (
148  IOobject
149  (
150  "cloudSUu",
151  runTime.timeName(),
152  mesh
153  ),
154  mesh,
156  zeroGradientFvPatchVectorField::typeName
157  );
158  volScalarField cloudSUp
159  (
160  IOobject
161  (
162  "cloudSUp",
163  runTime.timeName(),
164  mesh
165  ),
166  mesh,
168  zeroGradientFvPatchVectorField::typeName
169  );
170 
171  const cloudForceSplit cloudSUSplit =
172  pimple.dict().found("cloudForceSplit")
173  ? cloudForceSplitNames.read(pimple.dict().lookup("cloudForceSplit"))
174  : cloudForceSplit::faceExplicitCellImplicit;
175 
176  switch (cloudSUSplit)
177  {
178  case cloudForceSplit::faceExplicitCellImplicit:
179  cloudSUu.primitiveFieldRef() = -cloudSU.source()/mesh.V();
180  cloudSUu.correctBoundaryConditions();
181  cloudSUp.primitiveFieldRef() = Zero;
182  cloudSUp.correctBoundaryConditions();
183  //cloudSU.diag() = cloudSU.diag();
184  cloudSU.source() = Zero;
185  break;
186 
187  case cloudForceSplit::faceExplicitCellLagged:
188  cloudSUu.primitiveFieldRef() =
189  (cloudSU.diag()*Uc() - cloudSU.source())/mesh.V();
190  cloudSUu.correctBoundaryConditions();
191  cloudSUp.primitiveFieldRef() = Zero;
192  cloudSUp.correctBoundaryConditions();
193  //cloudSU.diag() = cloudSU.diag();
194  cloudSU.source() = cloudSU.diag()*Uc();
195  break;
196 
197  case cloudForceSplit::faceImplicit:
198  cloudSUu.primitiveFieldRef() = -cloudSU.source()/mesh.V();
199  cloudSUu.correctBoundaryConditions();
200  cloudSUp.primitiveFieldRef() = cloudSU.diag()/mesh.V();
201  cloudSUp.correctBoundaryConditions();
202  cloudSU.diag() = Zero;
203  cloudSU.source() = Zero;
204  break;
205  }
206 
207  // --- Pressure-velocity PIMPLE corrector loop
208  while (pimple.loop())
209  {
210  fvModels.correct();
211 
212  #include "UcEqn.H"
213 
214  // --- PISO loop
215  while (pimple.correct())
216  {
217  #include "pEqn.H"
218  }
219 
220  if (pimple.turbCorr())
221  {
222  continuousPhaseTurbulence->correct();
223  }
224  }
225 
226  runTime.write();
227 
228  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
229  << " ClockTime = " << runTime.elapsedClockTime() << " s"
230  << nl << endl;
231  }
232 
233  Info<< "End\n" << endl;
234 
235  return 0;
236 }
237 
238 
239 // ************************************************************************* //
pimpleNoLoopControl & pimple
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void correct()
Correct the fvModels.
Definition: fvModels.C:316
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
correctPhi
checkMeshCourantNo
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
Info<< "Reading field U\"<< endl;volVectorField Uc(IOobject(IOobject::groupName("U", continuousPhaseName), runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field p\"<< endl;volScalarField p(IOobject("p", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading/calculating continuous-phase face flux field phic\"<< endl;surfaceScalarField phic(IOobject(IOobject::groupName("phi", continuousPhaseName), runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), linearInterpolate(Uc) &mesh.Sf());pressureReference pressureReference(p, pimple.dict());mesh.setFluxRequired(p.name());Info<< "Creating turbulence model\"<< endl;singlePhaseTransportModel continuousPhaseTransport(Uc, phic);dimensionedScalar rhocValue(IOobject::groupName("rho", continuousPhaseName), dimDensity, continuousPhaseTransport.lookup(IOobject::groupName("rho", continuousPhaseName)));volScalarField rhoc(IOobject(rhocValue.name(), runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, rhocValue);volScalarField muc(IOobject(IOobject::groupName("mu", continuousPhaseName), runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), rhoc *continuousPhaseTransport.nu());Info<< "Creating field alphac\"<< endl;volScalarField alphac(IOobject(IOobject::groupName("alpha", continuousPhaseName), runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimless, 0));Info<< "Constructing clouds"<< endl;parcelCloudList clouds(rhoc, Uc, muc, g);scalar alphacMin(1 - mesh.solverDict(alphac.name()).lookup< scalar >"max"));alphac=max(1.0 - clouds.theta(), alphacMin);alphac.correctBoundaryConditions();surfaceScalarField alphacf("alphacf", fvc::interpolate(alphac));surfaceScalarField alphaPhic(IOobject::groupName("alphaPhi", continuousPhaseName), alphacf *phic);autoPtr< phaseIncompressible::momentumTransportModel > continuousPhaseTurbulence(phaseIncompressible::momentumTransportModel::New(alphac, Uc, alphaPhic, phic, continuousPhaseTransport))
const char * names[]
Definition: EDC.C:32
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Creates and initialises the continuous phase face velocity field Ufc if required. ...
dynamicFvMesh & mesh
phic
Definition: correctPhic.H:2
static const zero Zero
Definition: zero.H:97
const dimensionSet dimForce
static const char nl
Definition: Ostream.H:260
const dimensionSet dimVelocity
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Foam::fvModels & fvModels
Calculates and outputs the mean and maximum Courant Numbers.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
const dimensionSet dimVolume
messageStream Info
Execute application functionObjects to post-process existing results.
Namespace for OpenFOAM.