correctUphiBCs.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) 2015-2018 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 "CorrectPhi.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
31 (
32  volVectorField& U,
34 )
35 {
36  const fvMesh& mesh = U.mesh();
37 
38  if (mesh.changing())
39  {
40  volVectorField::Boundary& Ubf = U.boundaryFieldRef();
41  surfaceScalarField::Boundary& phibf =
42  phi.boundaryFieldRef();
43 
44  forAll(Ubf, patchi)
45  {
46  if (Ubf[patchi].fixesValue())
47  {
48  Ubf[patchi].initEvaluate();
49  }
50  }
51 
52  forAll(Ubf, patchi)
53  {
54  if (Ubf[patchi].fixesValue())
55  {
56  Ubf[patchi].evaluate();
57 
58  phibf[patchi] = Ubf[patchi] & mesh.Sf().boundaryField()[patchi];
59  }
60  }
61  }
62 }
63 
64 
66 (
67  const volScalarField& rho,
68  volVectorField& U,
70 )
71 {
72  const fvMesh& mesh = U.mesh();
73 
74  if (mesh.changing())
75  {
76  volVectorField::Boundary& Ubf = U.boundaryFieldRef();
77  surfaceScalarField::Boundary& phibf =
78  phi.boundaryFieldRef();
79 
80  forAll(Ubf, patchi)
81  {
82  if (Ubf[patchi].fixesValue())
83  {
84  Ubf[patchi].initEvaluate();
85  }
86  }
87 
88  forAll(Ubf, patchi)
89  {
90  if (Ubf[patchi].fixesValue())
91  {
92  Ubf[patchi].evaluate();
93 
94  phibf[patchi] =
95  rho.boundaryField()[patchi]
96  *(
97  Ubf[patchi]
98  & mesh.Sf().boundaryField()[patchi]
99  );
100  }
101  }
102  }
103 }
104 
105 
106 // ************************************************************************* //
bool changing() const
Is mesh changing (topology changing and/or moving)
Definition: polyMesh.H:530
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const surfaceVectorField & Sf() const
Return cell face area vectors.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void correctUphiBCs(volVectorField &U, surfaceScalarField &phi)
If the mesh is moving correct the velocity BCs on the moving walls to.
dynamicFvMesh & mesh
const Mesh & mesh() const
Return mesh.
label patchi
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78