FreeFOAM The Cross-Platform CFD Toolkit
sampledPatch.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 "sampledPatch.H"
27 #include <OpenFOAM/dictionary.H>
28 #include <OpenFOAM/polyMesh.H>
29 #include <OpenFOAM/polyPatch.H>
30 #include <finiteVolume/volFields.H>
31 
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(sampledPatch, 0);
39  addNamedToRunTimeSelectionTable(sampledSurface, sampledPatch, word, patch);
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  const word& name,
48  const polyMesh& mesh,
49  const word& patchName,
50  const bool triangulate
51 )
52 :
53  sampledSurface(name, mesh),
54  patchName_(patchName),
55  triangulate_(triangulate),
56  needsUpdate_(true),
57  patchFaceLabels_(0)
58 {}
59 
60 
62 (
63  const word& name,
64  const polyMesh& mesh,
65  const dictionary& dict
66 )
67 :
68  sampledSurface(name, mesh, dict),
69  patchName_(dict.lookup("patchName")),
70  triangulate_(dict.lookupOrDefault("triangulate", false)),
71  needsUpdate_(true),
72  patchFaceLabels_(0)
73 {}
74 
75 
76 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
77 
79 {}
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
85 {
86  return needsUpdate_;
87 }
88 
89 
91 {
92  // already marked as expired
93  if (needsUpdate_)
94  {
95  return false;
96  }
97 
100  patchFaceLabels_.clear();
101 
102  needsUpdate_ = true;
103  return true;
104 }
105 
106 
108 {
109  if (!needsUpdate_)
110  {
111  return false;
112  }
113 
114  const label patchI = mesh().boundaryMesh().findPatchID(patchName_);
115 
116  if (patchI != -1)
117  {
118  const polyPatch& p = mesh().boundaryMesh()[patchI];
119  this->storedPoints() = p.localPoints();
120  this->storedFaces() = p.localFaces();
121 
122  // an identity map
123  patchFaceLabels_.setSize(faces().size());
124  forAll(patchFaceLabels_, i)
125  {
126  patchFaceLabels_[i] = i;
127  }
128 
129  // triangulate uses remapFaces()
130  // - this is somewhat less efficient since it recopies the faces
131  // that we just created, but we probably don't want to do this
132  // too often anyhow.
133  if (triangulate_)
134  {
135  MeshStorage::triangulate();
136  }
137  }
138 
139  if (debug)
140  {
141  print(Pout);
142  Pout<< endl;
143  }
144 
145  needsUpdate_ = false;
146  return true;
147 }
148 
149 
150 // remap action on triangulation
151 void Foam::sampledPatch::remapFaces
152 (
153  const UList<label>& faceMap
154 )
155 {
156  // recalculate the cells cut
157  if (&faceMap && faceMap.size())
158  {
159  MeshStorage::remapFaces(faceMap);
160  patchFaceLabels_ = labelList
161  (
162  UIndirectList<label>(patchFaceLabels_, faceMap)
163  );
164  }
165 }
166 
167 
168 
171 (
172  const volScalarField& vField
173 ) const
174 {
175  return sampleField(vField);
176 }
177 
178 
181 (
182  const volVectorField& vField
183 ) const
184 {
185  return sampleField(vField);
186 }
187 
190 (
191  const volSphericalTensorField& vField
192 ) const
193 {
194  return sampleField(vField);
195 }
196 
197 
200 (
201  const volSymmTensorField& vField
202 ) const
203 {
204  return sampleField(vField);
205 }
206 
207 
210 (
211  const volTensorField& vField
212 ) const
213 {
214  return sampleField(vField);
215 }
216 
217 
220 (
221  const interpolation<scalar>& interpolator
222 ) const
223 {
224  return interpolateField(interpolator);
225 }
226 
227 
230 (
231  const interpolation<vector>& interpolator
232 ) const
233 {
234  return interpolateField(interpolator);
235 }
236 
239 (
240  const interpolation<sphericalTensor>& interpolator
241 ) const
242 {
243  return interpolateField(interpolator);
244 }
245 
246 
249 (
250  const interpolation<symmTensor>& interpolator
251 ) const
252 {
253  return interpolateField(interpolator);
254 }
255 
256 
259 (
260  const interpolation<tensor>& interpolator
261 ) const
262 {
263  return interpolateField(interpolator);
264 }
265 
266 
268 {
269  os << "sampledPatch: " << name() << " :"
270  << " patch:" << patchName()
271  << " faces:" << faces().size()
272  << " points:" << points().size();
273 }
274 
275 
276 // ************************ vim: set sw=4 sts=4 et: ************************ //