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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2015-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "IOobject.H"
30#include "dictionary.H"
31#include "fvMesh.H"
32#include "fvPatchFieldMapper.H"
33#include "volMesh.H"
34
35// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36
37template<class Type>
39(
40 const fvPatch& p,
42)
43:
44 Field<Type>(p.size()),
45 patch_(p),
46 internalField_(iF),
47 updated_(false),
48 manipulatedMatrix_(false),
49 useImplicit_(false),
50 patchType_()
51{}
52
53
54template<class Type>
56(
57 const fvPatch& p,
59 const Type& value
60)
61:
62 Field<Type>(p.size(), value),
63 patch_(p),
64 internalField_(iF),
65 updated_(false),
66 manipulatedMatrix_(false),
67 useImplicit_(false),
68 patchType_()
69{}
70
71
72template<class Type>
74(
75 const fvPatch& p,
77 const word& patchType
78)
79:
80 Field<Type>(p.size()),
81 patch_(p),
82 internalField_(iF),
83 updated_(false),
84 manipulatedMatrix_(false),
85 useImplicit_(false),
86 patchType_(patchType)
87{}
88
89
90template<class Type>
92(
93 const fvPatch& p,
95 const Field<Type>& f
96)
97:
98 Field<Type>(f),
99 patch_(p),
100 internalField_(iF),
101 updated_(false),
102 manipulatedMatrix_(false),
103 useImplicit_(false),
104 patchType_()
105{}
106
107
108template<class Type>
110(
111 const fvPatch& p,
113 const dictionary& dict,
114 const bool valueRequired
115)
116:
117 Field<Type>(p.size()),
118 patch_(p),
119 internalField_(iF),
120 updated_(false),
121 manipulatedMatrix_(false),
122 useImplicit_(false),
123 patchType_()
124{
125 dict.readIfPresent("useImplicit", useImplicit_, keyType::LITERAL);
126 dict.readIfPresent("patchType", patchType_, keyType::LITERAL);
127
128 if (valueRequired)
129 {
130 if (dict.found("value"))
131 {
133 (
134 Field<Type>("value", dict, p.size())
135 );
136 }
137 else
138 {
140 << "Essential entry 'value' missing on patch "
141 << p.name() << endl
142 << exit(FatalIOError);
143 }
144 }
145}
146
147
148template<class Type>
150(
151 const fvPatchField<Type>& ptf,
152 const fvPatch& p,
154 const fvPatchFieldMapper& mapper
155)
156:
157 Field<Type>(p.size()),
158 patch_(p),
159 internalField_(iF),
160 updated_(false),
161 manipulatedMatrix_(false),
162 useImplicit_(ptf.useImplicit_),
163 patchType_(ptf.patchType_)
164{
165 // For unmapped faces set to internal field value (zero-gradient)
166 if (notNull(iF) && mapper.hasUnmapped())
167 {
169 }
170 this->map(ptf, mapper);
171}
173
174template<class Type>
176(
177 const fvPatchField<Type>& ptf
178)
180 Field<Type>(ptf),
181 patch_(ptf.patch_),
182 internalField_(ptf.internalField_),
183 updated_(false),
184 manipulatedMatrix_(false),
185 useImplicit_(ptf.useImplicit_),
186 patchType_(ptf.patchType_)
188
189
190template<class Type>
192(
193 const fvPatchField<Type>& ptf,
196:
197 Field<Type>(ptf),
198 patch_(ptf.patch_),
199 internalField_(iF),
200 updated_(false),
201 manipulatedMatrix_(false),
202 useImplicit_(ptf.useImplicit_),
203 patchType_(ptf.patchType_)
204{}
205
206
207// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
208
209template<class Type>
211{
212 return patch_.boundaryMesh().mesh();
213}
214
215
216template<class Type>
218{
219 if (&patch_ != &(ptf.patch_))
220 {
222 << "different patches for fvPatchField<Type>s"
223 << abort(FatalError);
224 }
225}
226
227
228template<class Type>
231 return patch_.deltaCoeffs()*(*this - patchInternalField());
232}
233
234
235template<class Type>
238{
239 return patch_.patchInternalField(internalField_);
240}
241
242
243template<class Type>
245{
246 patch_.patchInternalField(internalField_, pif);
247}
248
249
250template<class Type>
252(
253 const fvPatchFieldMapper& mapper
254)
255{
256 Field<Type>& f = *this;
257
258 if (!this->size() && !mapper.distributed())
259 {
260 f.setSize(mapper.size());
261 if (f.size())
262 {
263 f = this->patchInternalField();
264 }
265 }
266 else
267 {
268 // Map all faces provided with mapping data
269 Field<Type>::autoMap(mapper);
270
271
272 // For unmapped faces set to internal field value (zero-gradient)
273 if (mapper.hasUnmapped())
274 {
275 Field<Type> pif(this->patchInternalField());
276
277 if
278 (
279 mapper.direct()
280 && notNull(mapper.directAddressing())
281 && mapper.directAddressing().size()
282 )
283 {
284 const labelList& mapAddressing = mapper.directAddressing();
285
286 forAll(mapAddressing, i)
287 {
288 if (mapAddressing[i] < 0)
289 {
290 f[i] = pif[i];
291 }
292 }
293 }
294 else if (!mapper.direct() && mapper.addressing().size())
295 {
296 const labelListList& mapAddressing = mapper.addressing();
297
298 forAll(mapAddressing, i)
299 {
300 const labelList& localAddrs = mapAddressing[i];
301
302 if (!localAddrs.size())
303 {
304 f[i] = pif[i];
305 }
306 }
307 }
308 }
309 }
310}
311
312
313template<class Type>
315(
316 const fvPatchField<Type>& ptf,
317 const labelList& addr
318)
319{
320 Field<Type>::rmap(ptf, addr);
321}
322
323
324template<class Type>
326{
327 updated_ = true;
328}
329
330
331template<class Type>
333{
334 // Default behaviour ignores the weights
335 if (!updated_)
336 {
337 updateCoeffs();
338
339 updated_ = true;
340 }
341}
342
343
344template<class Type>
346{
347 if (!updated_)
348 {
349 updateCoeffs();
350 }
351
352 updated_ = false;
353 manipulatedMatrix_ = false;
354}
355
356
357template<class Type>
360 manipulatedMatrix_ = true;
361}
362
363
364template<class Type>
366(
367 fvMatrix<Type>& matrix,
368 const scalarField& weights
369)
370{
371 manipulatedMatrix_ = true;
372}
373
374
375template<class Type>
377(
378 fvMatrix<Type>& matrix,
379 const label iMatrix,
380 const direction cmp
381)
382{
383 manipulatedMatrix_ = true;
384}
385
386
387template<class Type>
389{
390 os.writeEntry("type", type());
391
392 if (useImplicit_)
393 {
394 os.writeEntry("useImplicit", "true");
395 }
396
397 if (!patchType_.empty())
398 {
399 os.writeEntry("patchType", patchType_);
400 }
401}
402
403
404// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
405
406template<class Type>
408(
409 const UList<Type>& ul
410)
411{
414
415
416template<class Type>
418(
419 const fvPatchField<Type>& ptf
420)
421{
422 check(ptf);
424}
425
426
427template<class Type>
429(
430 const fvPatchField<Type>& ptf
431)
432{
433 check(ptf);
435}
436
437
438template<class Type>
440(
441 const fvPatchField<Type>& ptf
442)
443{
444 check(ptf);
446}
448
449template<class Type>
451(
452 const fvPatchField<scalar>& ptf
453)
454{
455 if (&patch_ != &ptf.patch())
456 {
458 << "incompatible patches for patch fields"
459 << abort(FatalError);
460 }
461
463}
464
465
466template<class Type>
469 const fvPatchField<scalar>& ptf
470)
471{
472 if (&patch_ != &ptf.patch())
473 {
475 << abort(FatalError);
476 }
477
479}
480
481
482template<class Type>
484(
485 const Field<Type>& tf
486)
487{
489}
490
491
492template<class Type>
494(
495 const Field<Type>& tf
496)
497{
499}
500
501
502template<class Type>
504(
505 const scalarField& tf
506)
507{
509}
510
511
512template<class Type>
514(
515 const scalarField& tf
516)
517{
519}
520
521
522template<class Type>
524(
525 const Type& t
526)
527{
529}
530
531
532template<class Type>
534(
535 const Type& t
536)
537{
540
541
542template<class Type>
544(
545 const Type& t
546)
547{
550
551
552template<class Type>
554(
555 const scalar s
556)
557{
559}
561
562template<class Type>
564(
565 const scalar s
567{
569}
570
572template<class Type>
581
582template<class Type>
584(
585 const Field<Type>& tf
590
591
592template<class Type>
595 const Type& t
597{
599}
600
601
602// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
603
604template<class Type>
606{
607 ptf.write(os);
608
610
611 return os;
612}
613
614
615// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
616
617#include "fvPatchFieldNew.C"
618
619// ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual bool direct() const =0
Is it a direct (non-interpolating) mapper?
virtual const labelUList & directAddressing() const
Return the direct addressing values.
Definition: FieldMapper.H:82
virtual label size() const =0
The size of the mapper.
virtual bool distributed() const
Does the mapper have remote contributions?
Definition: FieldMapper.H:72
virtual const labelListList & addressing() const
Return the interpolation addressing.
Definition: FieldMapper.H:102
virtual bool hasUnmapped() const =0
Any unmapped values?
Generic templated field type.
Definition: Field.H:82
void operator=(const Field< Type > &)
Copy assignment.
Definition: Field.C:641
void autoMap(const FieldMapper &map, const bool applyFlip=true)
Map from self.
Definition: Field.C:403
void operator+=(const UList< Type > &)
Definition: Field.C:693
void operator-=(const UList< Type > &)
Definition: Field.C:694
void operator*=(const UList< scalar > &)
Definition: Field.C:695
void operator/=(const UList< scalar > &)
Definition: Field.C:696
void map(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 map from the given field
Definition: Field.C:240
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:466
void evaluate()
Evaluate boundary conditions.
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:58
void setSize(const label n)
Alias for resize()
Definition: List.H:218
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:94
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
commsTypes
Types of communications.
Definition: UPstream.H:67
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionaryI.H:87
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
virtual bool write()
Write the output fields.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:82
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fvPatchField.C:252
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:388
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:237
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:229
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:210
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:408
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:325
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
Definition: fvPatchField.C:358
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:332
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:315
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:362
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
@ LITERAL
String literal.
Definition: keyType.H:81
Registry of regIOobjects.
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.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))
#define FUNCTION_NAME
static void check(const int retVal, const char *what)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
errorManip< error > abort(error &err)
Definition: errorManip.H:144
uint8_t direction
Definition: direction.H:56
bool notNull(const T *ptr)
True if ptr is not a pointer (of type T) to the nullObject.
Definition: nullObject.H:207
IOerror FatalIOError
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
labelList f(nPoints)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333