fvPatchField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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 Type& value
56 )
57 :
58  Field<Type>(p.size(), value),
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 word& patchType
73 )
74 :
75  Field<Type>(p.size()),
76  patch_(p),
77  internalField_(iF),
78  updated_(false),
79  manipulatedMatrix_(false),
80  patchType_(patchType)
81 {}
82 
83 
84 template<class Type>
86 (
87  const fvPatch& p,
89  const Field<Type>& f
90 )
91 :
92  Field<Type>(f),
93  patch_(p),
94  internalField_(iF),
95  updated_(false),
96  manipulatedMatrix_(false),
97  patchType_(word::null)
98 {}
99 
100 
101 template<class Type>
103 (
104  const fvPatch& p,
106  const dictionary& dict,
107  const bool valueRequired
108 )
109 :
110  Field<Type>(p.size()),
111  patch_(p),
112  internalField_(iF),
113  updated_(false),
114  manipulatedMatrix_(false),
115  patchType_(dict.lookupOrDefault<word>("patchType", word::null))
116 {
117  if (valueRequired)
118  {
119  if (dict.found("value"))
120  {
122  (
123  Field<Type>("value", dict, p.size())
124  );
125  }
126  else
127  {
129  (
130  dict
131  ) << "Essential entry 'value' missing"
132  << exit(FatalIOError);
133  }
134  }
135 }
136 
137 
138 template<class Type>
140 (
141  const fvPatchField<Type>& ptf,
142  const fvPatch& p,
144  const fvPatchFieldMapper& mapper
145 )
146 :
147  Field<Type>(p.size()),
148  patch_(p),
149  internalField_(iF),
150  updated_(false),
151  manipulatedMatrix_(false),
152  patchType_(ptf.patchType_)
153 {
154  // For unmapped faces set to internal field value (zero-gradient)
155  if (notNull(iF) && mapper.hasUnmapped())
156  {
157  fvPatchField<Type>::operator=(this->patchInternalField());
158  }
159  this->map(ptf, mapper);
160 }
161 
162 
163 template<class Type>
165 (
166  const fvPatchField<Type>& ptf
167 )
168 :
169  Field<Type>(ptf),
170  patch_(ptf.patch_),
171  internalField_(ptf.internalField_),
172  updated_(false),
173  manipulatedMatrix_(false),
174  patchType_(ptf.patchType_)
175 {}
176 
177 
178 template<class Type>
180 (
181  const fvPatchField<Type>& ptf,
183 )
184 :
185  Field<Type>(ptf),
186  patch_(ptf.patch_),
187  internalField_(iF),
188  updated_(false),
189  manipulatedMatrix_(false),
190  patchType_(ptf.patchType_)
191 {}
192 
193 
194 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
195 
196 template<class Type>
198 {
199  return patch_.boundaryMesh().mesh();
200 }
201 
202 
203 template<class Type>
205 {
206  if (&patch_ != &(ptf.patch_))
207  {
209  << "different patches for fvPatchField<Type>s"
210  << abort(FatalError);
211  }
212 }
213 
214 
215 template<class Type>
217 {
218  return patch_.deltaCoeffs()*(*this - patchInternalField());
219 }
220 
221 
222 template<class Type>
225 {
226  return patch_.patchInternalField(internalField_);
227 }
228 
229 
230 template<class Type>
232 {
233  patch_.patchInternalField(internalField_, pif);
234 }
235 
236 
237 template<class Type>
239 (
240  const fvPatchFieldMapper& mapper
241 )
242 {
243  Field<Type>& f = *this;
244 
245  if (!this->size() && !mapper.distributed())
246  {
247  f.setSize(mapper.size());
248  if (f.size())
249  {
250  f = this->patchInternalField();
251  }
252  }
253  else
254  {
255  // Map all faces provided with mapping data
256  Field<Type>::autoMap(mapper);
257 
258  // For unmapped faces set to internal field value (zero-gradient)
259  if (mapper.hasUnmapped())
260  {
261  Field<Type> pif(this->patchInternalField());
262 
263  if
264  (
265  mapper.direct()
266  && notNull(mapper.directAddressing())
267  && mapper.directAddressing().size()
268  )
269  {
270  const labelList& mapAddressing = mapper.directAddressing();
271 
272  forAll(mapAddressing, i)
273  {
274  if (mapAddressing[i] < 0)
275  {
276  f[i] = pif[i];
277  }
278  }
279  }
280  else if (!mapper.direct() && mapper.addressing().size())
281  {
282  const labelListList& mapAddressing = mapper.addressing();
283 
284  forAll(mapAddressing, i)
285  {
286  const labelList& localAddrs = mapAddressing[i];
287 
288  if (!localAddrs.size())
289  {
290  f[i] = pif[i];
291  }
292  }
293  }
294  }
295  }
296 }
297 
298 
299 template<class Type>
301 (
302  const fvPatchField<Type>& ptf,
303  const labelList& addr
304 )
305 {
306  Field<Type>::rmap(ptf, addr);
307 }
308 
309 
310 template<class Type>
312 {
313  updated_ = true;
314 }
315 
316 
317 template<class Type>
319 {
320  // Default behaviour ignores the weights
321  if (!updated_)
322  {
323  updateCoeffs();
324 
325  updated_ = true;
326  }
327 }
328 
329 
330 template<class Type>
332 {
333  if (!updated_)
334  {
335  updateCoeffs();
336  }
337 
338  updated_ = false;
339  manipulatedMatrix_ = false;
340 }
341 
342 
343 template<class Type>
345 {
346  manipulatedMatrix_ = true;
347 }
348 
349 
350 template<class Type>
352 (
353  fvMatrix<Type>& matrix,
354  const scalarField& weights
355 )
356 {
357  manipulatedMatrix_ = true;
358 }
359 
360 
361 template<class Type>
363 {
364  os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
365 
366  if (patchType_.size())
367  {
368  os.writeKeyword("patchType") << patchType_
369  << token::END_STATEMENT << nl;
370  }
371 }
372 
373 
374 template<class Type>
375 template<class EntryType>
377 (
378  Ostream& os,
379  const word& entryName,
380  const EntryType& value1,
381  const EntryType& value2
382 ) const
383 {
384  if (value1 != value2)
385  {
386  os.writeKeyword(entryName) << value2 << token::END_STATEMENT << nl;
387  }
388 }
389 
390 
391 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
392 
393 template<class Type>
394 void Foam::fvPatchField<Type>::operator=
395 (
396  const UList<Type>& ul
397 )
398 {
400 }
401 
402 
403 template<class Type>
404 void Foam::fvPatchField<Type>::operator=
405 (
406  const fvPatchField<Type>& ptf
407 )
408 {
409  check(ptf);
411 }
412 
413 
414 template<class Type>
415 void Foam::fvPatchField<Type>::operator+=
416 (
417  const fvPatchField<Type>& ptf
418 )
419 {
420  check(ptf);
422 }
423 
424 
425 template<class Type>
426 void Foam::fvPatchField<Type>::operator-=
427 (
428  const fvPatchField<Type>& ptf
429 )
430 {
431  check(ptf);
433 }
434 
435 
436 template<class Type>
437 void Foam::fvPatchField<Type>::operator*=
438 (
439  const fvPatchField<scalar>& ptf
440 )
441 {
442  if (&patch_ != &ptf.patch())
443  {
445  << "incompatible patches for patch fields"
446  << abort(FatalError);
447  }
448 
450 }
451 
452 
453 template<class Type>
454 void Foam::fvPatchField<Type>::operator/=
455 (
456  const fvPatchField<scalar>& ptf
457 )
458 {
459  if (&patch_ != &ptf.patch())
460  {
462  << abort(FatalError);
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 Field<Type>& tf
483 )
484 {
486 }
487 
488 
489 template<class Type>
490 void Foam::fvPatchField<Type>::operator*=
491 (
492  const scalarField& tf
493 )
494 {
496 }
497 
498 
499 template<class Type>
500 void Foam::fvPatchField<Type>::operator/=
501 (
502  const scalarField& tf
503 )
504 {
506 }
507 
508 
509 template<class Type>
510 void Foam::fvPatchField<Type>::operator=
511 (
512  const Type& t
513 )
514 {
516 }
517 
518 
519 template<class Type>
520 void Foam::fvPatchField<Type>::operator+=
521 (
522  const Type& t
523 )
524 {
526 }
527 
528 
529 template<class Type>
530 void Foam::fvPatchField<Type>::operator-=
531 (
532  const Type& t
533 )
534 {
536 }
537 
538 
539 template<class Type>
540 void Foam::fvPatchField<Type>::operator*=
541 (
542  const scalar s
543 )
544 {
546 }
547 
548 
549 template<class Type>
550 void Foam::fvPatchField<Type>::operator/=
551 (
552  const scalar s
553 )
554 {
556 }
557 
558 
559 template<class Type>
560 void Foam::fvPatchField<Type>::operator==
561 (
562  const fvPatchField<Type>& ptf
563 )
564 {
566 }
567 
568 
569 template<class Type>
570 void Foam::fvPatchField<Type>::operator==
571 (
572  const Field<Type>& tf
573 )
574 {
576 }
577 
578 
579 template<class Type>
580 void Foam::fvPatchField<Type>::operator==
581 (
582  const Type& t
583 )
584 {
586 }
587 
588 
589 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
590 
591 template<class Type>
592 Foam::Ostream& Foam::operator<<(Ostream& os, const fvPatchField<Type>& ptf)
593 {
594  ptf.write(os);
595 
596  os.check("Ostream& operator<<(Ostream&, const fvPatchField<Type>&");
597 
598  return os;
599 }
600 
601 
602 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
603 
604  #include "fvPatchFieldNew.C"
605 
606 // ************************************************************************* //
dictionary dict
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:431
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
type
Types of root.
Definition: Roots.H:52
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fvPatchField.C:239
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:216
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:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
commsTypes
Types of communications.
Definition: UPstream.H:64
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
virtual const labelListList & addressing() const
Definition: FieldMapper.H:94
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 bool direct() const =0
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
virtual const labelUList & directAddressing() const
Definition: FieldMapper.H:85
const tensorField & tf
virtual label size() const =0
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:57
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:301
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 void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
Definition: fvPatchField.C:344
void check(const fvPatchField< Type > &) const
Check fvPatchField<Type> against given fvPatchField<Type>
Definition: fvPatchField.C:204
errorManip< error > abort(error &err)
Definition: errorManip.H:131
virtual bool hasUnmapped() const =0
Are there unmapped values? I.e. do all size() elements get.
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:340
virtual label size() const
Return size.
Definition: fvPatch.H:161
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual bool distributed() const
Definition: FieldMapper.H:68
static const char nl
Definition: Ostream.H:262
void writeEntryIfDifferent(Ostream &os, const word &entryName, const EntryType &value1, const EntryType &value2) const
Helper function to write the keyword and entry only if the.
Definition: fvPatchField.C:377
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:224
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
labelList f(nPoints)
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:46
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:331
void setSize(const label)
Reset size of List.
Definition: List.C:281
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:311
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:197
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:318
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
fvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Definition: fvPatchField.C:36
IOerror FatalIOError