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-2018 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 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 int main(int argc, char *argv[])
46 {
47  #include "postProcess.H"
48 
49  #include "setRootCaseLists.H"
50  #include "createTime.H"
51  #include "createDynamicFvMesh.H"
52  #include "initContinuityErrs.H"
53  #include "createDyMControls.H"
54  #include "createFields.H"
55  #include "createUfIfPresent.H"
56  #include "CourantNo.H"
57  #include "setInitialDeltaT.H"
58 
59  turbulence->validate();
60 
61  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62 
63  Info<< "\nStarting time loop\n" << endl;
64 
65  while (runTime.run())
66  {
67  #include "readDyMControls.H"
68  #include "CourantNo.H"
69  #include "setDeltaT.H"
70 
71  runTime++;
72 
73  Info<< "Time = " << runTime.timeName() << nl << endl;
74 
75  // --- Pressure-velocity PIMPLE corrector loop
76  while (pimple.loop())
77  {
78  if (pimple.firstIter() || moveMeshOuterCorrectors)
79  {
80  mesh.update();
81 
82  if (mesh.changing())
83  {
84  MRF.update();
85 
86  if (correctPhi)
87  {
88  // Calculate absolute flux
89  // from the mapped surface velocity
90  phi = mesh.Sf() & Uf();
91 
92  #include "correctPhi.H"
93 
94  // Make the flux relative to the mesh motion
96  }
97 
99  {
100  #include "meshCourantNo.H"
101  }
102  }
103  }
104 
105  #include "UEqn.H"
106 
107  // --- Pressure corrector loop
108  while (pimple.correct())
109  {
110  #include "pEqn.H"
111  }
112 
113  if (pimple.turbCorr())
114  {
115  laminarTransport.correct();
116  turbulence->correct();
117  }
118  }
119 
120  runTime.write();
121 
122  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
123  << " ClockTime = " << runTime.elapsedClockTime() << " s"
124  << nl << endl;
125  }
126 
127  Info<< "End\n" << endl;
128 
129  return 0;
130 }
131 
132 
133 // ************************************************************************* //
pimpleNoLoopControl & pimple
surfaceScalarField & phi
IOMRFZoneList & MRF
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
correctPhi
checkMeshCourantNo
singlePhaseTransportModel laminarTransport(U, phi)
dynamicFvMesh & mesh
autoPtr< surfaceVectorField > Uf
static const char nl
Definition: Ostream.H:265
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
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
Creates and initialises the velocity field Uf if required.