overPimpleDyMFoam.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2016-2018 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Application
28 overPimpleDyMFoam
29
30Group
31 grpIncompressibleSolvers grpMovingMeshSolvers
32
33Description
34 Transient solver for incompressible flow of Newtonian fluids
35 on a moving mesh using the PIMPLE (merged PISO-SIMPLE) algorithm.
36
37 Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
38
39\*---------------------------------------------------------------------------*/
40
41#include "fvCFD.H"
42#include "dynamicFvMesh.H"
45#include "pimpleControl.H"
46#include "fvOptions.H"
47
50#include "localMin.H"
52#include "transform.H"
53#include "fvMeshSubset.H"
54#include "oversetAdjustPhi.H"
55
56// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57
58int main(int argc, char *argv[])
59{
60 argList::addNote
61 (
62 "Transient solver for incompressible, turbulent flow"
63 " on a moving mesh."
64 );
65
66 #include "postProcess.H"
67
68 #include "setRootCaseLists.H"
69 #include "createTime.H"
70 #include "createDynamicFvMesh.H"
71 #include "initContinuityErrs.H"
72
73 pimpleControl pimple(mesh);
74
75 #include "createFields.H"
76 #include "createUf.H"
77 #include "createMRF.H"
78 #include "createFvOptions.H"
79 #include "createControls.H"
80 #include "CourantNo.H"
81 #include "setInitialDeltaT.H"
82
83 turbulence->validate();
84
85 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
86
87 Info<< "\nStarting time loop\n" << endl;
88
89 while (runTime.run())
90 {
91 #include "readControls.H"
92 #include "CourantNo.H"
93
94 #include "setDeltaT.H"
95
96 ++runTime;
97
98 Info<< "Time = " << runTime.timeName() << nl << endl;
99
100 bool changed = mesh.update();
101
102 if (changed)
103 {
104 #include "setCellMask.H"
105 #include "setInterpolatedCells.H"
106
107 surfaceScalarField faceMaskOld
108 (
109 localMin<scalar>(mesh).interpolate(cellMask.oldTime())
110 );
111
112 // Zero Uf on old faceMask (H-I)
113 Uf *= faceMaskOld;
114 // Update Uf and phi on new C-I faces
115 Uf += (1-faceMaskOld)*fvc::interpolate(U);
116 phi = mesh.Sf() & Uf;
117
118 // Zero phi on current H-I
120 (
121 localMin<scalar>(mesh).interpolate(cellMask)
122 );
123 phi *= faceMask;
124 }
125
126
127 if (mesh.changing() && correctPhi)
128 {
129 // Calculate absolute flux from the mapped surface velocity
130 #include "correctPhi.H"
131 }
132
133 // Make the flux relative to the mesh motion
134 fvc::makeRelative(phi, U);
135
136 if (mesh.changing() && checkMeshCourantNo)
137 {
138 #include "meshCourantNo.H"
139 }
140
141 // --- Pressure-velocity PIMPLE corrector loop
142 while (pimple.loop())
143 {
144 #include "UEqn.H"
145
146 // --- Pressure corrector loop
147 while (pimple.correct())
148 {
149 #include "pEqn.H"
150 }
151
152 if (pimple.turbCorr())
153 {
154 laminarTransport.correct();
155 turbulence->correct();
156 }
157 }
158
159 runTime.write();
160
161 runTime.printExecutionTime(Info);
162 }
163
164 Info<< "End\n" << endl;
165
166 return 0;
167}
168
169
170// ************************************************************************* //
surfaceScalarField & phi
pimpleControl & pimple
U
Definition: pEqn.H:72
surfaceScalarField faceMask(localMin< scalar >(mesh).interpolate(cellMask))
mesh interpolate(rAU)
dynamicFvMesh & mesh
engineTime & runTime
autoPtr< surfaceVectorField > Uf
Creates and initialises the velocity velocity field Uf.
compressible::turbulenceModel & turbulence
Calculates and outputs the mean and maximum Courant Numbers.
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Adjust the balance of fluxes on the faces between interpolated and calculated to obey continuity.
Execute application functionObjects to post-process existing results.
checkMeshCourantNo
correctPhi
Sets blocked cells mask field.
Sets blocked cells mask field.
singlePhaseTransportModel laminarTransport(U, phi)
3D tensor transformation operations.