All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
fvPatchField.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-2019 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 "IOobject.H"
27 #include "dictionary.H"
28 #include "fvMesh.H"
29 #include "fvPatchFieldMapper.H"
30 #include "volMesh.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class Type>
36 (
37  const fvPatch& p,
39 )
40 :
41  Field<Type>(p.size()),
42  patch_(p),
43  internalField_(iF),
44  updated_(false),
45  manipulatedMatrix_(false),
46  patchType_(word::null)
47 {}
48 
49 
50 template<class Type>
52 (
53  const fvPatch& p,
55  const Field<Type>& f
56 )
57 :
58  Field<Type>(f),
59  patch_(p),
60  internalField_(iF),
61  updated_(false),
62  manipulatedMatrix_(false),
63  patchType_(word::null)
64 {}
65 
66 
67 template<class Type>
69 (
70  const fvPatch& p,
72  const dictionary& dict,
73  const bool valueRequired
74 )
75 :
76  Field<Type>(p.size()),
77  patch_(p),
78  internalField_(iF),
79  updated_(false),
80  manipulatedMatrix_(false),
81  patchType_(dict.lookupOrDefault<word>("patchType", word::null))
82 {
83  if (valueRequired)
84  {
85  if (dict.found("value"))
86  {
88  (
89  Field<Type>("value", dict, p.size())
90  );
91  }
92  else
93  {
95  (
96  dict
97  ) << "Essential entry 'value' missing"
98  << exit(FatalIOError);
99  }
100  }
101 }
102 
103 
104 template<class Type>
106 (
107  const fvPatchField<Type>& ptf,
108  const fvPatch& p,
110  const fvPatchFieldMapper& mapper,
111  const bool mappingRequired
112 )
113 :
114  Field<Type>(p.size()),
115  patch_(p),
116  internalField_(iF),
117  updated_(false),
118  manipulatedMatrix_(false),
119  patchType_(ptf.patchType_)
120 {
121  if (mappingRequired)
122  {
123  // For unmapped faces set to internal field value (zero-gradient)
124  if (notNull(iF) && mapper.hasUnmapped())
125  {
126  fvPatchField<Type>::operator=(this->patchInternalField());
127  }
128  mapper(*this, ptf);
129  }
130 }
131 
132 
133 template<class Type>
135 (
136  const fvPatchField<Type>& ptf
137 )
138 :
139  Field<Type>(ptf),
140  patch_(ptf.patch_),
141  internalField_(ptf.internalField_),
142  updated_(false),
143  manipulatedMatrix_(false),
144  patchType_(ptf.patchType_)
145 {}
146 
147 
148 template<class Type>
150 (
151  const fvPatchField<Type>& ptf,
153 )
154 :
155  Field<Type>(ptf),
156  patch_(ptf.patch_),
157  internalField_(iF),
158  updated_(false),
159  manipulatedMatrix_(false),
160  patchType_(ptf.patchType_)
161 {}
162 
163 
164 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
165 
166 template<class Type>
168 {
169  return patch_.boundaryMesh().mesh();
170 }
171 
172 
173 template<class Type>
175 {
176  if (&patch_ != &(ptf.patch_))
177  {
179  << "different patches for fvPatchField<Type>s"
180  << abort(FatalError);
181  }
182 }
183 
184 
185 template<class Type>
187 {
188  return patch_.deltaCoeffs()*(*this - patchInternalField());
189 }
190 
191 
192 template<class Type>
195 {
196  return patch_.patchInternalField(internalField_);
197 }
198 
199 
200 template<class Type>
202 {
203  patch_.patchInternalField(internalField_, pif);
204 }
205 
206 
207 template<class Type>
209 (
210  const fvPatchFieldMapper& mapper
211 )
212 {
213  mapper(*this, *this);
214 }
215 
216 
217 template<class Type>
219 (
220  const fvPatchField<Type>& ptf,
221  const labelList& addr
222 )
223 {
224  Field<Type>::rmap(ptf, addr);
225 }
226 
227 
228 template<class Type>
230 {
231  updated_ = true;
232 }
233 
234 
235 template<class Type>
237 {
238  // Default behaviour ignores the weights
239  if (!updated_)
240  {
241  updateCoeffs();
242 
243  updated_ = true;
244  }
245 }
246 
247 
248 template<class Type>
250 {
251  if (!updated_)
252  {
253  updateCoeffs();
254  }
255 
256  updated_ = false;
257  manipulatedMatrix_ = false;
258 }
259 
260 
261 template<class Type>
263 {
264  manipulatedMatrix_ = true;
265 }
266 
267 
268 template<class Type>
270 (
271  fvMatrix<Type>& matrix,
272  const scalarField& weights
273 )
274 {
275  manipulatedMatrix_ = true;
276 }
277 
278 
279 template<class Type>
281 {
282  writeEntry(os, "type", type());
283 
284  if (patchType_.size())
285  {
286  writeEntry(os, "patchType", patchType_);
287  }
288 }
289 
290 
291 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
292 
293 template<class Type>
294 void Foam::fvPatchField<Type>::operator=
295 (
296  const UList<Type>& ul
297 )
298 {
300 }
301 
302 
303 template<class Type>
304 void Foam::fvPatchField<Type>::operator=
305 (
306  const fvPatchField<Type>& ptf
307 )
308 {
309  check(ptf);
311 }
312 
313 
314 template<class Type>
315 void Foam::fvPatchField<Type>::operator+=
316 (
317  const fvPatchField<Type>& ptf
318 )
319 {
320  check(ptf);
322 }
323 
324 
325 template<class Type>
326 void Foam::fvPatchField<Type>::operator-=
327 (
328  const fvPatchField<Type>& ptf
329 )
330 {
331  check(ptf);
333 }
334 
335 
336 template<class Type>
337 void Foam::fvPatchField<Type>::operator*=
338 (
339  const fvPatchField<scalar>& ptf
340 )
341 {
342  if (&patch_ != &ptf.patch())
343  {
345  << "incompatible patches for patch fields"
346  << abort(FatalError);
347  }
348 
350 }
351 
352 
353 template<class Type>
354 void Foam::fvPatchField<Type>::operator/=
355 (
356  const fvPatchField<scalar>& ptf
357 )
358 {
359  if (&patch_ != &ptf.patch())
360  {
362  << abort(FatalError);
363  }
364 
366 }
367 
368 
369 template<class Type>
370 void Foam::fvPatchField<Type>::operator+=
371 (
372  const Field<Type>& tf
373 )
374 {
376 }
377 
378 
379 template<class Type>
380 void Foam::fvPatchField<Type>::operator-=
381 (
382  const Field<Type>& tf
383 )
384 {
386 }
387 
388 
389 template<class Type>
390 void Foam::fvPatchField<Type>::operator*=
391 (
392  const scalarField& tf
393 )
394 {
396 }
397 
398 
399 template<class Type>
400 void Foam::fvPatchField<Type>::operator/=
401 (
402  const scalarField& tf
403 )
404 {
406 }
407 
408 
409 template<class Type>
410 void Foam::fvPatchField<Type>::operator=
411 (
412  const Type& t
413 )
414 {
416 }
417 
418 
419 template<class Type>
420 void Foam::fvPatchField<Type>::operator+=
421 (
422  const Type& t
423 )
424 {
426 }
427 
428 
429 template<class Type>
430 void Foam::fvPatchField<Type>::operator-=
431 (
432  const Type& t
433 )
434 {
436 }
437 
438 
439 template<class Type>
440 void Foam::fvPatchField<Type>::operator*=
441 (
442  const scalar s
443 )
444 {
446 }
447 
448 
449 template<class Type>
450 void Foam::fvPatchField<Type>::operator/=
451 (
452  const scalar s
453 )
454 {
456 }
457 
458 
459 template<class Type>
460 void Foam::fvPatchField<Type>::operator==
461 (
462  const fvPatchField<Type>& ptf
463 )
464 {
466 }
467 
468 
469 template<class Type>
470 void Foam::fvPatchField<Type>::operator==
471 (
472  const Field<Type>& tf
473 )
474 {
476 }
477 
478 
479 template<class Type>
480 void Foam::fvPatchField<Type>::operator==
481 (
482  const Type& t
483 )
484 {
486 }
487 
488 
489 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
490 
491 template<class Type>
492 Foam::Ostream& Foam::operator<<(Ostream& os, const fvPatchField<Type>& ptf)
493 {
494  ptf.write(os);
495 
496  os.check("Ostream& operator<<(Ostream&, const fvPatchField<Type>&");
497 
498  return os;
499 }
500 
501 
502 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
503 
504  #include "fvPatchFieldNew.C"
505 
506 // ************************************************************************* //
dictionary dict
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fvPatchField.C:209
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:186
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
commsTypes
Types of communications.
Definition: UPstream.H:64
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
const tensorField & tf
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:219
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
virtual bool hasUnmapped() const =0
Are there unmapped values? I.e. do all size() elements get.
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
Definition: fvPatchField.C:262
void check(const fvPatchField< Type > &) const
Check fvPatchField<Type> against given fvPatchField<Type>
Definition: fvPatchField.C:174
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:325
virtual label size() const
Return size.
Definition: fvPatch.H:155
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:194
labelList f(nPoints)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:52
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:249
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:229
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:167
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:236
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
fvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Definition: fvPatchField.C:36
IOerror FatalIOError