sampledSurface.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 "sampledSurface.H"
27 #include "polyMesh.H"
28 #include "demandDrivenData.H"
29 
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(sampledSurface, 0);
36  defineRunTimeSelectionTable(sampledSurface, word);
37 }
38 
39 
40 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41 
43 {
44  deleteDemandDrivenData(SfPtr_);
45  deleteDemandDrivenData(magSfPtr_);
46  deleteDemandDrivenData(CfPtr_);
47  area_ = -1;
48 }
49 
50 
51 void Foam::sampledSurface::makeSf() const
52 {
53  // It is an error to recalculate if the pointer is already set
54  if (SfPtr_)
55  {
56  FatalErrorIn("Foam::sampledSurface::makeSf()")
57  << "face area vectors already exist"
58  << abort(FatalError);
59  }
60 
61  const faceList& theFaces = faces();
62  SfPtr_ = new vectorField(theFaces.size());
63 
64  vectorField& values = *SfPtr_;
65  forAll(theFaces, faceI)
66  {
67  values[faceI] = theFaces[faceI].normal(points());
68  }
69 }
70 
71 
72 void Foam::sampledSurface::makeMagSf() const
73 {
74  // It is an error to recalculate if the pointer is already set
75  if (magSfPtr_)
76  {
77  FatalErrorIn("Foam::sampledSurface::makeMagSf()")
78  << "mag face areas already exist"
79  << abort(FatalError);
80  }
81 
82  const faceList& theFaces = faces();
83  magSfPtr_ = new scalarField(theFaces.size());
84 
85  scalarField& values = *magSfPtr_;
86  forAll(theFaces, faceI)
87  {
88  values[faceI] = theFaces[faceI].mag(points());
89  }
90 }
91 
92 
93 void Foam::sampledSurface::makeCf() const
94 {
95  // It is an error to recalculate if the pointer is already set
96  if (CfPtr_)
97  {
98  FatalErrorIn("Foam::sampledSurface::makeCf()")
99  << "face centres already exist"
100  << abort(FatalError);
101  }
102 
103  const faceList& theFaces = faces();
104  CfPtr_ = new vectorField(theFaces.size());
105 
106  vectorField& values = *CfPtr_;
107  forAll(theFaces, faceI)
108  {
109  values[faceI] = theFaces[faceI].centre(points());
110  }
111 }
112 
113 
114 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
115 
117 (
118  const word& name,
119  const polyMesh& mesh,
120  const dictionary& dict
121 )
122 {
123  const word sampleType(dict.lookup("type"));
124 
125  if (debug)
126  {
127  Info<< "Selecting sampledType " << sampleType << endl;
128  }
129 
130  wordConstructorTable::iterator cstrIter =
131  wordConstructorTablePtr_->find(sampleType);
132 
133  if (cstrIter == wordConstructorTablePtr_->end())
134  {
136  (
137  "sampledSurface::New"
138  "(const word&, const polyMesh&, const dictionary&)"
139  ) << "Unknown sample type "
140  << sampleType << nl << nl
141  << "Valid sample types : " << endl
142  << wordConstructorTablePtr_->sortedToc()
143  << exit(FatalError);
144  }
145 
146  return autoPtr<sampledSurface>(cstrIter()(name, mesh, dict));
147 }
148 
149 
150 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
151 
153 (
154  const word& name,
155  const polyMesh& mesh,
156  const bool interpolate
157 )
158 :
159  name_(name),
160  mesh_(mesh),
161  interpolate_(interpolate),
162  SfPtr_(NULL),
163  magSfPtr_(NULL),
164  CfPtr_(NULL),
165  area_(-1)
166 {}
167 
168 
170 (
171  const word& name,
172  const polyMesh& mesh,
173  const dictionary& dict
174 )
175 :
176  name_(name),
177  mesh_(mesh),
178  interpolate_(dict.lookupOrDefault("interpolate", false)),
179  SfPtr_(NULL),
180  magSfPtr_(NULL),
181  CfPtr_(NULL),
182  area_(-1)
183 {
184  dict.readIfPresent("name", name_);
185 }
186 
187 
188 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
189 
191 {
192  clearGeom();
193 }
194 
195 
196 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
197 
199 {
200  if (!SfPtr_)
201  {
202  makeSf();
203  }
204 
205  return *SfPtr_;
206 }
207 
208 
210 {
211  if (!magSfPtr_)
212  {
213  makeMagSf();
214  }
215 
216  return *magSfPtr_;
217 }
218 
219 
221 {
222  if (!CfPtr_)
223  {
224  makeCf();
225  }
226 
227  return *CfPtr_;
228 }
229 
230 
231 Foam::scalar Foam::sampledSurface::area() const
232 {
233  if (area_ < 0)
234  {
235  area_ = sum(magSf());
236  reduce(area_, sumOp<scalar>());
237  }
238 
239  return area_;
240 }
241 
242 
244 (
245  const surfaceScalarField& sField
246 ) const
247 {
248  notImplemented("tmp<Foam::scalarField> sampledSurface::sample");
249  return tmp<scalarField>(NULL);
250 }
251 
252 
254 (
255  const surfaceVectorField& sField
256 ) const
257 {
258  notImplemented("tmp<Foam::vectorField> sampledSurface::sample");
259  return tmp<vectorField>(NULL);
260 }
261 
262 
264 (
265  const surfaceSphericalTensorField& sField
266 ) const
267 {
268  notImplemented("tmp<Foam::sphericalTensorField> sampledSurface::sample");
269  return tmp<sphericalTensorField>(NULL);
270 }
271 
272 
274 (
275  const surfaceSymmTensorField& sField
276 ) const
277 {
278  notImplemented("tmp<Foam::symmTensorField> sampledSurface::sample");
279  return tmp<symmTensorField>(NULL);
280 }
281 
282 
284 (
285  const surfaceTensorField& sField
286 ) const
287 {
288  notImplemented("tmp<Foam::tensorField> sampledSurface::sample");
289  return tmp<tensorField>(NULL);
290 }
291 
292 
294 Foam::sampledSurface::project(const Field<scalar>& field) const
295 {
296  tmp<Field<scalar> > tRes(new Field<scalar>(faces().size()));
297  Field<scalar>& res = tRes();
298 
299  forAll(faces(), faceI)
300  {
301  res[faceI] = field[faceI];
302  }
303 
304  return tRes;
305 }
306 
307 
309 Foam::sampledSurface::project(const Field<vector>& field) const
310 {
311  tmp<Field<scalar> > tRes(new Field<scalar>(faces().size()));
312  project(tRes(), field);
313  return tRes;
314 }
315 
316 
318 Foam::sampledSurface::project(const Field<sphericalTensor>& field) const
319 {
320  tmp<Field<vector> > tRes(new Field<vector>(faces().size()));
321  project(tRes(), field);
322  return tRes;
323 }
324 
325 
327 Foam::sampledSurface::project(const Field<symmTensor>& field) const
328 {
329  tmp<Field<vector> > tRes(new Field<vector>(faces().size()));
330  project(tRes(), field);
331  return tRes;
332 }
333 
334 
336 Foam::sampledSurface::project(const Field<tensor>& field) const
337 {
338  tmp<Field<vector> > tRes(new Field<vector>(faces().size()));
339  project(tRes(), field);
340  return tRes;
341 }
342 
343 
345 {
346  os << type();
347 }
348 
349 
350 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
351 
353 {
354  s.print(os);
355  os.check("Ostream& operator<<(Ostream&, const sampledSurface&");
356  return os;
357 }
358 
359 
360 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual tmp< scalarField > sample(const volScalarField &) const =0
Sample field on surface.
const pointField & points
virtual const vectorField & Cf() const
Return face centres as vectorField.
An abstract class for surfaces with sampling.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void deleteDemandDrivenData(DataPtr &dataPtr)
A class for handling words, derived from string.
Definition: word.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
sampledSurface(const word &name, const polyMesh &, const bool interpolate=false)
Construct from name, mesh.
virtual const scalarField & magSf() const
Return face area magnitudes.
messageStream Info
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Ostream & operator<<(Ostream &, const edgeMesh &)
Definition: edgeMeshIO.C:133
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
dictionary dict
virtual const vectorField & Sf() const
Return face area vectors.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Field< vector > vectorField
Specialisation of Field<T> for vector.
#define forAll(list, i)
Definition: UList.H:421
virtual void clearGeom() const
Template functions to aid in the implementation of demand driven data.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
scalar area() const
The total surface area.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static autoPtr< sampledSurface > New(const word &name, const polyMesh &, const dictionary &)
Return a reference to the selected surface.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
virtual void print(Ostream &) const
Write.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
#define notImplemented(functionName)
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:356
A class for managing temporary objects.
Definition: PtrList.H:118
virtual ~sampledSurface()
Destructor.
defineTypeNameAndDebug(combustionModel, 0)