FreeFOAM The Cross-Platform CFD Toolkit
sampledPlane.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 "sampledPlane.H"
27 #include <OpenFOAM/dictionary.H>
28 #include <OpenFOAM/polyMesh.H>
29 #include <finiteVolume/volFields.H>
30 
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(sampledPlane, 0);
38  addNamedToRunTimeSelectionTable(sampledSurface, sampledPlane, word, plane);
39 }
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const word& name,
46  const polyMesh& mesh,
47  const plane& planeDesc,
48  const word& zoneName
49 )
50 :
51  sampledSurface(name, mesh),
52  cuttingPlane(planeDesc),
53  zoneName_(zoneName),
54  needsUpdate_(true)
55 {
56  if (debug && zoneName_.size())
57  {
58  if (mesh.cellZones().findZoneID(zoneName_) < 0)
59  {
60  Info<< "cellZone \"" << zoneName_
61  << "\" not found - using entire mesh" << endl;
62  }
63  }
64 }
65 
66 
68 (
69  const word& name,
70  const polyMesh& mesh,
71  const dictionary& dict
72 )
73 :
74  sampledSurface(name, mesh, dict),
75  cuttingPlane(plane(dict.lookup("basePoint"), dict.lookup("normalVector"))),
76  zoneName_(word::null),
77  needsUpdate_(true)
78 {
79  // make plane relative to the coordinateSystem (Cartesian)
80  // allow lookup from global coordinate systems
81  if (dict.found("coordinateSystem"))
82  {
83  coordinateSystem cs(dict, mesh);
84 
85  point base = cs.globalPosition(planeDesc().refPoint());
86  vector norm = cs.globalVector(planeDesc().normal());
87 
88  // assign the plane description
89  static_cast<plane&>(*this) = plane(base, norm);
90  }
91 
92  dict.readIfPresent("zone", zoneName_);
93 
94  if (debug && zoneName_.size())
95  {
96  if (mesh.cellZones().findZoneID(zoneName_) < 0)
97  {
98  Info<< "cellZone \"" << zoneName_
99  << "\" not found - using entire mesh" << endl;
100  }
101  }
102 
103 }
104 
105 
106 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
107 
109 {}
110 
111 
112 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113 
115 {
116  return needsUpdate_;
117 }
118 
119 
121 {
122  // already marked as expired
123  if (needsUpdate_)
124  {
125  return false;
126  }
127 
129 
130  needsUpdate_ = true;
131  return true;
132 }
133 
134 
136 {
137  if (!needsUpdate_)
138  {
139  return false;
140  }
141 
143 
144  label zoneId = -1;
145  if (zoneName_.size())
146  {
147  zoneId = mesh().cellZones().findZoneID(zoneName_);
148  }
149 
150  if (zoneId < 0)
151  {
152  reCut(mesh());
153  }
154  else
155  {
156  reCut(mesh(), mesh().cellZones()[zoneId]);
157  }
158 
159  if (debug)
160  {
161  print(Pout);
162  Pout << endl;
163  }
164 
165  needsUpdate_ = false;
166  return true;
167 }
168 
169 
172 (
173  const volScalarField& vField
174 ) const
175 {
176  return sampleField(vField);
177 }
178 
179 
182 (
183  const volVectorField& vField
184 ) const
185 {
186  return sampleField(vField);
187 }
188 
189 
192 (
193  const volSphericalTensorField& vField
194 ) const
195 {
196  return sampleField(vField);
197 }
198 
199 
202 (
203  const volSymmTensorField& vField
204 ) const
205 {
206  return sampleField(vField);
207 }
208 
209 
212 (
213  const volTensorField& vField
214 ) const
215 {
216  return sampleField(vField);
217 }
218 
219 
222 (
223  const interpolation<scalar>& interpolator
224 ) const
225 {
226  return interpolateField(interpolator);
227 }
228 
229 
232 (
233  const interpolation<vector>& interpolator
234 ) const
235 {
236  return interpolateField(interpolator);
237 }
238 
241 (
242  const interpolation<sphericalTensor>& interpolator
243 ) const
244 {
245  return interpolateField(interpolator);
246 }
247 
248 
251 (
252  const interpolation<symmTensor>& interpolator
253 ) const
254 {
255  return interpolateField(interpolator);
256 }
257 
258 
261 (
262  const interpolation<tensor>& interpolator
263 ) const
264 {
265  return interpolateField(interpolator);
266 }
267 
268 
270 {
271  os << "sampledPlane: " << name() << " :"
272  << " base:" << refPoint()
273  << " normal:" << normal()
274  << " faces:" << faces().size()
275  << " points:" << points().size();
276 }
277 
278 
279 // ************************ vim: set sw=4 sts=4 et: ************************ //