OpenFOAM: API Guide
v2112
The open source CFD toolbox
nonNewtonianIcoFoam.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
-------------------------------------------------------------------------------
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
nonNewtonianIcoFoam
28
29
Group
30
grpIncompressibleSolvers
31
32
Description
33
Transient solver for incompressible, laminar flow of non-Newtonian fluids.
34
35
\*---------------------------------------------------------------------------*/
36
37
#include "
fvCFD.H
"
38
#include "
singlePhaseTransportModel.H
"
39
#include "
pisoControl.H
"
40
41
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
43
int
main(
int
argc,
char
*argv[])
44
{
45
argList::addNote
46
(
47
"Transient solver for incompressible laminar flow"
48
" of non-Newtonian fluids."
49
);
50
51
#include "
postProcess.H
"
52
53
#include "
addCheckCaseOptions.H
"
54
#include "
setRootCaseLists.H
"
55
#include "
createTime.H
"
56
#include "
createMeshNoClear.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 "CourantNo.H"
70
71
fluid
.correct();
72
73
// Momentum predictor
74
75
fvVectorMatrix
UEqn
76
(
77
fvm::ddt(
U
)
78
+ fvm::div(
phi
,
U
)
79
- fvm::laplacian(
fluid
.nu(),
U
)
80
- (fvc::grad(
U
) & fvc::grad(
fluid
.nu()))
81
);
82
83
if
(
piso
.momentumPredictor())
84
{
85
solve
(
UEqn
== -fvc::grad(
p
));
86
}
87
88
// --- PISO loop
89
while
(
piso
.correct())
90
{
91
volScalarField
rAU
(1.0/
UEqn
.A());
92
volVectorField
HbyA
(
constrainHbyA
(
rAU
*
UEqn
.H(),
U
,
p
));
93
surfaceScalarField
phiHbyA
94
(
95
"phiHbyA"
,
96
fvc::flux(
HbyA
)
97
+ fvc::interpolate(
rAU
)*fvc::ddtCorr(
U
,
phi
)
98
);
99
100
adjustPhi
(
phiHbyA
,
U
,
p
);
101
102
// Update the pressure BCs to ensure flux consistency
103
constrainPressure
(
p
,
U
,
phiHbyA
,
rAU
);
104
105
// Non-orthogonal pressure corrector loop
106
while
(
piso
.correctNonOrthogonal())
107
{
108
// Pressure corrector
109
110
fvScalarMatrix
pEqn
111
(
112
fvm::laplacian(
rAU
,
p
) == fvc::div(
phiHbyA
)
113
);
114
115
pEqn.setReference(
pRefCell
,
pRefValue
);
116
117
pEqn.solve(
mesh
.solver(
p
.select(
piso
.finalInnerIter())));
118
119
if
(
piso
.finalNonOrthogonalIter())
120
{
121
phi
=
phiHbyA
- pEqn.flux();
122
}
123
}
124
125
#include "continuityErrs.H"
126
127
U
=
HbyA
-
rAU
*fvc::grad(
p
);
128
U
.correctBoundaryConditions();
129
}
130
131
runTime
.write();
132
133
runTime
.printExecutionTime(Info);
134
}
135
136
Info
<<
"End\n"
<<
endl
;
137
138
return
0;
139
}
140
141
142
// ************************************************************************* //
addCheckCaseOptions.H
Required Classes.
phi
surfaceScalarField & phi
Definition:
setRegionFluidFields.H:8
pRefValue
const scalar pRefValue
Definition:
setRegionFluidFields.H:37
pRefCell
const label pRefCell
Definition:
setRegionFluidFields.H:36
fluid
twoPhaseSystem & fluid
Definition:
setRegionFluidFields.H:3
constrainHbyA
constrainPressure
U
U
Definition:
pEqn.H:72
p
volScalarField & p
Definition:
createFieldRefs.H:8
UEqn
fvVectorMatrix & UEqn
Definition:
UEqn.H:13
phiHbyA
phiHbyA
Definition:
pcEqn.H:73
HbyA
HbyA
Definition:
pcEqn.H:74
mesh
dynamicFvMesh & mesh
Definition:
createDynamicFvMesh.H:6
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createMeshNoClear.H
Required Variables.
piso
pisoControl piso(mesh)
createTime.H
fvCFD.H
adjustPhi
adjustPhi(phiHbyA, U, p_rgh)
rAU
tmp< volScalarField > rAU
Definition:
initCorrectPhi.H:1
Foam::volVectorField
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition:
volFieldsFwd.H:83
Foam::fvScalarMatrix
fvMatrix< scalar > fvScalarMatrix
Definition:
fvMatricesFwd.H:45
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition:
volFieldsFwd.H:82
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition:
surfaceFieldsFwd.H:66
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition:
Ostream.H:372
Foam::fvVectorMatrix
fvMatrix< vector > fvVectorMatrix
Definition:
fvMatricesFwd.H:46
Foam::nl
constexpr char nl
The newline '\n' character (0x0a)
Definition:
Ostream.H:53
pisoControl.H
postProcess.H
Execute application functionObjects to post-process existing results.
setRootCaseLists.H
singlePhaseTransportModel.H
solve
CEqn solve()
applications
solvers
incompressible
nonNewtonianIcoFoam
nonNewtonianIcoFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.