correctPhi.H
Go to the documentation of this file.
1if (mesh.changing())
2{
3 volVectorField::Boundary& bfld = U.boundaryFieldRef();
4 forAll(bfld, patchi)
5 {
6 if (bfld[patchi].fixesValue())
7 {
8 bfld[patchi].initEvaluate();
9 }
10 }
11
12 surfaceScalarField::Boundary& phiBfld = phi.boundaryFieldRef();
13 forAll(bfld, patchi)
14 {
15 if (bfld[patchi].fixesValue())
16 {
17 bfld[patchi].evaluate();
18
19 phiBfld[patchi] =
20 rho.boundaryField()[patchi]
21 * (
22 bfld[patchi]
23 & mesh.Sf().boundaryField()[patchi]
24 );
25 }
26 }
27}
28 // Initialize BCs list for pcorr to zero-gradient
29 wordList pcorrTypes
30 (
31 p.boundaryField().size(),
33 );
34
35 // Set BCs of pcorr to fixed-value for patches at which p is fixed
36 forAll(p.boundaryField(), patchi)
37 {
38 if (p.boundaryField()[patchi].fixesValue())
39 {
41 }
42 }
43
44 volScalarField pcorr
45 (
46 IOobject
47 (
48 "pcorr",
49 runTime.timeName(),
50 mesh,
51 IOobject::NO_READ,
52 IOobject::NO_WRITE
53 ),
54 mesh,
55 dimensionedScalar(p.dimensions(), Zero),
57 );
58
59 mesh.setFluxRequired(pcorr.name());
60
61{
62 dimensionedScalar rAUf("rAUf", dimTime, 1.0);
63
64 while (pimple.correctNonOrthogonal())
65 {
66 fvScalarMatrix pcorrEqn
67 (
68 fvm::ddt(psi, pcorr)
69 + fvc::div(phi)
70 - fvm::laplacian(rAUf, pcorr)
71 ==
72 divrhoU()
73 );
74
75 //pcorrEqn.solve(mesh.solver(pcorr.select(pimple.finalInnerIter())));
76 //Bypass virtual layer
77 const dictionary& d = mesh.solver
78 (
79 pcorr.select
80 (
81 pimple.finalInnerIter()
82 )
83 );
84 mesh.fvMesh::solve(pcorrEqn, d);
85
86 if (pimple.finalNonOrthogonalIter())
87 {
88 phi += pcorrEqn.flux();
89 }
90 }
91}
surfaceScalarField & phi
pimpleControl & pimple
U
Definition: pEqn.H:72
volScalarField & p
const volScalarField & psi
volScalarField pcorr(IOobject("pcorr", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), mesh, dimensionedScalar(p.dimensions(), Zero), pcorrTypes)
wordList pcorrTypes(p.boundaryField().size(), zeroGradientFvPatchScalarField::typeName)
dynamicFvMesh & mesh
engineTime & runTime
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU))
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
static const char *const typeName
The type name used in ensight case files.