fvOptionListTemplates.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) 2011-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 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
27 
28 template<class Type>
30 (
32  const word& fieldName,
33  const dimensionSet& ds
34 )
35 {
36  checkApplied();
37 
38  tmp<fvMatrix<Type>> tmtx(new fvMatrix<Type>(field, ds));
39  fvMatrix<Type>& mtx = tmtx.ref();
40 
41  forAll(*this, i)
42  {
43  option& source = this->operator[](i);
44 
45  label fieldi = source.applyToField(fieldName);
46 
47  if (fieldi != -1)
48  {
49  source.setApplied(fieldi);
50 
51  if (source.isActive())
52  {
53  if (debug)
54  {
55  Info<< "Applying source " << source.name() << " to field "
56  << fieldName << endl;
57  }
58 
59  source.addSup(mtx, fieldi);
60  }
61  }
62  }
63 
64  return tmtx;
65 }
66 
67 
68 template<class Type>
70 (
72 )
73 {
74  return this->operator()(field, field.name());
75 }
76 
77 
78 template<class Type>
80 (
82  const word& fieldName
83 )
84 {
85  return source(field, fieldName, field.dimensions()/dimTime*dimVolume);
86 }
87 
88 
89 template<class Type>
91 (
92  const volScalarField& rho,
94 )
95 {
96  return this->operator()(rho, field, field.name());
97 }
98 
99 
100 template<class Type>
102 (
103  const volScalarField& rho,
105  const word& fieldName
106 )
107 {
108  checkApplied();
109 
110  const dimensionSet ds
111  (
112  rho.dimensions()*field.dimensions()/dimTime*dimVolume
113  );
114 
115  tmp<fvMatrix<Type>> tmtx(new fvMatrix<Type>(field, ds));
116  fvMatrix<Type>& mtx = tmtx.ref();
117 
118  forAll(*this, i)
119  {
120  option& source = this->operator[](i);
121 
122  label fieldi = source.applyToField(fieldName);
123 
124  if (fieldi != -1)
125  {
126  source.setApplied(fieldi);
127 
128  if (source.isActive())
129  {
130  if (debug)
131  {
132  Info<< "Applying source " << source.name() << " to field "
133  << fieldName << endl;
134  }
135 
136  source.addSup(rho, mtx, fieldi);
137  }
138  }
139  }
140 
141  return tmtx;
142 }
143 
144 
145 template<class Type>
147 (
148  const volScalarField& alpha,
149  const volScalarField& rho,
151 )
152 {
153  return this->operator()(alpha, rho, field, field.name());
154 }
155 
156 
157 template<class Type>
159 (
160  const volScalarField& alpha,
161  const volScalarField& rho,
163  const word& fieldName
164 )
165 {
166  checkApplied();
167 
168  const dimensionSet ds
169  (
170  alpha.dimensions()*rho.dimensions()*field.dimensions()
172  );
173 
174  tmp<fvMatrix<Type>> tmtx(new fvMatrix<Type>(field, ds));
175  fvMatrix<Type>& mtx = tmtx.ref();
176 
177  forAll(*this, i)
178  {
179  option& source = this->operator[](i);
180 
181  label fieldi = source.applyToField(fieldName);
182 
183  if (fieldi != -1)
184  {
185  source.setApplied(fieldi);
186 
187  if (source.isActive())
188  {
189  if (debug)
190  {
191  Info<< "Applying source " << source.name() << " to field "
192  << fieldName << endl;
193  }
194 
195  source.addSup(alpha, rho, mtx, fieldi);
196  }
197  }
198  }
199 
200  return tmtx;
201 }
202 
203 
204 template<class Type>
206 (
207  const geometricOneField& alpha,
208  const geometricOneField& rho,
210 )
211 {
212  return this->operator()(field, field.name());
213 }
214 
215 
216 template<class Type>
218 (
219  const volScalarField& alpha,
220  const geometricOneField& rho,
222 )
223 {
225  (
226  IOobject
227  (
228  "one",
229  this->mesh_.time().timeName(),
230  this->mesh_,
231  IOobject::NO_READ,
232  IOobject::NO_WRITE,
233  false
234  ),
235  this->mesh_,
237  );
238 
239  return this->operator()(alpha, one, field, field.name());
240 }
241 
242 
243 template<class Type>
245 (
246  const geometricOneField& alpha,
247  const volScalarField& rho,
249 )
250 {
251  return this->operator()(rho, field, field.name());
252 }
253 
254 
255 template<class Type>
257 (
259 )
260 {
261  return this->d2dt2(field, field.name());
262 }
263 
264 
265 template<class Type>
267 (
269  const word& fieldName
270 )
271 {
272  return source(field, fieldName, field.dimensions()/sqr(dimTime)*dimVolume);
273 }
274 
275 
276 template<class Type>
278 {
279  checkApplied();
280 
281  forAll(*this, i)
282  {
283  option& source = this->operator[](i);
284 
285  label fieldi = source.applyToField(eqn.psi().name());
286 
287  if (fieldi != -1)
288  {
289  source.setApplied(fieldi);
290 
291  if (source.isActive())
292  {
293  if (debug)
294  {
295  Info<< "Applying constraint " << source.name()
296  << " to field " << eqn.psi().name() << endl;
297  }
298 
299  source.constrain(eqn, fieldi);
300  }
301  }
302  }
303 }
304 
305 
306 template<class Type>
308 (
310 )
311 {
312  const word& fieldName = field.name();
313 
314  forAll(*this, i)
315  {
316  option& source = this->operator[](i);
317 
318  label fieldi = source.applyToField(fieldName);
319 
320  if (fieldi != -1)
321  {
322  source.setApplied(fieldi);
323 
324  if (source.isActive())
325  {
326  if (debug)
327  {
328  Info<< "Correcting source " << source.name()
329  << " for field " << fieldName << endl;
330  }
331 
332  source.correct(field);
333  }
334  }
335  }
336 }
337 
338 
339 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
tmp< GeometricField< Type, fvPatchField, volMesh > > d2dt2(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcD2dt2.C:45
virtual void correct(volScalarField &field)
Definition: fvOption.C:306
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
const word & name() const
Return name.
Definition: IOobject.H:295
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:282
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const T & operator[](const label) const
Return element const reference.
Definition: UPtrListI.H:96
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Generic GeometricField class.
tmp< fvMatrix< Type > > source(GeometricField< Type, fvPatchField, volMesh > &field, const word &fieldName, const dimensionSet &ds)
Return source for equation with specified name and dimensions.
void setApplied(const label fieldi)
Set the applied flag to true for field index fieldi.
Definition: fvOptionI.H:52
friend Ostream & operator(Ostream &, const UPtrList< T > &)
Write UPtrList to Ostream.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
Dimension set for the base types.
Definition: dimensionSet.H:120
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Definition: fvOption.C:134
A class for handling words, derived from string.
Definition: word.H:59
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
const word & name() const
Return const access to the source name.
Definition: fvOptionI.H:28
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
void checkApplied() const
Check that all sources have been applied.
Definition: fvOptionList.C:73
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
virtual void constrain(fvMatrix< scalar > &eqn, const label fieldi)
Definition: fvOption.C:278
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
virtual bool isActive()
Is the source active?
Definition: fvOption.C:107
virtual label applyToField(const word &fieldName) const
Return index of field name if found in fieldNames list.
Definition: fvOption.C:113
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< fvMatrix< Type > > d2dt2(GeometricField< Type, fvPatchField, volMesh > &field)
Return source for equation with second time derivative.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
A class representing the concept of 1 (scalar(1)) used to avoid unnecessary manipulations for objects...
Definition: one.H:50
Finite volume options abstract base class. Provides a base set of controls, e.g.: ...
Definition: fvOption.H:66