correctPressure.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2023 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
27 #include "constrainHbyA.H"
28 #include "constrainPressure.H"
29 #include "adjustPhi.H"
30 #include "fvcMeshPhi.H"
31 #include "fvcFlux.H"
32 #include "fvcDdt.H"
33 #include "fvcSnGrad.H"
34 #include "fvcReconstruct.H"
35 #include "fvmLaplacian.H"
36 
37 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
38 
40 {
44 
45  fvVectorMatrix& UcEqn = tUcEqn.ref();
46 
47  const volScalarField rAUc(1.0/UcEqn.A());
48  const volScalarField r1ADUc(1/(1 + rAUc*Dc()));
49 
50  const surfaceScalarField rAUcf(fvc::interpolate(rAUc));
51  const surfaceScalarField r1ADUcf(1/(1 + rAUcf*Dcf()));
52  const surfaceScalarField rADUcf("Dp", r1ADUcf*rAUcf);
53 
54  volVectorField HbyA(constrainHbyA(rAUc*UcEqn.H(), Uc, p));
55 
56  surfaceScalarField phiHbyAD
57  (
58  "phiHbyAD",
59  (
60  r1ADUcf
61  *(
63  + alphacf*rAUcf*fvc::ddtCorr(Uc, phic, Ucf)
64  )
65  )
66  );
67 
68  if (p.needReference())
69  {
70  fvc::makeRelative(phiHbyAD, Uc);
71  adjustPhi(phiHbyAD, Uc, p);
72  fvc::makeAbsolute(phiHbyAD, Uc);
73  }
74 
75  // Face buoyancy force
76  const surfaceScalarField Fgf(g & mesh.Sf());
77 
78  phiHbyAD += rADUcf*(Fgf + Dcf()*phid());
79 
80  // Update the pressure BCs to ensure flux consistency
81  constrainPressure(p, Uc, phiHbyAD, rADUcf);
82 
83  // Non-orthogonal pressure corrector loop
85  {
86  fvScalarMatrix pEqn
87  (
88  fvm::laplacian(alphacf*rADUcf, p)
89  ==
91  + fvc::div(alphacf*phiHbyAD)
92  + fvModels().sourceProxy(alphac, p)
93  );
94 
95  pEqn.setReference
96  (
99  );
100 
101  pEqn.solve();
102 
104  {
105  phic = phiHbyAD - pEqn.flux()/alphacf;
106 
107  // Explicitly relax pressure for momentum corrector
108  p.relax();
109 
110  Uc =
111  r1ADUc
112  *(
113  HbyA
114  + rAUc
115  *(
117  (
118  Fgf - pEqn.flux()/alphacf/rADUcf
119  - Dcf()*(phic - phid())
120  )
121  + Dc()*fvc::reconstruct(phic - phid())
122  + Fd()
123  )
124  );
125 
128 
129  // Correct Ucf if the mesh is moving
131 
132  // Make the fluxes relative to the mesh motion
134  }
135  }
136 }
137 
138 
139 // ************************************************************************* //
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
Generic GeometricField class.
void relax(const scalar alpha)
Relax field (for steady-state solution).
bool needReference() const
Does the field need a reference level for solution.
void correctBoundaryConditions()
Correct boundary field.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:808
tmp< VolField< Type > > H() const
Return the H operation source.
Definition: fvMatrix.C:868
tmp< SurfaceField< Type > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:964
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:563
const surfaceVectorField & Sf() const
Return cell face area vectors.
bool finalNonOrthogonalIter() const
Flag to indicate the last non-orthogonal iteration.
bool correctNonOrthogonal()
Non-orthogonal corrector loop.
Provides controls for the pressure reference in closed-volume simulations.
scalar refValue() const
Return the pressure reference level.
label refCell() const
Return the cell in which the reference pressure is set.
Foam::fvModels & fvModels() const
Return the fvModels that are created on demand.
Definition: solver.C:96
pimpleNoLoopControl pimple
PIMPLE inner-loop controls.
Definition: solver.H:107
Foam::fvConstraints & fvConstraints() const
Return the fvConstraints that are created on demand.
Definition: solver.C:107
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
autoPtr< surfaceScalarField > phid
Effective volumetric flux of the dispersed phase.
uniformDimensionedVectorField g
Acceleration due to gravity.
const volVectorField & Uc
Reference to the continuous phase velocity field.
surfaceScalarField phic_
Continuous phase flux field.
surfaceScalarField alphacf
Interpolated continuous phase-fraction.
const volScalarField & alphac
Reference continuous phase-fraction.
volVectorField Uc_
Continuous phase velocity field.
autoPtr< surfaceVectorField > Ucf
Pointer to the surface momentum field.
tmp< fvVectorMatrix > tUcEqn
Cached momentum matrix.
autoPtr< volVectorField > Fd
Dispersed phase drag force.
virtual void correctPressure()
Construct the pressure equation.
autoPtr< surfaceScalarField > Dcf
Face continuous-dispersed phase drag coefficient.
const surfaceScalarField & phic
Reference to the continuous phase volumetric-flux field.
const volScalarField & p
Reference to the pressure field.
Calculate the first temporal derivative.
Calculate the face-flux of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Reconstruct volField from a face flux field.
Calculate the snGrad of the given volField.
Calculate the matrix for the laplacian of the field.
volVectorField & HbyA
Definition: pEqn.H:13
tmp< SurfaceField< typename innerProduct< vector, Type >::type > > flux(const VolField< Type > &vf)
Return the face-flux field obtained from the given volVectorField.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
void correctUf(autoPtr< surfaceVectorField > &Uf, const volVectorField &U, const surfaceScalarField &phi, const MRFType &MRF)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
void makeAbsolute(surfaceScalarField &phi, const volVectorField &U)
Make the given flux absolute.
Definition: fvcMeshPhi.C:128
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:89
tmp< SurfaceField< typename Foam::flux< Type >::type > > ddtCorr(const VolField< Type > &U, const SurfaceField< Type > &Uf)
Definition: fvcDdt.C:196
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
void constrainPressure(volScalarField &p, const RhoType &rho, const volVectorField &U, const surfaceScalarField &phiHbyA, const RAUType &rhorAU, const MRFType &MRF)
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
Definition: adjustPhi.C:35
void constrainHbyA(volVectorField &HbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:34