DAFoam
v4.0.2
Discrete Adjoint with OpenFOAM
dafoam
src
adjoint
DASolver
DAPimpleFoam
pEqnPimple.H
Go to the documentation of this file.
1
volScalarField
rAU
(1.0 /
UEqn
.A());
2
//***************** NOTE *******************
3
// constrainHbyA has been used since OpenFOAM-v1606; however, it may degrade the accuracy of derivatives
4
// because constraining variables will create discontinuity. Here we have a option to use the old
5
// implementation in OpenFOAM-3.0+ and before (no constraint for HbyA)
6
autoPtr<volVectorField>
HbyAPtr
=
nullptr
;
7
label
useConstrainHbyA
= daOptionPtr_->getOption<label>(
"useConstrainHbyA"
);
8
if
(
useConstrainHbyA
)
9
{
10
HbyAPtr
.reset(
new
volVectorField(constrainHbyA(
rAU
*
UEqn
.H(),
U
,
p
)));
11
}
12
else
13
{
14
HbyAPtr
.reset(
new
volVectorField(
"HbyA"
,
U
));
15
HbyAPtr
() =
rAU
*
UEqn
.H();
16
}
17
volVectorField&
HbyA
=
HbyAPtr
();
18
surfaceScalarField
phiHbyA
(
19
"phiHbyA"
,
20
fvc::flux(
HbyA
));
21
22
if
(
p
.needReference())
23
{
24
adjustPhi
(
phiHbyA
,
U
,
p
);
25
}
26
27
tmp<volScalarField>
rAtU
(
rAU
);
28
29
if
(
pimple
.consistent())
30
{
31
rAtU
= 1.0 / max(1.0 /
rAU
-
UEqn
.H1(), 0.1 /
rAU
);
32
phiHbyA
+=
33
fvc::interpolate(
rAtU
() -
rAU
) * fvc::snGrad(
p
) *
mesh
.magSf();
34
HbyA
-= (
rAU
-
rAtU
()) * fvc::grad(
p
);
35
}
36
37
if
(
pimple
.nCorrPISO() <= 1)
38
{
39
tUEqn
.clear();
40
}
41
42
// Non-orthogonal pressure corrector loop
43
while
(
pimple
.correctNonOrthogonal())
44
{
45
// Pressure corrector
46
fvScalarMatrix pEqn(
47
fvm::laplacian(
rAtU
(),
p
) == fvc::div(
phiHbyA
));
48
49
pEqn.setReference(
pRefCell
,
pRefValue
);
50
51
// get the solver performance info such as initial
52
// and final residuals
53
SolverPerformance<scalar> solverP = pEqn.solve(
mesh
.solver(
p
.select(
pimple
.finalInnerIter())));
54
55
DAUtility::primalResidualControl(solverP, pimplePrintToScreen,
"p"
, daGlobalVarPtr_->primalMaxRes);
56
57
if
(
pimple
.finalNonOrthogonalIter())
58
{
59
phi
=
phiHbyA
- pEqn.flux();
60
}
61
}
62
63
if
(pimplePrintToScreen)
64
{
65
#include "
continuityErrsPython.H
"
66
}
67
68
// Explicitly relax pressure for momentum corrector
69
p
.relax();
70
71
// bound p
72
DAUtility::boundVar(
allOptions
,
p
, pimplePrintToScreen);
73
74
U
=
HbyA
-
rAtU
() * fvc::grad(
p
);
75
// bound U
76
DAUtility::boundVar(
allOptions
,
U
, pimplePrintToScreen);
77
U
.correctBoundaryConditions();
rAtU
tmp< volScalarField > rAtU(rAU)
pRefValue
scalar & pRefValue
Definition:
createRefsPimpleDyM.H:12
pimple
pimpleControlDF & pimple
Definition:
createRefsPimpleDyM.H:5
useConstrainHbyA
label useConstrainHbyA
Definition:
pEqnPimple.H:7
UEqn
fvVectorMatrix & UEqn
Definition:
UEqnPimpleDyM.H:15
rAU
volScalarField rAU(1.0/UEqn.A())
HbyA
volVectorField & HbyA
Definition:
pEqnPimple.H:17
continuityErrsPython.H
pRefCell
label & pRefCell
Definition:
createRefsPimpleDyM.H:11
mesh
fvMesh & mesh
Definition:
createRefsHeatTransfer.H:4
p
volScalarField & p
Definition:
createRefsPimpleDyM.H:6
phi
surfaceScalarField & phi
Definition:
createRefsPimpleDyM.H:8
tUEqn
tmp< fvVectorMatrix > tUEqn(fvm::ddt(U)+fvm::div(phi, U)+turbulencePtr_->divDevReff(U) - fvSource)
allOptions
const dictionary & allOptions
Definition:
createRefsPimpleDyM.H:14
adjustPhi
adjustPhi(phiHbyA, U, p)
HbyAPtr
autoPtr< volVectorField > HbyAPtr
Definition:
pEqnPimple.H:6
U
U
Definition:
pEqnPimple.H:74
phiHbyA
surfaceScalarField phiHbyA("phiHbyA", fvc::flux(HbyA))
Generated by
1.8.17