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