MRFZone.H
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) 2011-2022 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 Class
25  Foam::MRFZone
26 
27 Description
28  MRF zone definition based on cell zone and parameters
29  obtained from a control dictionary constructed from the given stream.
30 
31  The rotation of the MRF region is defined by an origin and axis of
32  rotation and an angular speed.
33 
34 SourceFiles
35  MRFZone.C
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #ifndef MRFZone_H
40 #define MRFZone_H
41 
42 #include "dictionary.H"
43 #include "wordList.H"
44 #include "labelList.H"
45 #include "dimensionedScalar.H"
46 #include "dimensionedVector.H"
47 #include "volFieldsFwd.H"
48 #include "surfaceFields.H"
49 #include "fvMatricesFwd.H"
50 #include "polyTopoChangeMap.H"
51 #include "Function1.H"
52 #include "autoPtr.H"
53 
54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56 namespace Foam
57 {
58 
59 /*---------------------------------------------------------------------------*\
60  Class MRFZone Declaration
61 \*---------------------------------------------------------------------------*/
62 
63 class MRFZone
64 {
65  // Private Data
66 
67  //- Reference to the mesh database
68  const fvMesh& mesh_;
69 
70  //- Name of the MRF region
71  const word name_;
72 
73  //- Coefficients dictionary
74  dictionary coeffs_;
75 
76  //- Name of cell zone
77  word cellZoneName_;
78 
79  //- Cell zone ID
80  label cellZoneID_;
81 
82  const wordReList excludedPatchNames_;
83 
84  labelList excludedPatchLabels_;
85 
86  //- Internal faces that are part of MRF
87  labelList internalFaces_;
88 
89  //- Outside faces (per patch) that move with the MRF
90  labelListList includedFaces_;
91 
92  //- Excluded faces (per patch) that do not move with the MRF
93  labelListList excludedFaces_;
94 
95  //- Origin of the axis
96  const vector origin_;
97 
98  //- Axis vector
99  vector axis_;
100 
101  //- Angular velocity (rad/sec)
103 
104 
105  // Private Member Functions
106 
107  //- Divide faces in frame according to patch
108  void setMRFFaces();
109 
110  //- Make the given absolute mass/vol flux relative within the MRF region
111  template<class RhoFieldType>
112  void makeRelativeRhoFlux
113  (
114  const RhoFieldType& rho,
116  ) const;
117 
118  //- Make the given absolute mass/vol flux relative within the MRF region
119  template<class RhoFieldType>
120  void makeRelativeRhoFlux
121  (
122  const RhoFieldType& rho,
124  ) const;
125 
126  //- Make the given absolute mass/vol flux relative within the MRF region
127  template<class RhoFieldType>
128  void makeRelativeRhoFlux
129  (
130  const RhoFieldType& rho,
131  Field<scalar>& phi,
132  const label patchi
133  ) const;
134 
135  //- Make the given relative mass/vol flux absolute within the MRF region
136  template<class RhoFieldType>
137  void makeAbsoluteRhoFlux
138  (
139  const RhoFieldType& rho,
140  surfaceScalarField& phi
141  ) const;
142 
143 
144 public:
145 
146  // Declare name of the class and its debug switch
147  ClassName("MRFZone");
148 
149 
150  // Constructors
151 
152  //- Construct from fvMesh
153  MRFZone
154  (
155  const word& name,
156  const fvMesh& mesh,
157  const dictionary& dict,
158  const word& cellZoneName = word::null
159  );
160 
161  //- Disallow default bitwise copy construction
162  MRFZone(const MRFZone&) = delete;
163 
164  //- Return clone
165  autoPtr<MRFZone> clone() const
166  {
168  return autoPtr<MRFZone>(nullptr);
169  }
170 
171 
172  // Member Functions
173 
174  //- Return const access to the MRF region name
175  inline const word& name() const;
176 
177  //- Return the current Omega vector
178  vector Omega() const;
179 
180  //- Add the Coriolis force contribution to the acceleration field
181  void addCoriolis
182  (
183  const volVectorField& U,
184  volVectorField& ddtU
185  ) const;
186 
187  //- Add the Coriolis force contribution to the momentum equation
188  // Adds to the lhs of the equation; optionally add to rhs
189  void addCoriolis
190  (
192  const bool rhs = false
193  ) const;
194 
195  //- Add the Coriolis force contribution to the momentum equation
196  // Adds to the lhs of the equation; optionally add to rhs
197  void addCoriolis
198  (
199  const volScalarField& rho,
200  fvVectorMatrix& UEqn,
201  const bool rhs = false
202  ) const;
203 
204  //- Add the centrifugal acceleration
206  (
207  volVectorField& centrifugalAcceleration
208  ) const;
209 
210  //- Make the given absolute velocity relative within the MRF region
211  void makeRelative(volVectorField& U) const;
212 
213  //- Make the given absolute flux relative within the MRF region
214  void makeRelative(surfaceScalarField& phi) const;
215 
216  //- Make the given absolute boundary flux relative
217  // within the MRF region
219 
220  //- Make the given absolute patch flux relative
221  // within the MRF region
222  void makeRelative(Field<scalar>& phi, const label patchi) const;
223 
224  //- Make the given absolute mass-flux relative within the MRF region
225  void makeRelative
226  (
227  const surfaceScalarField& rho,
228  surfaceScalarField& phi
229  ) const;
230 
231  //- Make the given relative velocity absolute within the MRF region
232  void makeAbsolute(volVectorField& U) const;
233 
234  //- Make the given relative flux absolute within the MRF region
235  void makeAbsolute(surfaceScalarField& phi) const;
236 
237  //- Make the given relative mass-flux absolute within the MRF region
238  void makeAbsolute
239  (
240  const surfaceScalarField& rho,
241  surfaceScalarField& phi
242  ) const;
243 
244  //- Correct the boundary velocity for the rotation of the MRF region
246 
247  //- Zero the MRF region of the given field
248  template<class Type>
250 
251  //- Update MRFZone faces if the mesh topology changes
252  void update();
253 
254 
255  // I-O
256 
257  //- Write
258  void writeData(Ostream& os) const;
259 
260  //- Read MRF dictionary
261  bool read(const dictionary& dict);
262 
263 
264  // Member Operators
265 
266  //- Disallow default bitwise assignment
267  void operator=(const MRFZone&) = delete;
268 };
269 
270 
271 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
272 
273 } // End namespace Foam
274 
275 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
276 
277 #ifdef NoRepository
278  #include "MRFZoneTemplates.C"
279 #endif
280 
281 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
282 
283 #include "MRFZoneI.H"
284 
285 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
286 
287 #endif
288 
289 // ************************************************************************* //
bool read(const dictionary &dict)
Read MRF dictionary.
Definition: MRFZone.C:621
Foam::surfaceFields.
void addCoriolis(const volVectorField &U, volVectorField &ddtU) const
Add the Coriolis force contribution to the acceleration field.
Definition: MRFZone.C:310
void makeRelative(volVectorField &U) const
Make the given absolute velocity relative within the MRF region.
Definition: MRFZone.C:444
dictionary dict
U
Definition: pEqn.H:72
autoPtr< MRFZone > clone() const
Return clone.
Definition: MRFZone.H:164
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
MRFZone(const word &name, const fvMesh &mesh, const dictionary &dict, const word &cellZoneName=word::null)
Construct from fvMesh.
Definition: MRFZone.C:244
MRF zone definition based on cell zone and parameters obtained from a control dictionary constructed ...
Definition: MRFZone.H:62
void makeAbsolute(volVectorField &U) const
Make the given relative velocity absolute within the MRF region.
Definition: MRFZone.C:517
void addCentrifugalAcceleration(volVectorField &centrifugalAcceleration) const
Add the centrifugal acceleration.
Definition: MRFZone.C:406
void update()
Update MRFZone faces if the mesh topology changes.
Definition: MRFZone.C:632
fvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
const word & name() const
Return const access to the MRF region name.
Definition: MRFZoneI.H:26
void zero(GeometricField< Type, fvsPatchField, surfaceMesh > &phi) const
Zero the MRF region of the given field.
static const word null
An empty word.
Definition: word.H:77
phi
Definition: correctPhi.H:3
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
void operator=(const MRFZone &)=delete
Disallow default bitwise assignment.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
void correctBoundaryVelocity(volVectorField &U) const
Correct the boundary velocity for the rotation of the MRF region.
Definition: MRFZone.C:577
vector Omega() const
Return the current Omega vector.
Definition: MRFZone.C:303
Forward declarations of fvMatrix specialisations.
label patchi
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
void writeData(Ostream &os) const
Write.
Definition: MRFZone.C:602
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:353
ClassName("MRFZone")
Namespace for OpenFOAM.
fvVectorMatrix & UEqn
Definition: UEqn.H:13