sonicLiquidFoam.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-------------------------------------------------------------------------------
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
26Group
27 grpCompressibleSolvers
28
29Application
30 sonicLiquidFoam
31
32Description
33 Transient solver for trans-sonic/supersonic, laminar flow of a
34 compressible liquid.
35
36\*---------------------------------------------------------------------------*/
37
38#include "fvCFD.H"
39#include "pimpleControl.H"
40
41// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
43int main(int argc, char *argv[])
44{
45 argList::addNote
46 (
47 "Transient solver for trans-sonic/supersonic, laminar flow"
48 " of a compressible liquid."
49 );
50
51 #include "postProcess.H"
52
53 #include "addCheckCaseOptions.H"
54 #include "setRootCaseLists.H"
55 #include "createTime.H"
56 #include "createMesh.H"
57 #include "createControl.H"
58 #include "createFields.H"
59 #include "initContinuityErrs.H"
60
61 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62
63 Info<< "\nStarting time loop\n" << endl;
64
65 while (runTime.loop())
66 {
67 Info<< "Time = " << runTime.timeName() << nl << endl;
68
69 #include "compressibleCourantNo.H"
70
71 solve(fvm::ddt(rho) + fvc::div(phi));
72
73 // --- Pressure-velocity PIMPLE corrector loop
74 while (pimple.loop())
75 {
77 (
78 fvm::ddt(rho, U)
79 + fvm::div(phi, U)
80 - fvm::laplacian(mu, U)
81 );
82
83 solve(UEqn == -fvc::grad(p));
84
85 // --- Pressure corrector loop
86 while (pimple.correct())
87 {
88 volScalarField rAU("rAU", 1.0/UEqn.A());
90 (
91 "rhorAUf",
92 fvc::interpolate(rho*rAU)
93 );
94
95 U = rAU*UEqn.H();
96
98 (
99 "phid",
100 psi
101 *(
102 fvc::flux(U)
103 + rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
104 )
105 );
106
107 phi = (rhoO/psi)*phid;
108
109 fvScalarMatrix pEqn
110 (
111 fvm::ddt(psi, p)
112 + fvc::div(phi)
113 + fvm::div(phid, p)
114 - fvm::laplacian(rhorAUf, p)
115 );
116
117 pEqn.solve();
118
119 phi += pEqn.flux();
120
121 solve(fvm::ddt(rho) + fvc::div(phi));
122 #include "compressibleContinuityErrs.H"
123
124 U -= rAU*fvc::grad(p);
125 U.correctBoundaryConditions();
126 }
127 }
128
129 rho = rhoO + psi*p;
130
131 runTime.write();
132
133 runTime.printExecutionTime(Info);
134 }
135
136 Info<< "End\n" << endl;
137
138 return 0;
139}
140
141
142// ************************************************************************* //
Required Classes.
surfaceScalarField & phi
pimpleControl & pimple
U
Definition: pEqn.H:72
volScalarField & p
const volScalarField & psi
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
fvVectorMatrix & UEqn
Definition: UEqn.H:13
const volScalarField & mu
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
engineTime & runTime
Required Variables.
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
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
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:46
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Execute application functionObjects to post-process existing results.
CEqn solve()