FreeFOAM The Cross-Platform CFD Toolkit
faceSource.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) 2009-2011 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 "faceSource.H"
27 #include <finiteVolume/fvMesh.H>
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  namespace fieldValues
38  {
39  defineTypeNameAndDebug(faceSource, 0);
40  }
41 
42  template<>
43  const char* NamedEnum<fieldValues::faceSource::sourceType, 3>::
44  names[] = {"faceZone", "patch", "sampledSurface"};
45 
46  const NamedEnum<fieldValues::faceSource::sourceType, 3>
48 
49  template<>
50  const char* NamedEnum<fieldValues::faceSource::operationType, 7>::
51  names[] =
52  {
53  "none", "sum", "areaAverage",
54  "areaIntegrate", "weightedAverage", "min", "max"
55  };
56 
57  const NamedEnum<fieldValues::faceSource::operationType, 7>
59 
60 }
61 
62 
63 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
64 
65 void Foam::fieldValues::faceSource::setFaceZoneFaces()
66 {
67  label zoneId = mesh().faceZones().findZoneID(sourceName_);
68 
69  if (zoneId < 0)
70  {
71  FatalErrorIn("faceSource::faceSource::setFaceZoneFaces()")
72  << type() << " " << name_ << ": "
73  << sourceTypeNames_[source_] << "(" << sourceName_ << "):" << nl
74  << " Unknown face zone name: " << sourceName_
75  << ". Valid face zones are: " << mesh().faceZones().names()
76  << nl << exit(FatalError);
77  }
78 
79  const faceZone& fZone = mesh().faceZones()[zoneId];
80 
81  faceId_.setSize(fZone.size());
82  facePatchId_.setSize(fZone.size());
83  faceSign_.setSize(fZone.size());
84 
85  label count = 0;
86  forAll(fZone, i)
87  {
88  label faceI = fZone[i];
89 
90  label faceId = -1;
91  label facePatchId = -1;
92  if (mesh().isInternalFace(faceI))
93  {
94  faceId = faceI;
95  facePatchId = -1;
96  }
97  else
98  {
99  facePatchId = mesh().boundaryMesh().whichPatch(faceI);
100  const polyPatch& pp = mesh().boundaryMesh()[facePatchId];
101  if (isA<processorPolyPatch>(pp))
102  {
103  if (refCast<const processorPolyPatch>(pp).owner())
104  {
105  faceId = pp.whichFace(faceI);
106  }
107  else
108  {
109  faceId = -1;
110  }
111  }
112  else if (isA<cyclicPolyPatch>(pp))
113  {
114  label patchFaceI = faceI - pp.start();
115  if (patchFaceI < pp.size()/2)
116  {
117  faceId = patchFaceI;
118  }
119  else
120  {
121  faceId = -1;
122  }
123  }
124  else if (!isA<emptyPolyPatch>(pp))
125  {
126  faceId = faceI - pp.start();
127  }
128  else
129  {
130  faceId = -1;
131  facePatchId = -1;
132  }
133  }
134 
135  if (faceId >= 0)
136  {
137  if (fZone.flipMap()[i])
138  {
139  faceSign_[count] = -1;
140  }
141  else
142  {
143  faceSign_[count] = 1;
144  }
145  faceId_[count] = faceId;
146  facePatchId_[count] = facePatchId;
147  count++;
148  }
149  }
150 
151  faceId_.setSize(count);
152  facePatchId_.setSize(count);
153  faceSign_.setSize(count);
154  nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
155 
156  if (debug)
157  {
158  Pout<< "Original face zone size = " << fZone.size()
159  << ", new size = " << count << endl;
160  }
161 }
162 
163 
164 void Foam::fieldValues::faceSource::setPatchFaces()
165 {
166  label patchId = mesh().boundaryMesh().findPatchID(sourceName_);
167 
168  if (patchId < 0)
169  {
170  FatalErrorIn("faceSource::constructFaceAddressing()")
171  << type() << " " << name_ << ": "
172  << sourceTypeNames_[source_] << "(" << sourceName_ << "):" << nl
173  << " Unknown patch name: " << sourceName_
174  << ". Valid patch names are: "
175  << mesh().boundaryMesh().names() << nl
176  << exit(FatalError);
177  }
178 
179  const polyPatch& pp = mesh().boundaryMesh()[patchId];
180 
181  label nFaces = pp.size();
182  if (isA<cyclicPolyPatch>(pp))
183  {
184  nFaces /= 2;
185  }
186  else if (isA<emptyPolyPatch>(pp))
187  {
188  nFaces = 0;
189  }
190 
191  faceId_.setSize(nFaces);
192  facePatchId_.setSize(nFaces);
193  faceSign_.setSize(nFaces);
194  nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
195 
196  forAll(faceId_, faceI)
197  {
198  faceId_[faceI] = faceI;
199  facePatchId_[faceI] = patchId;
200  faceSign_[faceI] = 1;
201  }
202 }
203 
204 
205 void Foam::fieldValues::faceSource::sampledSurfaceFaces(const dictionary& dict)
206 {
207  surfacePtr_ = sampledSurface::New
208  (
209  name_,
210  mesh(),
211  dict.subDict("sampledSurfaceDict")
212  );
213  surfacePtr_().update();
214  nFaces_ = returnReduce(surfacePtr_().faces().size(), sumOp<label>());
215 }
216 
217 
218 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
219 
221 {
222  switch (source_)
223  {
224  case stFaceZone:
225  {
226  setFaceZoneFaces();
227  break;
228  }
229  case stPatch:
230  {
231  setPatchFaces();
232  break;
233  }
234  case stSampledSurface:
235  {
236  sampledSurfaceFaces(dict);
237  break;
238  }
239  default:
240  {
241  FatalErrorIn("faceSource::initialise()")
242  << type() << " " << name_ << ": "
243  << sourceTypeNames_[source_] << "(" << sourceName_ << "):"
244  << nl << " Unknown source type. Valid source types are:"
245  << sourceTypeNames_ << nl << exit(FatalError);
246  }
247  }
248 
249  scalar totalArea;
250 
251  if (surfacePtr_.valid())
252  {
253  surfacePtr_().update();
254  totalArea = gSum(surfacePtr_().magSf());
255  }
256  else
257  {
258  totalArea = gSum(filterField(mesh().magSf(), false));
259  }
260 
261  Info<< type() << " " << name_ << ":" << nl
262  << " total faces = " << nFaces_
263  << nl
264  << " total area = " << totalArea
265  << nl;
266 
267  if (operation_ == opWeightedAverage)
268  {
269  dict.lookup("weightField") >> weightFieldName_;
270  if
271  (
272  obr().foundObject<volScalarField>(weightFieldName_)
273  || obr().foundObject<surfaceScalarField>(weightFieldName_)
274  )
275  {
276  Info<< " weight field = " << weightFieldName_;
277  }
278  else
279  {
280  FatalErrorIn("faceSource::initialise()")
281  << type() << " " << name_ << ": "
282  << sourceTypeNames_[source_] << "(" << sourceName_ << "):"
283  << nl << " Weight field " << weightFieldName_
284  << " must be either a " << volScalarField::typeName << " or "
286  }
287  }
288 
289  Info<< nl << endl;
290 }
291 
292 
294 {
295  if (outputFilePtr_.valid())
296  {
297  outputFilePtr_()
298  << "# Source : " << sourceTypeNames_[source_] << " "
299  << sourceName_ << nl << "# Faces : " << nFaces_ << nl
300  << "# Time" << tab << "sum(magSf)";
301 
302  forAll(fields_, i)
303  {
304  outputFilePtr_()
305  << tab << operationTypeNames_[operation_]
306  << "(" << fields_[i] << ")";
307  }
308 
309  outputFilePtr_() << endl;
310  }
311 }
312 
313 
314 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
315 
317 (
318  const word& name,
319  const objectRegistry& obr,
320  const dictionary& dict,
321  const bool loadFromFiles
322 )
323 :
324  fieldValue(name, obr, dict, loadFromFiles),
325  source_(sourceTypeNames_.read(dict.lookup("source"))),
326  operation_(operationTypeNames_.read(dict.lookup("operation"))),
327  weightFieldName_("undefinedWeightedFieldName"),
328  nFaces_(0),
329  faceId_(),
330  facePatchId_(),
331  faceSign_()
332 {
333  read(dict);
334 }
335 
336 
337 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
338 
340 {}
341 
342 
343 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
344 
346 {
347  fieldValue::read(dict);
348 
349  if (active_)
350  {
351  initialise(dict);
352  }
353 }
354 
355 
357 {
359 
360  if (active_)
361  {
362  scalar totalArea;
363 
364  if (surfacePtr_.valid())
365  {
366  surfacePtr_().update();
367  totalArea = gSum(surfacePtr_().magSf());
368  }
369  else
370  {
371  totalArea = gSum(filterField(mesh().magSf(), false));
372  }
373 
374  if (Pstream::master())
375  {
376  outputFilePtr_() << obr_.time().value() << tab << totalArea;
377  }
378 
379  forAll(fields_, i)
380  {
381  writeValues<scalar>(fields_[i]);
382  writeValues<vector>(fields_[i]);
383  writeValues<sphericalTensor>(fields_[i]);
384  writeValues<symmTensor>(fields_[i]);
385  writeValues<tensor>(fields_[i]);
386  }
387 
388  if (Pstream::master())
389  {
390  outputFilePtr_()<< endl;
391  }
392 
393  if (log_)
394  {
395  Info<< endl;
396  }
397  }
398 }
399 
400 
401 // ************************ vim: set sw=4 sts=4 et: ************************ //