MRFZoneListDF.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2 
3  DAFoam : Discrete Adjoint with OpenFOAM
4  Version : v3
5 
6  Description:
7  A modified version of MRF that allows changing the rotation speed at
8  the runtime to enable derivatives wrt the rotation speed
9 
10  This class is modified from OpenFOAM's source code
11  src/finiteVolume/cfdTools/general/MRF
12 
13  OpenFOAM: The Open Source CFD Toolbox
14 
15  Copyright (C): 2011-2016 OpenFOAM Foundation
16 
17  OpenFOAM License:
18 
19  OpenFOAM is free software: you can redistribute it and/or modify it
20  under the terms of the GNU General Public License as published by
21  the Free Software Foundation, either version 3 of the License, or
22  (at your option) any later version.
23 
24  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
25  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
26  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
27  for more details.
28 
29  You should have received a copy of the GNU General Public License
30  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #include "MRFZoneListDF.H"
35 #include "volFields.H"
36 #include "fixedValueFvsPatchFields.H"
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 Foam::MRFZoneListDF::MRFZoneListDF(
41  const fvMesh& mesh,
42  const dictionary& dict)
43  : PtrList<MRFZoneDF>(),
44  mesh_(mesh)
45 {
46  reset(dict);
47 
48  active(true);
49 }
50 
51 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
52 
54 {
55 }
56 
57 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
58 
59 bool Foam::MRFZoneListDF::active(const bool warn) const
60 {
61  bool a = false;
62  forAll(*this, i)
63  {
64  a = a || this->operator[](i).active();
65  }
66 
67  if (warn && this->size() && !a)
68  {
69  Info << " No MRF zones active" << endl;
70  }
71 
72  return a;
73 }
74 
75 void Foam::MRFZoneListDF::reset(const dictionary& dict)
76 {
77  label count = 0;
78  for (const entry& dEntry : dict)
79  {
80  if (dEntry.isDict())
81  {
82  ++count;
83  }
84  }
85 
86  this->resize(count);
87 
88  count = 0;
89  for (const entry& dEntry : dict)
90  {
91  if (dEntry.isDict())
92  {
93  const word& name = dEntry.keyword();
94  const dictionary& modelDict = dEntry.dict();
95 
96  Info << " creating MRF zone: " << name << endl;
97 
98  this->set(
99  count++,
100  new MRFZoneDF(name, mesh_, modelDict));
101  }
102  }
103 }
104 
105 /*
106 bool Foam::MRFZoneListDF::read(const dictionary& dict)
107 {
108  // do nothing
109  return true;
110 }
111 
112 bool Foam::MRFZoneListDF::writeData(Ostream& os) const
113 {
114  // do nothing
115 
116  return true;
117 }
118 */
119 
120 Foam::tmp<Foam::volVectorField> Foam::MRFZoneListDF::DDt(
121  const volVectorField& U) const
122 {
123  tmp<volVectorField> tacceleration(
124  new volVectorField(
125  IOobject(
126  "MRFZoneListDF:acceleration",
127  U.mesh().time().timeName(),
128  U.mesh()),
129  U.mesh(),
130  dimensionedVector(U.dimensions() / dimTime, Zero)));
131  volVectorField& acceleration = tacceleration.ref();
132 
133  forAll(*this, i)
134  {
135  operator[](i).addCoriolis(U, acceleration);
136  }
137 
138  return tacceleration;
139 }
140 
141 Foam::tmp<Foam::volVectorField> Foam::MRFZoneListDF::DDt(
142  const volScalarField& rho,
143  const volVectorField& U) const
144 {
145  return rho * DDt(U);
146 }
147 
148 void Foam::MRFZoneListDF::makeRelative(volVectorField& U) const
149 {
150  forAll(*this, i)
151  {
152  operator[](i).makeRelative(U);
153  }
154 }
155 
156 void Foam::MRFZoneListDF::makeRelative(surfaceScalarField& phi) const
157 {
158  forAll(*this, i)
159  {
160  operator[](i).makeRelative(phi);
161  }
162 }
163 
164 Foam::tmp<Foam::surfaceScalarField> Foam::MRFZoneListDF::relative(
165  const tmp<surfaceScalarField>& tphi) const
166 {
167  // this function should not be called, if yes, return an error
168  FatalErrorInFunction
169  << "This is unexpected!"
170  << exit(FatalError);
171  return tmp<surfaceScalarField>(tphi, true);
172 }
173 
174 Foam::tmp<Foam::FieldField<Foam::fvsPatchField, Foam::scalar>>
176  const tmp<FieldField<fvsPatchField, scalar>>& tphi) const
177 {
178  // this function should not be called, if yes, return an error
179  FatalErrorInFunction
180  << "This is unexpected!"
181  << exit(FatalError);
182  return tmp<FieldField<fvsPatchField, scalar>>(tphi, true);
183 }
184 
185 Foam::tmp<Foam::Field<Foam::scalar>>
187  const tmp<Field<scalar>>& tphi,
188  const label patchi) const
189 {
190  // this function should not be called, if yes, return an error
191  FatalErrorInFunction
192  << "This is unexpected!"
193  << exit(FatalError);
194  return tmp<Field<scalar>>(tphi, true);
195 }
196 
198  const surfaceScalarField& rho,
199  surfaceScalarField& phi) const
200 {
201  forAll(*this, i)
202  {
203  operator[](i).makeRelative(rho, phi);
204  }
205 }
206 
208 {
209  forAll(*this, i)
210  {
211  operator[](i).correctBoundaryVelocity(U);
212  }
213 }
214 
215 const Foam::scalar& Foam::MRFZoneListDF::getOmegaRef() const
216 {
217  label nObjs = 0;
218  forAll(*this, i)
219  {
220  nObjs++;
221  }
222 
223  if (nObjs > 1)
224  {
225  FatalErrorInFunction
226  << "Do not support more than one MRF zones!"
227  << exit(FatalError);
228  }
229 
230  return operator[](0).getOmegaRef();
231 }
232 
233 /*
234 void Foam::MRFZoneListDF::update()
235 {
236  if (mesh_.topoChanging())
237  {
238  forAll(*this, i)
239  {
240  operator[](i).update();
241  }
242  }
243 }
244 */
245 
246 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
247 
248 Foam::Ostream& Foam::operator<<(
249  Ostream& os,
250  const MRFZoneListDF& models)
251 {
252  //models.writeData(os);
253  return os;
254 }
255 
256 // ************************************************************************* //
Foam::MRFZoneListDF::DDt
tmp< volVectorField > DDt(const volVectorField &U) const
Definition: MRFZoneListDF.C:120
Foam::MRFZoneListDF::~MRFZoneListDF
~MRFZoneListDF()
Definition: MRFZoneListDF.C:53
U
U
Definition: pEqnPimpleDyM.H:60
forAll
forAll(pseudoP.boundaryField(), patchI)
Definition: solvePseudoPEqn.H:10
Foam::MRFZoneListDF::relative
tmp< surfaceScalarField > relative(const tmp< surfaceScalarField > &phi) const
Definition: MRFZoneListDF.C:164
Foam::MRFZoneListDF
Definition: MRFZoneListDF.H:55
Foam::MRFZoneDF
Definition: MRFZoneDF.H:61
Foam::operator<<
Ostream & operator<<(Ostream &os, const MRFZoneListDF &models)
Definition: MRFZoneListDF.C:248
mesh
fvMesh & mesh
Definition: createRefsHeatTransfer.H:4
phi
surfaceScalarField & phi
Definition: createRefsPimple.H:8
MRFZoneListDF.H
Foam::MRFZoneListDF::reset
void reset(const dictionary &dict)
Definition: MRFZoneListDF.C:75
Foam::MRFZoneListDF::makeRelative
void makeRelative(volVectorField &U) const
Definition: MRFZoneListDF.C:148
Foam::MRFZoneListDF::getOmegaRef
const scalar & getOmegaRef() const
Definition: MRFZoneListDF.C:215
Foam::MRFZoneListDF::correctBoundaryVelocity
void correctBoundaryVelocity(volVectorField &U) const
Definition: MRFZoneListDF.C:207
Foam::MRFZoneListDF::active
bool active(const bool warn=false) const
Definition: MRFZoneListDF.C:59
rho
volScalarField & rho
Definition: createRefsRhoSimpleC.H:8