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-2024 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 "fvCellSet.H"
43 #include "volFieldsFwd.H"
44 #include "surfaceFields.H"
45 #include "omega1.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 /*---------------------------------------------------------------------------*\
53  Class MRFZone Declaration
54 \*---------------------------------------------------------------------------*/
55 
56 class MRFZone
57 {
58  // Private Data
59 
60  //- Reference to the mesh database
61  const fvMesh& mesh_;
62 
63  //- Name of the MRF region
64  const word name_;
65 
66  //- Coefficients dictionary
67  dictionary coeffs_;
68 
69  //- MRF zone cell set
70  fvCellSet cellSet_;
71 
72  //- Internal faces that are in the MRF region
73  labelList internalFaces_;
74 
75  //- Patch faces that are in the MRF region
76  labelListList patchFaces_;
77 
78  //- Origin of the axis
79  const vector origin_;
80 
81  //- Axis vector
82  vector axis_;
83 
84  //- Angular velocity
85  Function1s::omega omega_;
86 
87 
88  // Private Member Functions
89 
90  //- Divide faces in frame according to patch
91  void setMRFFaces();
92 
93  //- Check that the case has been updated with correct MRF BCs
94  void checkMRFBCs(const volVectorField& U) const;
95 
96  //- Make the given absolute mass/vol flux relative within the MRF region
97  template<class RhoFieldType>
98  void makeRelativeRhoFlux
99  (
100  const RhoFieldType& rho,
101  surfaceScalarField& phi
102  ) const;
103 
104  //- Make the given absolute mass/vol flux relative within the MRF region
105  template<class RhoFieldType>
106  void makeRelativeRhoFlux
107  (
108  const RhoFieldType& rho,
110  ) const;
111 
112  //- Make the given absolute mass/vol flux relative within the MRF region
113  template<class RhoFieldType>
114  void makeRelativeRhoFlux
115  (
116  const RhoFieldType& rho,
117  Field<scalar>& phi,
118  const label patchi
119  ) const;
120 
121  //- Make the given relative mass/vol flux absolute within the MRF region
122  template<class RhoFieldType>
123  void makeAbsoluteRhoFlux
124  (
125  const RhoFieldType& rho,
126  surfaceScalarField& phi
127  ) const;
128 
129 
130 public:
131 
132  // Declare name of the class and its debug switch
133  ClassName("MRFZone");
134 
135 
136  // Constructors
137 
138  //- Construct from fvMesh
139  MRFZone
140  (
141  const word& name,
142  const fvMesh& mesh,
143  const dictionary& dict
144  );
145 
146  //- Disallow default bitwise copy construction
147  MRFZone(const MRFZone&) = delete;
148 
149  //- Return clone
150  autoPtr<MRFZone> clone() const
151  {
153  return autoPtr<MRFZone>(nullptr);
154  }
155 
156 
157  // Member Functions
158 
159  //- Return const access to the MRF region name
160  inline const word& name() const;
161 
162  //- Return the current Omega vector
163  vector Omega() const;
164 
165  //- Add the Coriolis force contribution to the acceleration field
166  void addCoriolis
167  (
168  const volVectorField& U,
169  volVectorField& ddtU
170  ) const;
171 
172  //- Add the centrifugal acceleration
174  (
175  volVectorField& centrifugalAcceleration
176  ) const;
177 
178  //- Make the given absolute velocity relative within the MRF region
179  void makeRelative(volVectorField& U) const;
180 
181  //- Make the given absolute flux relative within the MRF region
182  void makeRelative(surfaceScalarField& phi) const;
183 
184  //- Make the given absolute boundary flux relative
185  // within the MRF region
187 
188  //- Make the given absolute patch flux relative
189  // within the MRF region
190  void makeRelative(Field<scalar>& phi, const label patchi) const;
191 
192  //- Make the given absolute mass-flux relative within the MRF region
193  void makeRelative
194  (
195  const surfaceScalarField& rho,
196  surfaceScalarField& phi
197  ) const;
198 
199  //- Make the given relative patch velocity relative
200  // within the MRF region
201  void makeRelative(Field<vector>& Up, const label patchi) const;
202 
203  //- Make the given relative velocity absolute within the MRF region
204  void makeAbsolute(volVectorField& U) const;
205 
206  //- Make the given relative flux absolute within the MRF region
207  void makeAbsolute(surfaceScalarField& phi) const;
208 
209  //- Make the given relative mass-flux absolute within the MRF region
210  void makeAbsolute
211  (
212  const surfaceScalarField& rho,
213  surfaceScalarField& phi
214  ) const;
215 
216  //- Make the given relative patch velocity absolute
217  // within the MRF region
218  void makeAbsolute(Field<vector>& Up, const label patchi) const;
219 
220  //- Zero the MRF region of the given field
221  template<class Type>
222  void zero(SurfaceField<Type>& phi) const;
223 
224  //- Update MRFZone faces if the mesh topology changes
225  void update();
226 
227 
228  // I-O
229 
230  //- Read MRF dictionary
231  bool read(const dictionary& dict);
232 
233 
234  // Member Operators
235 
236  //- Disallow default bitwise assignment
237  void operator=(const MRFZone&) = delete;
238 };
239 
240 
241 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
242 
243 } // End namespace Foam
244 
245 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
246 
247 #ifdef NoRepository
248  #include "MRFZoneTemplates.C"
249 #endif
250 
251 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
252 
253 #include "MRFZoneI.H"
254 
255 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
256 
257 #endif
258 
259 // ************************************************************************* //
Generic field type.
Definition: FieldField.H:77
Convenience class to handle the input of time-varying rotational speed. Reads an omega Function1 entr...
Definition: omega1.H:132
Generic GeometricField class.
MRF zone definition based on cell zone and parameters obtained from a control dictionary constructed ...
Definition: MRFZone.H:56
vector Omega() const
Return the current Omega vector.
Definition: MRFZone.C:187
void zero(SurfaceField< Type > &phi) const
Zero the MRF region of the given field.
ClassName("MRFZone")
autoPtr< MRFZone > clone() const
Return clone.
Definition: MRFZone.H:149
bool read(const dictionary &dict)
Read MRF dictionary.
Definition: MRFZone.C:362
void makeAbsolute(volVectorField &U) const
Make the given relative velocity absolute within the MRF region.
Definition: MRFZone.C:311
MRFZone(const word &name, const fvMesh &mesh, const dictionary &dict)
Construct from fvMesh.
Definition: MRFZone.C:166
void makeRelative(volVectorField &U) const
Make the given absolute velocity relative within the MRF region.
Definition: MRFZone.C:247
void operator=(const MRFZone &)=delete
Disallow default bitwise assignment.
void addCentrifugalAcceleration(volVectorField &centrifugalAcceleration) const
Add the centrifugal acceleration.
Definition: MRFZone.C:216
void update()
Update MRFZone faces if the mesh topology changes.
Definition: MRFZone.C:372
void addCoriolis(const volVectorField &U, volVectorField &ddtU) const
Add the Coriolis force contribution to the acceleration field.
Definition: MRFZone.C:194
const word & name() const
Return const access to the MRF region name.
Definition: MRFZoneI.H:26
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
General run-time selected cell set selection class for fvMesh.
Definition: fvCellSet.H:88
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
label patchi
U
Definition: pEqn.H:72
Namespace for OpenFOAM.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dictionary dict
Foam::surfaceFields.