pointPatchField.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-2022 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 "pointPatchField.H"
27 #include "pointMesh.H"
28 #include "dictionary.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class Type>
34 (
35  const pointPatch& p,
37 )
38 :
39  patch_(p),
40  internalField_(iF),
41  updated_(false)
42 {}
43 
44 
45 template<class Type>
47 (
48  const pointPatch& p,
50  const dictionary& dict
51 )
52 :
53  patch_(p),
54  internalField_(iF),
55  updated_(false)
56 {}
57 
58 
59 template<class Type>
61 (
62  const pointPatchField<Type>& ptf,
63  const pointPatch& p,
66 )
67 :
68  patch_(p),
69  internalField_(iF),
70  updated_(false)
71 {}
72 
73 
74 template<class Type>
76 (
77  const pointPatchField<Type>& ptf,
79 )
80 :
81  patch_(ptf.patch_),
82  internalField_(iF),
83  updated_(false)
84 {}
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
89 template<class Type>
91 {
92  return patch_.boundaryMesh().mesh()();
93 }
94 
95 
96 template<class Type>
98 {
99  writeEntry(os, "type", type());
100 
101  if (overridesConstraint())
102  {
103  writeEntry(os, "patchType", patch().type());
104  }
105 }
106 
107 
108 template<class Type>
111 {
112  return patchInternalField(primitiveField());
113 }
114 
115 
116 template<class Type>
117 template<class Type1>
120 (
121  const Field<Type1>& iF,
122  const labelList& meshPoints
123 ) const
124 {
125  // Check size
126  if (iF.size() != primitiveField().size())
127  {
129  << "given internal field does not correspond to the mesh. "
130  << "Field size: " << iF.size()
131  << " mesh size: " << primitiveField().size()
132  << abort(FatalError);
133  }
134 
135  return tmp<Field<Type1>>(new Field<Type1>(iF, meshPoints));
136 }
137 
138 
139 template<class Type>
140 template<class Type1>
143 (
144  const Field<Type1>& iF
145 ) const
146 {
147  return patchInternalField(iF, patch().meshPoints());
148 }
149 
150 
151 template<class Type>
152 template<class Type1>
154 (
155  Field<Type1>& iF,
156  const Field<Type1>& pF
157 ) const
158 {
159  // Check size
160  if (iF.size() != primitiveField().size())
161  {
163  << "given internal field does not correspond to the mesh. "
164  << "Field size: " << iF.size()
165  << " mesh size: " << primitiveField().size()
166  << abort(FatalError);
167  }
168 
169  if (pF.size() != size())
170  {
172  << "given patch field does not correspond to the mesh. "
173  << "Field size: " << pF.size()
174  << " mesh size: " << size()
175  << abort(FatalError);
176  }
177 
178  // Get the addressing
179  const labelList& mp = patch().meshPoints();
180 
181  forAll(mp, pointi)
182  {
183  iF[mp[pointi]] += pF[pointi];
184  }
185 }
186 
187 
188 template<class Type>
189 template<class Type1>
191 (
192  Field<Type1>& iF,
193  const Field<Type1>& pF,
194  const labelList& points
195 ) const
196 {
197  // Check size
198  if (iF.size() != primitiveField().size())
199  {
201  << "given internal field does not correspond to the mesh. "
202  << "Field size: " << iF.size()
203  << " mesh size: " << primitiveField().size()
204  << abort(FatalError);
205  }
206 
207  if (pF.size() != size())
208  {
210  << "given patch field does not correspond to the mesh. "
211  << "Field size: " << pF.size()
212  << " mesh size: " << size()
213  << abort(FatalError);
214  }
215 
216  // Get the addressing
217  const labelList& mp = patch().meshPoints();
218 
219  forAll(points, i)
220  {
221  label pointi = points[i];
222  iF[mp[pointi]] += pF[pointi];
223  }
224 }
225 
226 
227 template<class Type>
228 template<class Type1>
230 (
231  Field<Type1>& iF,
232  const Field<Type1>& pF,
233  const labelList& meshPoints
234 ) const
235 {
236  // Check size
237  if (iF.size() != primitiveField().size())
238  {
240  << "given internal field does not correspond to the mesh. "
241  << "Field size: " << iF.size()
242  << " mesh size: " << primitiveField().size()
243  << abort(FatalError);
244  }
245 
246  if (pF.size() != meshPoints.size())
247  {
249  << "given patch field does not correspond to the meshPoints. "
250  << "Field size: " << pF.size()
251  << " meshPoints size: " << size()
252  << abort(FatalError);
253  }
254 
255  forAll(meshPoints, pointi)
256  {
257  iF[meshPoints[pointi]] = pF[pointi];
258  }
259 }
260 
261 
262 template<class Type>
263 template<class Type1>
265 (
266  Field<Type1>& iF,
267  const Field<Type1>& pF
268 ) const
269 {
270  setInternalField(iF, pF, patch().meshPoints());
271 }
272 
273 
274 template<class Type>
276 {
277  if (!updated_)
278  {
279  updateCoeffs();
280  }
281 
282  updated_ = false;
283 }
284 
285 
286 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
287 
288 template<class Type>
289 Foam::Ostream& Foam::operator<<
290 (
291  Ostream& os,
292  const pointPatchField<Type>& ptf
293 )
294 {
295  ptf.write(os);
296 
297  os.check("Ostream& operator<<(Ostream&, const pointPatchField<Type>&)");
298 
299  return os;
300 }
301 
302 
303 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
304 
305 #include "pointPatchFieldNew.C"
306 
307 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Pre-declare SubField and related Field type.
Definition: Field.H:82
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
virtual Ostream & write(const char)=0
Write character.
commsTypes
Types of communications.
Definition: UPstream.H:65
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
Registry of regIOobjects.
Foam::pointPatchFieldMapper.
Abstract base class for point-mesh patch fields.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking)
Evaluate the patch field.
virtual void write(Ostream &) const
Write.
tmp< Field< Type > > patchInternalField() const
Return field created from appropriate internal field values.
const objectRegistry & db() const
Return local objectRegistry.
void addToInternalField(Field< Type1 > &iF, const Field< Type1 > &pF) const
Given the internal field and a patch field,.
void setInternalField(Field< Type1 > &iF, const Field< Type1 > &pF, const labelList &meshPoints) const
Given the internal field and a patch field,.
pointPatchField(const pointPatch &, const DimensionedField< Type, pointMesh > &)
Construct from patch and internal field.
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:57
A class for managing temporary objects.
Definition: tmp.H:55
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
const pointField & points
const dimensionedScalar mp
Proton mass.
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
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
error FatalError
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
volScalarField & p