OpenFOAM: API Guide
v2112
The open source CFD toolbox
sprayDyMFoam.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) 2015-2017 OpenFOAM Foundation
9
-------------------------------------------------------------------------------
10
License
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
26
Application
27
sprayDyMFoam
28
29
Group
30
grpLagrangianSolvers grpMovingMeshSolvers
31
32
Description
33
Transient solver for compressible, turbulent flow with a spray particle
34
cloud, with optional mesh motion and mesh topology changes.
35
36
\*---------------------------------------------------------------------------*/
37
38
#include "
fvCFD.H
"
39
#include "
dynamicFvMesh.H
"
40
#include "
turbulenceModel.H
"
41
#include "
basicSprayCloud.H
"
42
#include "
psiReactionThermo.H
"
43
#include "
CombustionModel.H
"
44
#include "
radiationModel.H
"
45
#include "
SLGThermo.H
"
46
#include "
pimpleControl.H
"
47
#include "
CorrectPhi.H
"
48
#include "
fvOptions.H
"
49
50
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51
52
int
main(
int
argc,
char
*argv[])
53
{
54
argList::addNote
55
(
56
"Transient solver for compressible, turbulent flow"
57
" with a spray particle cloud.\n"
58
"With optional mesh motion and mesh topology changes.\n"
59
);
60
61
#include "
postProcess.H
"
62
63
#include "
setRootCaseLists.H
"
64
#include "
createTime.H
"
65
#include "
createDynamicFvMesh.H
"
66
#include "createDyMControls.H"
67
#include "createFields.H"
68
#include "createFieldRefs.H"
69
#include "
createRhoUf.H
"
70
#include "compressibleCourantNo.H"
71
#include "setInitialDeltaT.H"
72
#include "initContinuityErrs.H"
73
74
turbulence
->validate();
75
76
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77
78
Info
<<
"\nStarting time loop\n"
<<
endl
;
79
80
while
(
runTime
.run())
81
{
82
#include "
readDyMControls.H
"
83
84
{
85
// Store divrhoU from the previous time-step/mesh for the correctPhi
86
volScalarField
divrhoU
87
(
88
"divrhoU"
,
89
fvc::div(fvc::absolute(
phi
,
rho
,
U
))
90
);
91
92
#include "compressibleCourantNo.H"
93
#include "setDeltaT.H"
94
95
++
runTime
;
96
97
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
98
99
// Store momentum to set rhoUf for introduced faces.
100
volVectorField
rhoU(
"rhoU"
,
rho
*
U
);
101
102
// Store the particle positions
103
parcels.storeGlobalPositions();
104
105
// Do any mesh changes
106
mesh
.update();
107
108
if
(
mesh
.changing())
109
{
110
MRF
.update();
111
112
if
(
correctPhi
)
113
{
114
// Calculate absolute flux from the mapped surface velocity
115
phi
=
mesh
.Sf() &
rhoUf
;
116
117
#include "correctPhi.H"
118
119
// Make the fluxes relative to the mesh-motion
120
fvc::makeRelative(
phi
,
rho
,
U
);
121
}
122
123
if
(
checkMeshCourantNo
)
124
{
125
#include "
meshCourantNo.H
"
126
}
127
}
128
}
129
130
parcels.evolve();
131
132
#include "rhoEqn.H"
133
134
// --- Pressure-velocity PIMPLE corrector loop
135
while
(
pimple
.loop())
136
{
137
#include "UEqn.H"
138
#include "YEqn.H"
139
#include "EEqn.H"
140
141
// --- Pressure corrector loop
142
while
(
pimple
.correct())
143
{
144
#include "pEqn.H"
145
}
146
147
if
(
pimple
.turbCorr())
148
{
149
turbulence
->correct();
150
}
151
}
152
153
rho
=
thermo
.rho();
154
155
if
(
runTime
.write())
156
{
157
combustion
->Qdot()().
write
();
158
}
159
160
runTime
.printExecutionTime(Info);
161
}
162
163
Info
<<
"End\n"
<<
endl
;
164
165
return
0;
166
}
167
168
169
// ************************************************************************* //
CombustionModel.H
CorrectPhi.H
SLGThermo.H
basicSprayCloud.H
phi
surfaceScalarField & phi
Definition:
setRegionFluidFields.H:8
MRF
IOMRFZoneList & MRF
Definition:
setRegionFluidFields.H:22
pimple
pimpleControl & pimple
Definition:
setRegionFluidFields.H:56
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
U
Definition:
pEqn.H:72
createDynamicFvMesh.H
mesh
dynamicFvMesh & mesh
Definition:
createDynamicFvMesh.H:6
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
rhoUf
autoPtr< surfaceVectorField > rhoUf
Definition:
createRhoUfIfPresent.H:33
createRhoUf.H
Creates and initialises the velocity velocity field rhoUf.
createTime.H
dynamicFvMesh.H
turbulence
compressible::turbulenceModel & turbulence
Definition:
setRegionFluidFields.H:30
fvCFD.H
fvOptions.H
meshCourantNo.H
Calculates and outputs the mean and maximum Courant Numbers.
Foam::volVectorField
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition:
volFieldsFwd.H:83
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition:
volFieldsFwd.H:82
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition:
Ostream.H:372
Foam::nl
constexpr char nl
The newline '\n' character (0x0a)
Definition:
Ostream.H:53
write
runTime write()
pimpleControl.H
postProcess.H
Execute application functionObjects to post-process existing results.
psiReactionThermo.H
radiationModel.H
readDyMControls.H
checkMeshCourantNo
checkMeshCourantNo
Definition:
readDyMControls.H:9
correctPhi
correctPhi
Definition:
readDyMControls.H:3
rho
rho
Definition:
readInitialConditions.H:88
setRootCaseLists.H
combustion
Info<< "Creating combustion model\n"<< endl;autoPtr< CombustionModel< psiReactionThermo > > combustion(CombustionModel< psiReactionThermo >::New(thermo, turbulence()))
turbulenceModel.H
applications
solvers
lagrangian
sprayFoam
sprayDyMFoam
sprayDyMFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.