constrainPressure.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 | Copyright (C) 2016 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 
26 #include "constrainPressure.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "geometricOneField.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 template<class RhoType, class RAUType, class MRFType>
36 (
37  volScalarField& p,
38  const RhoType& rho,
39  const volVectorField& U,
40  const surfaceScalarField& phiHbyA,
41  const RAUType& rhorAU,
42  const MRFType& MRF
43 )
44 {
45  const fvMesh& mesh = p.mesh();
46 
47  volScalarField::Boundary& pBf = p.boundaryFieldRef();
48 
49  const volVectorField::Boundary& UBf = U.boundaryField();
50  const surfaceScalarField::Boundary& phiHbyABf =
51  phiHbyA.boundaryField();
52  const typename RAUType::Boundary& rhorAUBf =
53  rhorAU.boundaryField();
54  const surfaceVectorField::Boundary& SfBf =
55  mesh.Sf().boundaryField();
56  const surfaceScalarField::Boundary& magSfBf =
57  mesh.magSf().boundaryField();
58 
59  forAll(pBf, patchi)
60  {
61  if (isA<fixedFluxPressureFvPatchScalarField>(pBf[patchi]))
62  {
63  refCast<fixedFluxPressureFvPatchScalarField>
64  (
65  pBf[patchi]
66  ).updateCoeffs
67  (
68  (
69  phiHbyABf[patchi]
70  - rho.boundaryField()[patchi]
71  *MRF.relative(SfBf[patchi] & UBf[patchi], patchi)
72  )
73  /(magSfBf[patchi]*rhorAUBf[patchi])
74  );
75  }
76  }
77 }
78 
79 
80 template<class RAUType>
82 (
83  volScalarField& p,
84  const volScalarField& rho,
85  const volVectorField& U,
86  const surfaceScalarField& phiHbyA,
87  const RAUType& rAU
88 )
89 {
90  constrainPressure(p, rho, U, phiHbyA, rAU, NullMRF());
91 }
92 
93 
94 template<class RAUType, class MRFType>
96 (
97  volScalarField& p,
98  const volVectorField& U,
99  const surfaceScalarField& phiHbyA,
100  const RAUType& rAU,
101  const MRFType& MRF
102 )
103 {
104  constrainPressure(p, geometricOneField(), U, phiHbyA, rAU, MRF);
105 }
106 
107 
108 template<class RAUType>
110 (
111  volScalarField& p,
112  const volVectorField& U,
113  const surfaceScalarField& phiHbyA,
114  const RAUType& rAU
115 )
116 {
117  constrainPressure(p, U, phiHbyA, rAU, NullMRF());
118 }
119 
120 
121 // ************************************************************************* //
Foam::surfaceFields.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const surfaceVectorField & Sf() const
Return cell face area vectors.
void constrainPressure(volScalarField &p, const RhoType &rho, const volVectorField &U, const surfaceScalarField &phiHbyA, const RAUType &rhorAU, const MRFType &MRF)
dynamicFvMesh & mesh
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
label patchi
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const Mesh & mesh() const
Return mesh.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78