potentialFreeSurfaceDyMFoam.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) 2014-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Application
27 potentialFreeSurfaceDyMFoam
28
29Group
30 grpMultiphaseSolvers grpMovingMeshSolvers
31
32Description
33 Incompressible Navier-Stokes solver with inclusion of a wave height field
34 to enable single-phase free-surface approximations, with optional mesh
35 motion and mesh topology changes.
36
37 Wave height field, zeta, used by pressure boundary conditions.
38
39 Optional mesh motion and mesh topology changes including adaptive
40 re-meshing.
41
42 Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
43
44\*---------------------------------------------------------------------------*/
45
46#include "fvCFD.H"
47#include "dynamicFvMesh.H"
50#include "pimpleControl.H"
51#include "fvOptions.H"
52#include "CorrectPhi.H"
53
54// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55
56int main(int argc, char *argv[])
57{
58 argList::addNote
59 (
60 "Incompressible Navier-Stokes solver with inclusion of a wave height"
61 " field to enable single-phase free-surface approximations.\n"
62 "With optional mesh motion and mesh topology changes."
63 );
64
65 #include "postProcess.H"
66
67 #include "setRootCaseLists.H"
68 #include "createTime.H"
69 #include "createDynamicFvMesh.H"
70 #include "initContinuityErrs.H"
71 #include "createDyMControls.H"
72 #include "createFields.H"
73
75 (
76 IOobject
77 (
78 "rAU",
79 runTime.timeName(),
80 mesh,
81 IOobject::READ_IF_PRESENT,
82 IOobject::AUTO_WRITE
83 ),
84 mesh,
85 dimensionedScalar("rAUf", dimTime, 1.0)
86 );
87
88 #include "correctPhi.H"
89 #include "createUf.H"
90
91 turbulence->validate();
92
93 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
94
95 Info<< "\nStarting time loop\n" << endl;
96
97 while (runTime.run())
98 {
99 #include "readDyMControls.H"
100 #include "CourantNo.H"
101 #include "setDeltaT.H"
102
103 ++runTime;
104
105 Info<< "Time = " << runTime.timeName() << nl << endl;
106
107 // --- Pressure-velocity PIMPLE corrector loop
108 while (pimple.loop())
109 {
110 if (pimple.firstIter() || moveMeshOuterCorrectors)
111 {
112 scalar timeBeforeMeshUpdate = runTime.elapsedCpuTime();
113
114 mesh.update();
115
116 if (mesh.changing())
117 {
118 Info<< "Execution time for mesh.update() = "
119 << runTime.elapsedCpuTime() - timeBeforeMeshUpdate
120 << " s" << endl;
121
122 MRF.update();
123
124 if (correctPhi)
125 {
126 // Calculate absolute flux
127 // from the mapped surface velocity
128 phi = mesh.Sf() & Uf;
129
130 #include "correctPhi.H"
131
132 // Make the flux relative to the mesh motion
133 fvc::makeRelative(phi, U);
134 }
135
137 {
138 #include "meshCourantNo.H"
139 }
140 }
141 }
142
143 #include "UEqn.H"
144
145 // --- Pressure corrector loop
146 while (pimple.correct())
147 {
148 #include "pEqn.H"
149 }
150
151 if (pimple.turbCorr())
152 {
153 laminarTransport.correct();
154 turbulence->correct();
155 }
156 }
157
158 runTime.write();
159
160 runTime.printExecutionTime(Info);
161 }
162
163 Info<< "End\n" << endl;
164
165 return 0;
166}
167
168
169// ************************************************************************* //
surfaceScalarField & phi
IOMRFZoneList & MRF
pimpleControl & pimple
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
engineTime & runTime
autoPtr< surfaceVectorField > Uf
Creates and initialises the velocity velocity field Uf.
compressible::turbulenceModel & turbulence
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
Calculates and outputs the mean and maximum Courant Numbers.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
messageStream Info
Information stream (stdout output on master, null elsewhere)
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
Execute application functionObjects to post-process existing results.
checkMeshCourantNo
correctPhi
moveMeshOuterCorrectors
singlePhaseTransportModel laminarTransport(U, phi)