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-2013 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 word& patchType
56 )
57 :
58  Field<Type>(p.size()),
59  patch_(p),
60  internalField_(iF),
61  updated_(false),
62  manipulatedMatrix_(false),
63  patchType_(patchType)
64 {}
65 
66 
67 template<class Type>
69 (
70  const fvPatch& p,
72  const Field<Type>& f
73 )
74 :
75  Field<Type>(f),
76  patch_(p),
77  internalField_(iF),
78  updated_(false),
79  manipulatedMatrix_(false),
80  patchType_(word::null)
81 {}
82 
83 
84 template<class Type>
86 (
87  const fvPatchField<Type>& ptf,
88  const fvPatch& p,
90  const fvPatchFieldMapper& mapper
91 )
92 :
93  Field<Type>(p.size()),
94  patch_(p),
95  internalField_(iF),
96  updated_(false),
97  manipulatedMatrix_(false),
98  patchType_(ptf.patchType_)
99 {
100  // For unmapped faces set to internal field value (zero-gradient)
101  if (notNull(iF) && iF.size())
102  {
103  fvPatchField<Type>::operator=(this->patchInternalField());
104  }
105  this->map(ptf, mapper);
106 }
107 
108 
109 template<class Type>
111 (
112  const fvPatch& p,
114  const dictionary& dict,
115  const bool valueRequired
116 )
117 :
118  Field<Type>(p.size()),
119  patch_(p),
120  internalField_(iF),
121  updated_(false),
122  manipulatedMatrix_(false),
123  patchType_(dict.lookupOrDefault<word>("patchType", word::null))
124 {
125  if (dict.found("value"))
126  {
128  (
129  Field<Type>("value", dict, p.size())
130  );
131  }
132  else if (!valueRequired)
133  {
135  }
136  else
137  {
139  (
140  "fvPatchField<Type>::fvPatchField"
141  "("
142  "const fvPatch& p,"
143  "const DimensionedField<Type, volMesh>& iF,"
144  "const dictionary& dict,"
145  "const bool valueRequired"
146  ")",
147  dict
148  ) << "Essential entry 'value' missing"
149  << exit(FatalIOError);
150  }
151 }
152 
153 
154 template<class Type>
156 (
157  const fvPatchField<Type>& ptf
158 )
159 :
160  Field<Type>(ptf),
161  patch_(ptf.patch_),
162  internalField_(ptf.internalField_),
163  updated_(false),
164  manipulatedMatrix_(false),
165  patchType_(ptf.patchType_)
166 {}
167 
168 
169 template<class Type>
171 (
172  const fvPatchField<Type>& ptf,
174 )
175 :
176  Field<Type>(ptf),
177  patch_(ptf.patch_),
178  internalField_(iF),
179  updated_(false),
180  manipulatedMatrix_(false),
181  patchType_(ptf.patchType_)
182 {}
183 
184 
185 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
186 
187 template<class Type>
189 {
190  return patch_.boundaryMesh().mesh();
191 }
192 
193 
194 template<class Type>
196 {
197  if (&patch_ != &(ptf.patch_))
198  {
199  FatalErrorIn("PatchField<Type>::check(const fvPatchField<Type>&)")
200  << "different patches for fvPatchField<Type>s"
201  << abort(FatalError);
202  }
203 }
204 
205 
206 template<class Type>
208 {
209  return patch_.deltaCoeffs()*(*this - patchInternalField());
210 }
211 
212 
213 template<class Type>
216 {
217  return patch_.patchInternalField(internalField_);
218 }
219 
220 
221 template<class Type>
223 {
224  patch_.patchInternalField(internalField_, pif);
225 }
226 
227 
228 template<class Type>
230 (
231  const fvPatchFieldMapper& mapper
232 )
233 {
234  Field<Type>& f = *this;
235 
236  if (!this->size())
237  {
238  f.setSize(mapper.size());
239  if (f.size())
240  {
241  f = this->patchInternalField();
242  }
243  }
244  else
245  {
246  // Map all faces provided with mapping data
247  Field<Type>::autoMap(mapper);
248 
249  // For unmapped faces set to internal field value (zero-gradient)
250  if
251  (
252  mapper.direct()
253  && notNull(mapper.directAddressing())
254  && mapper.directAddressing().size()
255  )
256  {
257  Field<Type> pif(this->patchInternalField());
258 
259  const labelList& mapAddressing = mapper.directAddressing();
260 
261  forAll(mapAddressing, i)
262  {
263  if (mapAddressing[i] < 0)
264  {
265  f[i] = pif[i];
266  }
267  }
268  }
269  else if (!mapper.direct() && mapper.addressing().size())
270  {
271  Field<Type> pif(this->patchInternalField());
272 
273  const labelListList& mapAddressing = mapper.addressing();
274 
275  forAll(mapAddressing, i)
276  {
277  const labelList& localAddrs = mapAddressing[i];
278 
279  if (!localAddrs.size())
280  {
281  f[i] = pif[i];
282  }
283  }
284  }
285  }
286 }
287 
288 
289 template<class Type>
291 (
292  const fvPatchField<Type>& ptf,
293  const labelList& addr
294 )
295 {
296  Field<Type>::rmap(ptf, addr);
297 }
298 
299 
300 template<class Type>
302 {
303  updated_ = true;
304 }
305 
306 
307 template<class Type>
309 {
310  if (!updated_)
311  {
312  updateCoeffs();
313 
314  Field<Type>& fld = *this;
315  fld *= weights;
316 
317  updated_ = true;
318  }
319 }
320 
321 
322 template<class Type>
324 {
325  if (!updated_)
326  {
327  updateCoeffs();
328  }
329 
330  updated_ = false;
331  manipulatedMatrix_ = false;
332 }
333 
334 
335 template<class Type>
337 {
338  manipulatedMatrix_ = true;
339 }
340 
341 
342 template<class Type>
344 (
345  fvMatrix<Type>& matrix,
346  const scalarField& weights
347 )
348 {
349  manipulatedMatrix_ = true;
350 }
351 
352 
353 template<class Type>
355 {
356  os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
357 
358  if (patchType_.size())
359  {
360  os.writeKeyword("patchType") << patchType_
361  << token::END_STATEMENT << nl;
362  }
363 }
364 
365 
366 template<class Type>
367 template<class EntryType>
369 (
370  Ostream& os,
371  const word& entryName,
372  const EntryType& value1,
373  const EntryType& value2
374 ) const
375 {
376  if (value1 != value2)
377  {
378  os.writeKeyword(entryName) << value2 << token::END_STATEMENT << nl;
379  }
380 }
381 
382 
383 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
384 
385 template<class Type>
386 void Foam::fvPatchField<Type>::operator=
387 (
388  const UList<Type>& ul
389 )
390 {
392 }
393 
394 
395 template<class Type>
396 void Foam::fvPatchField<Type>::operator=
397 (
398  const fvPatchField<Type>& ptf
399 )
400 {
401  check(ptf);
403 }
404 
405 
406 template<class Type>
407 void Foam::fvPatchField<Type>::operator+=
408 (
409  const fvPatchField<Type>& ptf
410 )
411 {
412  check(ptf);
414 }
415 
416 
417 template<class Type>
418 void Foam::fvPatchField<Type>::operator-=
419 (
420  const fvPatchField<Type>& ptf
421 )
422 {
423  check(ptf);
425 }
426 
427 
428 template<class Type>
429 void Foam::fvPatchField<Type>::operator*=
430 (
431  const fvPatchField<scalar>& ptf
432 )
433 {
434  if (&patch_ != &ptf.patch())
435  {
437  (
438  "PatchField<Type>::operator*=(const fvPatchField<scalar>& ptf)"
439  ) << "incompatible patches for patch fields"
440  << abort(FatalError);
441  }
442 
444 }
445 
446 
447 template<class Type>
448 void Foam::fvPatchField<Type>::operator/=
449 (
450  const fvPatchField<scalar>& ptf
451 )
452 {
453  if (&patch_ != &ptf.patch())
454  {
456  (
457  "PatchField<Type>::operator/=(const fvPatchField<scalar>& ptf)"
458  ) << " incompatible patches for patch fields"
459  << abort(FatalError);
460  }
461 
463 }
464 
465 
466 template<class Type>
467 void Foam::fvPatchField<Type>::operator+=
468 (
469  const Field<Type>& tf
470 )
471 {
473 }
474 
475 
476 template<class Type>
477 void Foam::fvPatchField<Type>::operator-=
478 (
479  const Field<Type>& tf
480 )
481 {
483 }
484 
485 
486 template<class Type>
487 void Foam::fvPatchField<Type>::operator*=
488 (
489  const scalarField& tf
490 )
491 {
493 }
494 
495 
496 template<class Type>
497 void Foam::fvPatchField<Type>::operator/=
498 (
499  const scalarField& tf
500 )
501 {
503 }
504 
505 
506 template<class Type>
507 void Foam::fvPatchField<Type>::operator=
508 (
509  const Type& t
510 )
511 {
513 }
514 
515 
516 template<class Type>
517 void Foam::fvPatchField<Type>::operator+=
518 (
519  const Type& t
520 )
521 {
523 }
524 
525 
526 template<class Type>
527 void Foam::fvPatchField<Type>::operator-=
528 (
529  const Type& t
530 )
531 {
533 }
534 
535 
536 template<class Type>
537 void Foam::fvPatchField<Type>::operator*=
538 (
539  const scalar s
540 )
541 {
543 }
544 
545 
546 template<class Type>
547 void Foam::fvPatchField<Type>::operator/=
548 (
549  const scalar s
550 )
551 {
553 }
554 
555 
556 // Force an assignment, overriding fixedValue status
557 template<class Type>
558 void Foam::fvPatchField<Type>::operator==
559 (
560  const fvPatchField<Type>& ptf
561 )
562 {
564 }
565 
566 
567 template<class Type>
568 void Foam::fvPatchField<Type>::operator==
569 (
570  const Field<Type>& tf
571 )
572 {
574 }
575 
576 
577 template<class Type>
578 void Foam::fvPatchField<Type>::operator==
579 (
580  const Type& t
581 )
582 {
584 }
585 
586 
587 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
588 
589 template<class Type>
590 Foam::Ostream& Foam::operator<<(Ostream& os, const fvPatchField<Type>& ptf)
591 {
592  ptf.write(os);
593 
594  os.check("Ostream& operator<<(Ostream&, const fvPatchField<Type>&");
595 
596  return os;
597 }
598 
599 
600 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
601 
602 # include "fvPatchFieldNew.C"
603 
604 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
Definition: fvPatchField.C:336
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:207
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
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 ))
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:301
labelList f(nPoints)
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:369
A class for handling words, derived from string.
Definition: word.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Foam::fvPatchFieldMapper.
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:68
bool notNull(const T &t)
Return true if t is not a reference to the nullObject of type T.
Definition: nullObjectI.H:46
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:354
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
dictionary dict
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
IOerror FatalIOError
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:291
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:215
commsTypes
Types of communications.
Definition: UPstream.H:64
#define forAll(list, i)
Definition: UList.H:421
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
errorManip< error > abort(error &err)
Definition: errorManip.H:131
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:323
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
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){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
virtual label size() const =0
const tensorField & tf
virtual const labelUList & directAddressing() const
Definition: FieldMapper.H:70
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fvPatchField.C:230
void check(const fvPatchField< Type > &) const
Check fvPatchField<Type> against given fvPatchField<Type>
Definition: fvPatchField.C:195
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Registry of regIOobjects.
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:188
virtual const labelListList & addressing() const
Definition: FieldMapper.H:79
fvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Definition: fvPatchField.C:36
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual bool direct() const =0
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual label size() const
Return size.
Definition: fvPatch.H:161
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:300
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
A class for managing temporary objects.
Definition: PtrList.H:118