rhoPimpleDyMFoam.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-2017 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  rhoPimpleDyMFoam
26 
27 Description
28  Transient solver for turbulent flow of compressible fluids for HVAC and
29  similar applications, with optional mesh motion and mesh topology changes.
30 
31  Uses the flexible PIMPLE (PISO-SIMPLE) solution for time-resolved and
32  pseudo-transient simulations.
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #include "fvCFD.H"
37 #include "dynamicFvMesh.H"
38 #include "fluidThermo.H"
40 #include "bound.H"
41 #include "pimpleControl.H"
42 #include "pressureControl.H"
43 #include "CorrectPhi.H"
44 #include "fvOptions.H"
45 #include "localEulerDdtScheme.H"
46 #include "fvcSmooth.H"
47 
48 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49 
50 int main(int argc, char *argv[])
51 {
52  #include "postProcess.H"
53 
54  #include "setRootCase.H"
55  #include "createTime.H"
56  #include "createDynamicFvMesh.H"
57  #include "createControl.H"
58  #include "initContinuityErrs.H"
59  #include "createFields.H"
60  #include "createFieldRefs.H"
61  #include "createFvOptions.H"
62  #include "createRhoUf.H"
63  #include "createControls.H"
64 
65  turbulence->validate();
66 
67  if (!LTS)
68  {
69  #include "compressibleCourantNo.H"
70  #include "setInitialDeltaT.H"
71  }
72 
73  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
74 
75  Info<< "\nStarting time loop\n" << endl;
76 
77  while (runTime.run())
78  {
79  #include "readControls.H"
80 
81  {
82  // Store divrhoU from the previous mesh so that it can be mapped
83  // and used in correctPhi to ensure the corrected phi has the
84  // same divergence
85  volScalarField divrhoU
86  (
87  "divrhoU",
89  );
90 
91  if (LTS)
92  {
93  #include "setRDeltaT.H"
94  }
95  else
96  {
97  #include "compressibleCourantNo.H"
98  #include "setDeltaT.H"
99  }
100 
101  runTime++;
102 
103  Info<< "Time = " << runTime.timeName() << nl << endl;
104 
105  // Store momentum to set rhoUf for introduced faces.
106  volVectorField rhoU("rhoU", rho*U);
107 
108  // Do any mesh changes
109  mesh.update();
110 
111  if (mesh.changing() && correctPhi)
112  {
113  // Calculate absolute flux from the mapped surface velocity
114  phi = mesh.Sf() & rhoUf;
115 
116  #include "correctPhi.H"
117 
118  // Make the fluxes relative to the mesh-motion
120  }
121  }
122 
123  if (mesh.changing() && checkMeshCourantNo)
124  {
125  #include "meshCourantNo.H"
126  }
127 
128  #include "rhoEqn.H"
129  Info<< "rhoEqn max/min : " << max(rho).value()
130  << " " << min(rho).value() << endl;
131 
132  // --- Pressure-velocity PIMPLE corrector loop
133  while (pimple.loop())
134  {
135  #include "UEqn.H"
136  #include "EEqn.H"
137 
138  // --- Pressure corrector loop
139  while (pimple.correct())
140  {
141  #include "pEqn.H"
142  }
143 
144  if (pimple.turbCorr())
145  {
146  turbulence->correct();
147  }
148  }
149 
150  runTime.write();
151 
152  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
153  << " ClockTime = " << runTime.elapsedClockTime() << " s"
154  << nl << endl;
155  }
156 
157  Info<< "End\n" << endl;
158 
159  return 0;
160 }
161 
162 
163 // ************************************************************************* //
surfaceScalarField & phi
U
Definition: pEqn.H:83
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
const dictionary & pimple
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
static const char nl
Definition: Ostream.H:262
Bound the given scalar field if it has gone unbounded.
bool LTS
Definition: createRDeltaT.H:1
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
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Calculates and outputs the mean and maximum Courant Numbers.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
rhoUf
Definition: pEqn.H:91
messageStream Info
checkMeshCourantNo
Definition: readControls.H:5
correctPhi
Definition: readControls.H:3
Execute application functionObjects to post-process existing results.
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:75
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
Creates and initialises the velocity velocity field rhoUf.