All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
pimpleFoam.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) 2011-2020 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  pimpleFoam
26 
27 Description
28  Transient solver for incompressible, turbulent flow of Newtonian fluids,
29  with optional mesh motion and mesh topology changes.
30 
31  Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #include "fvCFD.H"
36 #include "dynamicFvMesh.H"
39 #include "pimpleControl.H"
40 #include "CorrectPhi.H"
41 #include "fvOptions.H"
42 #include "localEulerDdtScheme.H"
43 #include "fvcSmooth.H"
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 int main(int argc, char *argv[])
48 {
49  #include "postProcess.H"
50 
51  #include "setRootCaseLists.H"
52  #include "createTime.H"
53  #include "createDynamicFvMesh.H"
54  #include "initContinuityErrs.H"
55  #include "createDyMControls.H"
56  #include "createFields.H"
57  #include "createUfIfPresent.H"
58 
59  turbulence->validate();
60 
61  if (!LTS)
62  {
63  #include "CourantNo.H"
64  #include "setInitialDeltaT.H"
65  }
66 
67  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
68 
69  Info<< "\nStarting time loop\n" << endl;
70 
71  while (pimple.run(runTime))
72  {
73  #include "readDyMControls.H"
74 
75  if (LTS)
76  {
77  #include "setRDeltaT.H"
78  }
79  else
80  {
81  #include "CourantNo.H"
82  #include "setDeltaT.H"
83  }
84 
85  runTime++;
86 
87  Info<< "Time = " << runTime.timeName() << nl << endl;
88 
89  // --- Pressure-velocity PIMPLE corrector loop
90  while (pimple.loop())
91  {
92  if (pimple.firstPimpleIter() || moveMeshOuterCorrectors)
93  {
94  mesh.update();
95 
96  if (mesh.changing())
97  {
98  MRF.update();
99 
100  if (correctPhi)
101  {
102  // Calculate absolute flux
103  // from the mapped surface velocity
104  phi = mesh.Sf() & Uf();
105 
106  #include "correctPhi.H"
107 
108  // Make the flux relative to the mesh motion
110  }
111 
112  if (checkMeshCourantNo)
113  {
114  #include "meshCourantNo.H"
115  }
116  }
117  }
118 
119  #include "UEqn.H"
120 
121  // --- Pressure corrector loop
122  while (pimple.correct())
123  {
124  #include "pEqn.H"
125  }
126 
127  if (pimple.turbCorr())
128  {
129  laminarTransport.correct();
130  turbulence->correct();
131  }
132  }
133 
134  runTime.write();
135 
136  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
137  << " ClockTime = " << runTime.elapsedClockTime() << " s"
138  << nl << endl;
139  }
140 
141  Info<< "End\n" << endl;
142 
143  return 0;
144 }
145 
146 
147 // ************************************************************************* //
pimpleNoLoopControl & pimple
IOMRFZoneList & MRF
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
correctPhi
checkMeshCourantNo
phi
Definition: pEqn.H:104
dynamicFvMesh & mesh
autoPtr< surfaceVectorField > Uf
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(dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar(dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::momentumTransportModel > turbulence(compressible::momentumTransportModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
static const char nl
Definition: Ostream.H:260
bool LTS
Definition: createRDeltaT.H:1
moveMeshOuterCorrectors
Calculates and outputs the mean and maximum Courant Numbers.
U
Definition: pEqn.H:72
messageStream Info
Execute application functionObjects to post-process existing results.
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:75
singlePhaseTransportModel laminarTransport(U, phi)
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
Creates and initialises the velocity field Uf if required.