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 | Copyright (C) 2014-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  potentialFreeSurfaceDyMFoam
26 
27 Description
28  Incompressible Navier-Stokes solver with inclusion of a wave height field
29  to enable single-phase free-surface approximations.
30 
31  Wave height field, zeta, used by pressure boundary conditions.
32 
33  Optional mesh motion and mesh topology changes including adaptive
34  re-meshing.
35 
36  Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
37 
38 \*---------------------------------------------------------------------------*/
39 
40 #include "fvCFD.H"
41 #include "dynamicFvMesh.H"
44 #include "pimpleControl.H"
45 #include "fvIOoptionList.H"
46 #include "CorrectPhi.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 int main(int argc, char *argv[])
52 {
53  #include "setRootCase.H"
54  #include "createTime.H"
55  #include "createDynamicFvMesh.H"
56  #include "initContinuityErrs.H"
57 
58  pimpleControl pimple(mesh);
59 
60  #include "createControls.H"
61  #include "createFields.H"
62  #include "createMRF.H"
63  #include "createFvOptions.H"
64 
66  (
67  IOobject
68  (
69  "rAU",
70  runTime.timeName(),
71  mesh,
72  IOobject::READ_IF_PRESENT,
73  IOobject::AUTO_WRITE
74  ),
75  mesh,
76  dimensionedScalar("rAUf", dimTime, 1.0)
77  );
78 
79  #include "correctPhi.H"
80  #include "createUf.H"
81 
82  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
83 
84  Info<< "\nStarting time loop\n" << endl;
85 
86  while (runTime.run())
87  {
88  #include "readControls.H"
89  #include "CourantNo.H"
90  #include "setDeltaT.H"
91 
92  runTime++;
93 
94  Info<< "Time = " << runTime.timeName() << nl << endl;
95 
96  // --- Pressure-velocity PIMPLE corrector loop
97  while (pimple.loop())
98  {
99  if (pimple.firstIter() || moveMeshOuterCorrectors)
100  {
101  scalar timeBeforeMeshUpdate = runTime.elapsedCpuTime();
102 
103  mesh.update();
104 
105  if (mesh.changing())
106  {
107  Info<< "Execution time for mesh.update() = "
108  << runTime.elapsedCpuTime() - timeBeforeMeshUpdate
109  << " s" << endl;
110  }
111 
112  if (mesh.changing() && correctPhi)
113  {
114  // Calculate absolute flux from the mapped surface velocity
115  phi = mesh.Sf() & Uf;
116 
117  #include "correctPhi.H"
118 
119  // Make the flux relative to the mesh motion
121  }
122 
123  if (mesh.changing() && checkMeshCourantNo)
124  {
125  #include "meshCourantNo.H"
126  }
127  }
128 
129  #include "UEqn.H"
130 
131  // --- Pressure corrector loop
132  while (pimple.correct())
133  {
134  #include "pEqn.H"
135  }
136 
137  if (pimple.turbCorr())
138  {
139  laminarTransport.correct();
140  turbulence->correct();
141  }
142  }
143 
144  runTime.write();
145 
146  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
147  << " ClockTime = " << runTime.elapsedClockTime() << " s"
148  << nl << endl;
149  }
150 
151  Info<< "End\n" << endl;
152 
153  return 0;
154 }
155 
156 
157 // ************************************************************************* //
Uf
Definition: pEqn.H:78
singlePhaseTransportModel laminarTransport(U, phi)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:74
surfaceScalarField & phi
volScalarField rAU(1.0/UEqn.A())
correctPhi
Definition: readControls.H:3
moveMeshOuterCorrectors
Definition: readControls.H:8
Calculates and outputs the mean and maximum Courant Numbers.
int main(int argc, char *argv[])
Definition: postCalc.C:54
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
U
Definition: pEqn.H:82
Creates and initialises the velocity velocity field Uf.