FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/polyMesh.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 
116 
119 (
120  const word& name,
121  const polyMesh& mesh,
122  const dictionary& dict
123 )
124 {
125  word sampleType(dict.lookup("type"));
126  if (debug)
127  {
128  Info<< "Selecting sampledType " << sampleType << endl;
129  }
130 
131  wordConstructorTable::iterator cstrIter =
132  wordConstructorTablePtr_->find(sampleType);
133 
134  if (cstrIter == wordConstructorTablePtr_->end())
135  {
137  (
138  "sampledSurface::New"
139  "(const word&, const polyMesh&, const dictionary&)"
140  ) << "Unknown sample type " << sampleType
141  << endl << endl
142  << "Valid sample types : " << endl
143  << wordConstructorTablePtr_->sortedToc()
144  << exit(FatalError);
145  }
146 
147  return autoPtr<sampledSurface>(cstrIter()(name, mesh, dict));
148 }
149 
150 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
151 
153 (
154  const word& name,
155  const polyMesh& mesh
156 )
157 :
158  name_(name),
159  mesh_(mesh),
160  interpolate_(false),
161  SfPtr_(NULL),
162  magSfPtr_(NULL),
163  CfPtr_(NULL),
164  area_(-1)
165 {}
166 
167 
168 // Construct from dictionary
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 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
196 
198 {
199  if (!SfPtr_)
200  {
201  makeSf();
202  }
203 
204  return *SfPtr_;
205 }
206 
207 
209 {
210  if (!magSfPtr_)
211  {
212  makeMagSf();
213  }
214 
215  return *magSfPtr_;
216 }
217 
218 
220 {
221  if (!CfPtr_)
222  {
223  makeCf();
224  }
225 
226  return *CfPtr_;
227 }
228 
229 
230 Foam::scalar Foam::sampledSurface::area() const
231 {
232  if (area_ < 0)
233  {
234  area_ = sum(magSf());
235  reduce(area_, sumOp<scalar>());
236  }
237 
238  return area_;
239 }
240 
241 
242 // do not project scalar - just copy values
244 Foam::sampledSurface::project(const Field<scalar>& field) const
245 {
246  tmp<Field<scalar> > tRes(new Field<scalar>(faces().size()));
247  Field<scalar>& res = tRes();
248 
249  forAll(faces(), faceI)
250  {
251  res[faceI] = field[faceI];
252  }
253 
254  return tRes;
255 }
256 
257 
259 Foam::sampledSurface::project(const Field<vector>& field) const
260 {
261  tmp<Field<scalar> > tRes(new Field<scalar>(faces().size()));
262  project(tRes(), field);
263  return tRes;
264 }
265 
266 
268 Foam::sampledSurface::project(const Field<sphericalTensor>& field) const
269 {
270  tmp<Field<vector> > tRes(new Field<vector>(faces().size()));
271  project(tRes(), field);
272  return tRes;
273 }
274 
275 
277 Foam::sampledSurface::project(const Field<symmTensor>& field) const
278 {
279  tmp<Field<vector> > tRes(new Field<vector>(faces().size()));
280  project(tRes(), field);
281  return tRes;
282 }
283 
284 
286 Foam::sampledSurface::project(const Field<tensor>& field) const
287 {
288  tmp<Field<vector> > tRes(new Field<vector>(faces().size()));
289  project(tRes(), field);
290  return tRes;
291 }
292 
293 
295 {
296  os << type();
297 }
298 
299 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
300 
302 {
303  s.print(os);
304  os.check("Ostream& operator<<(Ostream&, const sampledSurface&");
305  return os;
306 }
307 
308 // ************************ vim: set sw=4 sts=4 et: ************************ //