Skip to content

Commit e8888c1

Browse files
committed
legacy::cavitationFoam: recover suspended cavitatingFoam
See `cavitatingFoam` solver in parent of OpenFOAM/OpenFOAM-dev@776ecc9
1 parent 8f82089 commit e8888c1

File tree

19 files changed

+1062
-4
lines changed

19 files changed

+1062
-4
lines changed
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,58 @@
1+
/*---------------------------------------------------------------------------*\
2+
========= |
3+
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4+
\\ / O peration | Website: https://openfoam.org
5+
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
6+
\\/ M anipulation |
7+
-------------------------------------------------------------------------------
8+
License
9+
This file is originating from OpenFOAM and modified by the authors
10+
described below.
11+
12+
OpenFOAM is free software: you can redistribute it and/or modify it
13+
under the terms of the GNU General Public License as published by
14+
the Free Software Foundation, either version 3 of the License, or
15+
(at your option) any later version.
16+
17+
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
18+
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
19+
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
20+
for more details.
21+
22+
You should have received a copy of the GNU General Public License
23+
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
24+
25+
Global
26+
CourantNo
27+
28+
Description
29+
Calculates and outputs the mean and maximum Courant Numbers.
30+
31+
\*---------------------------------------------------------------------------*/
32+
33+
scalar CoNum = 0;
34+
scalar acousticCoNum = 0;
35+
36+
{
37+
const scalarField sumPhi(fvc::surfaceSum(mag(phi))().primitiveField());
38+
39+
CoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue();
40+
41+
const scalar meanCoNum =
42+
0.5*(gSum(sumPhi)/gSum(mesh.V().field()))*runTime.deltaTValue();
43+
44+
acousticCoNum = 0.5*gMax
45+
(
46+
fvc::surfaceSum
47+
(
48+
fvc::interpolate(scalar(1)/sqrt(mixture.psi()))*mesh.magSf()
49+
)().primitiveField()/mesh.V().field()
50+
)*runTime.deltaTValue();
51+
52+
Info<< "phi Courant Number mean: " << meanCoNum
53+
<< " max: " << CoNum
54+
<< " acoustic max: " << acousticCoNum
55+
<< endl;
56+
}
57+
58+
// ************************************************************************* //
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,4 @@
1+
cavitatingTwoPhaseMixture/cavitatingTwoPhaseMixture.C
2+
cavitationFoam.C
3+
4+
EXE = $(FOAM_USER_APPBIN)/cavitationFoam
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,21 @@
1+
EXE_INC = \
2+
-IcavitatingTwoPhaseMixture \
3+
-I$(LIB_SRC)/physicalProperties/lnInclude \
4+
-I$(LIB_SRC)/twoPhaseModels/VoF \
5+
-I$(LIB_SRC)/twoPhaseModels/twoPhaseMixture/lnInclude \
6+
-I$(LIB_SRC)/thermophysicalModels/barotropicCompressibilityModel/lnInclude \
7+
-I$(LIB_SRC)/MomentumTransportModels/momentumTransportModels/lnInclude \
8+
-I$(LIB_SRC)/MomentumTransportModels/incompressible/lnInclude \
9+
-I$(LIB_SRC)/finiteVolume/lnInclude \
10+
-I$(LIB_SRC)/meshTools/lnInclude
11+
12+
EXE_LIBS = \
13+
-lphysicalProperties \
14+
-ltwoPhaseMixture \
15+
-lbarotropicCompressibilityModel \
16+
-lmomentumTransportModels \
17+
-lincompressibleMomentumTransportModels \
18+
-lfiniteVolume \
19+
-lmeshTools \
20+
-lfvModels \
21+
-lfvConstraints
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
fvVectorMatrix UEqn
2+
(
3+
fvm::ddt(rho, U)
4+
+ fvm::div(rhoPhi, U)
5+
+ turbulence->divDevTau(rho, U)
6+
);
7+
8+
UEqn.relax();
9+
10+
if (pimple.momentumPredictor())
11+
{
12+
solve(UEqn == -fvc::grad(p));
13+
}
14+
15+
Info<< "max(U) " << max(mag(U)).value() << endl;
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,251 @@
1+
/*---------------------------------------------------------------------------*\
2+
========= |
3+
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4+
\\ / O peration | Website: https://openfoam.org
5+
\\ / A nd | Copyright (C) 2022 OpenFOAM Foundation
6+
\\/ M anipulation |
7+
-------------------------------------------------------------------------------
8+
License
9+
This file is originating from OpenFOAM but modified by authors described
10+
below.
11+
12+
OpenFOAM is free software: you can redistribute it and/or modify it
13+
under the terms of the GNU General Public License as published by
14+
the Free Software Foundation, either version 3 of the License, or
15+
(at your option) any later version.
16+
17+
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
18+
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
19+
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
20+
for more details.
21+
22+
You should have received a copy of the GNU General Public License
23+
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
24+
25+
Authors
26+
Henry Weller, CFD-Direct, 2023.
27+
28+
\*---------------------------------------------------------------------------*/
29+
30+
#include "cavitatingTwoPhaseMixture.H"
31+
#include "barotropicCompressibilityModel.H"
32+
33+
34+
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35+
36+
namespace Foam
37+
{
38+
defineTypeNameAndDebug(cavitatingTwoPhaseMixture, 0);
39+
}
40+
41+
42+
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43+
44+
Foam::cavitatingTwoPhaseMixture::cavitatingTwoPhaseMixture
45+
(
46+
const fvMesh& mesh
47+
)
48+
:
49+
twoPhaseMixture(mesh),
50+
51+
alphav_(alpha1()),
52+
alphal_(alpha2()),
53+
54+
nuModelv_(viscosityModel::New(mesh, phase1Name())),
55+
nuModell_(viscosityModel::New(mesh, phase2Name())),
56+
57+
rhov_("rho", dimDensity, nuModelv_()),
58+
rhol_("rho", dimDensity, nuModell_()),
59+
60+
nu_
61+
(
62+
IOobject
63+
(
64+
"nu",
65+
mesh.time().name(),
66+
mesh
67+
),
68+
mesh,
69+
dimensionedScalar(dimViscosity, 0),
70+
calculatedFvPatchScalarField::typeName
71+
),
72+
73+
thermodynamicProperties_
74+
(
75+
IOobject
76+
(
77+
"thermodynamicProperties",
78+
mesh.time().constant(),
79+
mesh,
80+
IOobject::MUST_READ_IF_MODIFIED,
81+
IOobject::NO_WRITE
82+
)
83+
),
84+
85+
psil_
86+
(
87+
"psil",
88+
dimCompressibility,
89+
thermodynamicProperties_
90+
),
91+
92+
rholSat_
93+
(
94+
"rholSat",
95+
dimDensity,
96+
thermodynamicProperties_
97+
),
98+
99+
psiv_
100+
(
101+
"psiv",
102+
dimCompressibility,
103+
thermodynamicProperties_
104+
),
105+
106+
pSat_
107+
(
108+
"pSat",
109+
dimPressure,
110+
thermodynamicProperties_
111+
),
112+
113+
rhovSat_("rhovSat", psiv_*pSat_),
114+
115+
rhol0_("rhol0", rholSat_ - pSat_*psil_),
116+
117+
rhoMin_
118+
(
119+
"rhoMin",
120+
dimDensity,
121+
thermodynamicProperties_
122+
),
123+
124+
p_
125+
(
126+
IOobject
127+
(
128+
"p",
129+
mesh.time().name(),
130+
mesh,
131+
IOobject::MUST_READ,
132+
IOobject::AUTO_WRITE
133+
),
134+
mesh
135+
),
136+
137+
psiModel_
138+
(
139+
barotropicCompressibilityModel::New
140+
(
141+
thermodynamicProperties_,
142+
alphav_
143+
)
144+
),
145+
146+
psi_(psiModel_->psi()),
147+
148+
rho_
149+
(
150+
IOobject
151+
(
152+
"rho",
153+
mesh.time().name(),
154+
mesh,
155+
IOobject::MUST_READ,
156+
IOobject::AUTO_WRITE
157+
),
158+
mesh
159+
)
160+
{
161+
alphav_.oldTime();
162+
163+
rho_ = max
164+
(
165+
psi_*p_
166+
+ alphal_*rhol0_
167+
+ ((alphav_*psiv_ + alphal_*psil_) - psi_)*pSat_,
168+
rhoMin_
169+
);
170+
171+
correct();
172+
}
173+
174+
175+
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
176+
177+
Foam::cavitatingTwoPhaseMixture::~cavitatingTwoPhaseMixture()
178+
{}
179+
180+
181+
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
182+
183+
bool Foam::cavitatingTwoPhaseMixture::read()
184+
{
185+
if (twoPhaseMixture::read())
186+
{
187+
nuModelv_->lookup("rho") >> rhov_;
188+
nuModell_->lookup("rho") >> rhol_;
189+
190+
return true;
191+
}
192+
else
193+
{
194+
return false;
195+
}
196+
}
197+
198+
199+
void Foam::cavitatingTwoPhaseMixture::correctPressure()
200+
{
201+
p_ =
202+
(
203+
rho_
204+
- alphal_*rhol0_
205+
- ((alphav_*psiv_ + alphal_*psil_) - psi_)*pSat_
206+
)/psi_;
207+
208+
p_.correctBoundaryConditions();
209+
}
210+
211+
212+
void Foam::cavitatingTwoPhaseMixture::correct()
213+
{
214+
rho_ == max(rho_, rhoMin_);
215+
216+
alphav_ =
217+
max
218+
(
219+
min
220+
(
221+
(rho_ - rholSat_)/(rhovSat_ - rholSat_),
222+
scalar(1)
223+
),
224+
scalar(0)
225+
);
226+
alphal_ = 1.0 - alphav_;
227+
228+
Info<< "max-min alphav: " << max(alphav_).value()
229+
<< " " << min(alphav_).value() << endl;
230+
231+
psiModel_->correct();
232+
233+
nuModelv_->correct();
234+
nuModell_->correct();
235+
236+
const volScalarField limitedAlphav
237+
(
238+
"limitedAlphav",
239+
min(max(alphav_, scalar(0)), scalar(1))
240+
);
241+
242+
// Mixture kinematic viscosity calculated from mixture dynamic viscosity
243+
nu_ =
244+
(
245+
limitedAlphav*rhov_*nuModelv_->nu()
246+
+ (scalar(1) - limitedAlphav)*rhol_*nuModell_->nu()
247+
)/(limitedAlphav*rhov_ + (scalar(1) - limitedAlphav)*rhol_);
248+
}
249+
250+
251+
// ************************************************************************* //

0 commit comments

Comments
 (0)