searchableSurfaceCollection.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 
28 #include "SortableList.H"
29 #include "Time.H"
30 #include "ListOps.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 defineTypeNameAndDebug(searchableSurfaceCollection, 0);
39 (
40  searchableSurface,
41  searchableSurfaceCollection,
42  dict
43 );
44 
45 }
46 
47 
48 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 
50 void Foam::searchableSurfaceCollection::findNearest
51 (
52  const pointField& samples,
53  scalarField& minDistSqr,
54  List<pointIndexHit>& nearestInfo,
55  labelList& nearestSurf
56 ) const
57 {
58  // Initialise
59  nearestInfo.setSize(samples.size());
60  nearestInfo = pointIndexHit();
61  nearestSurf.setSize(samples.size());
62  nearestSurf = -1;
63 
64  List<pointIndexHit> hitInfo(samples.size());
65 
66  const scalarField localMinDistSqr(samples.size(), GREAT);
67 
68  forAll(subGeom_, surfI)
69  {
70  subGeom_[surfI].findNearest
71  (
72  cmptDivide // Transform then divide
73  (
74  transform_[surfI].localPosition(samples),
75  scale_[surfI]
76  ),
77  localMinDistSqr,
78  hitInfo
79  );
80 
81  forAll(hitInfo, pointI)
82  {
83  if (hitInfo[pointI].hit())
84  {
85  // Rework back into global coordinate sys. Multiply then
86  // transform
87  point globalPt = transform_[surfI].globalPosition
88  (
90  (
91  hitInfo[pointI].rawPoint(),
92  scale_[surfI]
93  )
94  );
95 
96  scalar distSqr = magSqr(globalPt - samples[pointI]);
97 
98  if (distSqr < minDistSqr[pointI])
99  {
100  minDistSqr[pointI] = distSqr;
101  nearestInfo[pointI].setPoint(globalPt);
102  nearestInfo[pointI].setHit();
103  nearestInfo[pointI].setIndex
104  (
105  hitInfo[pointI].index()
106  + indexOffset_[surfI]
107  );
108  nearestSurf[pointI] = surfI;
109  }
110  }
111  }
112  }
113 }
114 
115 
116 // Sort hits into per-surface bins. Misses are rejected. Maintains map back
117 // to position
118 void Foam::searchableSurfaceCollection::sortHits
119 (
120  const List<pointIndexHit>& info,
121  List<List<pointIndexHit> >& surfInfo,
122  labelListList& infoMap
123 ) const
124 {
125  // Count hits per surface.
126  labelList nHits(subGeom_.size(), 0);
127 
128  forAll(info, pointI)
129  {
130  if (info[pointI].hit())
131  {
132  label index = info[pointI].index();
133  label surfI = findLower(indexOffset_, index+1);
134  nHits[surfI]++;
135  }
136  }
137 
138  // Per surface the hit
139  surfInfo.setSize(subGeom_.size());
140  // Per surface the original position
141  infoMap.setSize(subGeom_.size());
142 
143  forAll(surfInfo, surfI)
144  {
145  surfInfo[surfI].setSize(nHits[surfI]);
146  infoMap[surfI].setSize(nHits[surfI]);
147  }
148  nHits = 0;
149 
150  forAll(info, pointI)
151  {
152  if (info[pointI].hit())
153  {
154  label index = info[pointI].index();
155  label surfI = findLower(indexOffset_, index+1);
156 
157  // Store for correct surface and adapt indices back to local
158  // ones
159  label localI = nHits[surfI]++;
160  surfInfo[surfI][localI] = pointIndexHit
161  (
162  info[pointI].hit(),
163  info[pointI].rawPoint(),
164  index-indexOffset_[surfI]
165  );
166  infoMap[surfI][localI] = pointI;
167  }
168  }
169 }
170 
171 
172 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
173 
174 Foam::searchableSurfaceCollection::searchableSurfaceCollection
175 (
176  const IOobject& io,
177  const dictionary& dict
178 )
179 :
180  searchableSurface(io),
181  instance_(dict.size()),
182  scale_(dict.size()),
183  transform_(dict.size()),
184  subGeom_(dict.size()),
185  mergeSubRegions_(dict.lookup("mergeSubRegions")),
186  indexOffset_(dict.size()+1)
187 {
188  Info<< "SearchableCollection : " << name() << endl;
189 
190  label surfI = 0;
191  label startIndex = 0;
192  forAllConstIter(dictionary, dict, iter)
193  {
194  if (dict.isDict(iter().keyword()))
195  {
196  instance_[surfI] = iter().keyword();
197 
198  const dictionary& subDict = dict.subDict(instance_[surfI]);
199 
200  scale_[surfI] = subDict.lookup("scale");
201  transform_.set
202  (
203  surfI,
205  (
206  subDict.subDict("transform")
207  )
208  );
209 
210  const word subGeomName(subDict.lookup("surface"));
211  //Pout<< "Trying to find " << subGeomName << endl;
212 
213  const searchableSurface& s =
214  io.db().lookupObject<searchableSurface>(subGeomName);
215 
216  // I don't know yet how to handle the globalSize combined with
217  // regionOffset. Would cause non-consecutive indices locally
218  // if all indices offset by globalSize() of the local region...
219  if (s.size() != s.globalSize())
220  {
222  (
223  "searchableSurfaceCollection::searchableSurfaceCollection"
224  "(const IOobject&, const dictionary&)"
225  ) << "Cannot use a distributed surface in a collection."
226  << exit(FatalError);
227  }
228 
229  subGeom_.set(surfI, &const_cast<searchableSurface&>(s));
230 
231  indexOffset_[surfI] = startIndex;
232  startIndex += subGeom_[surfI].size();
233 
234  Info<< " instance : " << instance_[surfI] << endl;
235  Info<< " surface : " << s.name() << endl;
236  Info<< " scale : " << scale_[surfI] << endl;
237  Info<< " coordsys : " << transform_[surfI] << endl;
238 
239  surfI++;
240  }
241  }
242  indexOffset_[surfI] = startIndex;
243 
244  instance_.setSize(surfI);
245  scale_.setSize(surfI);
246  transform_.setSize(surfI);
247  subGeom_.setSize(surfI);
248  indexOffset_.setSize(surfI+1);
249 
250  // Bounds is the overall bounds
251  bounds() = boundBox(point::max, point::min);
252 
253  forAll(subGeom_, surfI)
254  {
255  const boundBox& surfBb = subGeom_[surfI].bounds();
256 
257  // Transform back to global coordinate sys.
258  const point surfBbMin = transform_[surfI].globalPosition
259  (
261  (
262  surfBb.min(),
263  scale_[surfI]
264  )
265  );
266  const point surfBbMax = transform_[surfI].globalPosition
267  (
269  (
270  surfBb.max(),
271  scale_[surfI]
272  )
273  );
274 
275  bounds().min() = min(bounds().min(), surfBbMin);
276  bounds().max() = max(bounds().max(), surfBbMax);
277  }
278 }
279 
280 
281 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
282 
284 {}
285 
286 
287 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
288 
290 {
291  if (regions_.size() == 0)
292  {
293  regionOffset_.setSize(subGeom_.size());
294 
295  DynamicList<word> allRegions;
296  forAll(subGeom_, surfI)
297  {
298  regionOffset_[surfI] = allRegions.size();
299 
300  if (mergeSubRegions_)
301  {
302  // Single name regardless how many regions subsurface has
303  allRegions.append(instance_[surfI] + "_" + Foam::name(surfI));
304  }
305  else
306  {
307  const wordList& subRegions = subGeom_[surfI].regions();
308 
309  forAll(subRegions, i)
310  {
311  allRegions.append(instance_[surfI] + "_" + subRegions[i]);
312  }
313  }
314  }
315  regions_.transfer(allRegions.shrink());
316  }
317  return regions_;
318 }
319 
320 
322 {
323  return indexOffset_.last();
324 }
325 
326 
329 {
330  tmp<pointField> tCtrs = tmp<pointField>(new pointField(size()));
331  pointField& ctrs = tCtrs();
332 
333  // Append individual coordinates
334  label coordI = 0;
335 
336  forAll(subGeom_, surfI)
337  {
338  const pointField subCoords(subGeom_[surfI].coordinates());
339 
340  forAll(subCoords, i)
341  {
342  ctrs[coordI++] = transform_[surfI].globalPosition
343  (
345  (
346  subCoords[i],
347  scale_[surfI]
348  )
349  );
350  }
351  }
352 
353  return tCtrs;
354 }
355 
356 
358 (
359  pointField& centres,
360  scalarField& radiusSqr
361 ) const
362 {
363  centres.setSize(size());
364  radiusSqr.setSize(centres.size());
365 
366  // Append individual coordinates
367  label coordI = 0;
368 
369  forAll(subGeom_, surfI)
370  {
371  scalar maxScale = cmptMax(scale_[surfI]);
372 
373  pointField subCentres;
374  scalarField subRadiusSqr;
375  subGeom_[surfI].boundingSpheres(subCentres, subRadiusSqr);
376 
377  forAll(subCentres, i)
378  {
379  centres[coordI] = transform_[surfI].globalPosition
380  (
382  (
383  subCentres[i],
384  scale_[surfI]
385  )
386  );
387  radiusSqr[coordI] = maxScale*subRadiusSqr[i];
388  coordI++;
389  }
390  }
391 }
392 
393 
396 {
397  // Get overall size
398  label nPoints = 0;
399 
400  forAll(subGeom_, surfI)
401  {
402  nPoints += subGeom_[surfI].points()().size();
403  }
404 
405  tmp<pointField> tPts(new pointField(nPoints));
406  pointField& pts = tPts();
407 
408  // Append individual coordinates
409  nPoints = 0;
410 
411  forAll(subGeom_, surfI)
412  {
413  const pointField subCoords(subGeom_[surfI].points());
414 
415  forAll(subCoords, i)
416  {
417  pts[nPoints++] = transform_[surfI].globalPosition
418  (
420  (
421  subCoords[i],
422  scale_[surfI]
423  )
424  );
425  }
426  }
427 
428  return tPts;
429 }
430 
431 
432 void Foam::searchableSurfaceCollection::findNearest
433 (
434  const pointField& samples,
435  const scalarField& nearestDistSqr,
436  List<pointIndexHit>& nearestInfo
437 ) const
438 {
439  // How to scale distance?
440  scalarField minDistSqr(nearestDistSqr);
441 
442  labelList nearestSurf;
443  findNearest
444  (
445  samples,
446  minDistSqr,
447  nearestInfo,
448  nearestSurf
449  );
450 }
451 
452 
454 (
455  const pointField& start,
456  const pointField& end,
457  List<pointIndexHit>& info
458 ) const
459 {
460  info.setSize(start.size());
461  info = pointIndexHit();
462 
463  // Current nearest (to start) intersection
464  pointField nearest(end);
465 
466  List<pointIndexHit> hitInfo(start.size());
467 
468  forAll(subGeom_, surfI)
469  {
470  // Starting point
472  (
473  transform_[surfI].localPosition
474  (
475  start
476  ),
477  scale_[surfI]
478  );
479 
480  // Current best end point
482  (
483  transform_[surfI].localPosition
484  (
485  nearest
486  ),
487  scale_[surfI]
488  );
489 
490  subGeom_[surfI].findLine(e0, e1, hitInfo);
491 
492  forAll(hitInfo, pointI)
493  {
494  if (hitInfo[pointI].hit())
495  {
496  // Transform back to global coordinate sys.
497  nearest[pointI] = transform_[surfI].globalPosition
498  (
500  (
501  hitInfo[pointI].rawPoint(),
502  scale_[surfI]
503  )
504  );
505  info[pointI] = hitInfo[pointI];
506  info[pointI].rawPoint() = nearest[pointI];
507  info[pointI].setIndex
508  (
509  hitInfo[pointI].index()
510  + indexOffset_[surfI]
511  );
512  }
513  }
514  }
515 
516 
517  // Debug check
518  if (false)
519  {
520  forAll(info, pointI)
521  {
522  if (info[pointI].hit())
523  {
524  vector n(end[pointI] - start[pointI]);
525  scalar magN = mag(n);
526 
527  if (magN > SMALL)
528  {
529  n /= mag(n);
530 
531  scalar s = ((info[pointI].rawPoint()-start[pointI])&n);
532 
533  if (s < 0 || s > 1)
534  {
536  (
537  "searchableSurfaceCollection::findLine(..)"
538  ) << "point:" << info[pointI]
539  << " s:" << s
540  << " outside vector "
541  << " start:" << start[pointI]
542  << " end:" << end[pointI]
543  << abort(FatalError);
544  }
545  }
546  }
547  }
548  }
549 }
550 
551 
553 (
554  const pointField& start,
555  const pointField& end,
556  List<pointIndexHit>& info
557 ) const
558 {
559  // To be done ...
560  findLine(start, end, info);
561 }
562 
563 
565 (
566  const pointField& start,
567  const pointField& end,
568  List<List<pointIndexHit> >& info
569 ) const
570 {
571  // To be done. Assume for now only one intersection.
572  List<pointIndexHit> nearestInfo;
573  findLine(start, end, nearestInfo);
574 
575  info.setSize(start.size());
576  forAll(info, pointI)
577  {
578  if (nearestInfo[pointI].hit())
579  {
580  info[pointI].setSize(1);
581  info[pointI][0] = nearestInfo[pointI];
582  }
583  else
584  {
585  info[pointI].clear();
586  }
587  }
588 }
589 
590 
592 (
593  const List<pointIndexHit>& info,
594  labelList& region
595 ) const
596 {
597  if (subGeom_.size() == 0)
598  {}
599  else if (subGeom_.size() == 1)
600  {
601  if (mergeSubRegions_)
602  {
603  region.setSize(info.size());
604  region = regionOffset_[0];
605  }
606  else
607  {
608  subGeom_[0].getRegion(info, region);
609  }
610  }
611  else
612  {
613  // Multiple surfaces. Sort by surface.
614 
615  // Per surface the hit
616  List<List<pointIndexHit> > surfInfo;
617  // Per surface the original position
618  List<List<label> > infoMap;
619  sortHits(info, surfInfo, infoMap);
620 
621  region.setSize(info.size());
622  region = -1;
623 
624  // Do region tests
625 
626  if (mergeSubRegions_)
627  {
628  // Actually no need for surfInfo. Just take region for surface.
629  forAll(infoMap, surfI)
630  {
631  const labelList& map = infoMap[surfI];
632  forAll(map, i)
633  {
634  region[map[i]] = regionOffset_[surfI];
635  }
636  }
637  }
638  else
639  {
640  forAll(infoMap, surfI)
641  {
642  labelList surfRegion;
643  subGeom_[surfI].getRegion(surfInfo[surfI], surfRegion);
644 
645  const labelList& map = infoMap[surfI];
646  forAll(map, i)
647  {
648  region[map[i]] = regionOffset_[surfI] + surfRegion[i];
649  }
650  }
651  }
652  }
653 }
654 
655 
657 (
658  const List<pointIndexHit>& info,
659  vectorField& normal
660 ) const
661 {
662  if (subGeom_.size() == 0)
663  {}
664  else if (subGeom_.size() == 1)
665  {
666  subGeom_[0].getNormal(info, normal);
667  }
668  else
669  {
670  // Multiple surfaces. Sort by surface.
671 
672  // Per surface the hit
673  List<List<pointIndexHit> > surfInfo;
674  // Per surface the original position
675  List<List<label> > infoMap;
676  sortHits(info, surfInfo, infoMap);
677 
678  normal.setSize(info.size());
679 
680  // Do region tests
681  forAll(surfInfo, surfI)
682  {
683  vectorField surfNormal;
684  subGeom_[surfI].getNormal(surfInfo[surfI], surfNormal);
685 
686  // Transform back to global coordinate sys.
687  surfNormal = transform_[surfI].globalVector(surfNormal);
688 
689  const labelList& map = infoMap[surfI];
690  forAll(map, i)
691  {
692  normal[map[i]] = surfNormal[i];
693  }
694  }
695  }
696 }
697 
698 
700 (
701  const pointField& points,
702  List<volumeType>& volType
703 ) const
704 {
706  (
707  "searchableSurfaceCollection::getVolumeType(const pointField&"
708  ", List<volumeType>&) const"
709  ) << "Volume type not supported for collection."
710  << exit(FatalError);
711 }
712 
713 
715 (
716  const List<treeBoundBox>& bbs,
717  const bool keepNonLocal,
718  autoPtr<mapDistribute>& faceMap,
719  autoPtr<mapDistribute>& pointMap
720 )
721 {
722  forAll(subGeom_, surfI)
723  {
724  // Note:Tranform the bounding boxes? Something like
725  // pointField bbPoints =
726  // cmptDivide
727  // (
728  // transform_[surfI].localPosition
729  // (
730  // bbs[i].points()
731  // ),
732  // scale_[surfI]
733  // );
734  // treeBoundBox newBb(bbPoints);
735 
736  // Note: what to do with faceMap, pointMap from multiple surfaces?
737  subGeom_[surfI].distribute
738  (
739  bbs,
740  keepNonLocal,
741  faceMap,
742  pointMap
743  );
744  }
745 }
746 
747 
749 {
750  forAll(subGeom_, surfI)
751  {
752  subGeom_[surfI].setField
753  (
754  static_cast<const labelList&>
755  (
757  (
758  values,
759  subGeom_[surfI].size(),
760  indexOffset_[surfI]
761  )
762  )
763  );
764  }
765 }
766 
767 
769 (
770  const List<pointIndexHit>& info,
771  labelList& values
772 ) const
773 {
774  if (subGeom_.size() == 0)
775  {}
776  else if (subGeom_.size() == 1)
777  {
778  subGeom_[0].getField(info, values);
779  }
780  else
781  {
782  // Multiple surfaces. Sort by surface.
783 
784  // Per surface the hit
785  List<List<pointIndexHit> > surfInfo;
786  // Per surface the original position
787  List<List<label> > infoMap;
788  sortHits(info, surfInfo, infoMap);
789 
790  // Do surface tests
791  forAll(surfInfo, surfI)
792  {
793  labelList surfValues;
794  subGeom_[surfI].getField(surfInfo[surfI], surfValues);
795 
796  if (surfValues.size())
797  {
798  // Size values only when we have a surface that supports it.
799  values.setSize(info.size());
800 
801  const labelList& map = infoMap[surfI];
802  forAll(map, i)
803  {
804  values[map[i]] = surfValues[i];
805  }
806  }
807  }
808  }
809 }
810 
811 
812 // ************************************************************************* //
const pointField & points
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
vector point
Point is a vector.
Definition: point.H:41
PtrList< coordinateSystem > coordinates(solidRegions.size())
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
virtual tmp< pointField > points() const
Get the points that define the surface.
Base class of (analytical or triangulated) surface. Encapsulates all the search routines. WIP.
dimensioned< scalar > mag(const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual void getVolumeType(const pointField &, List< volumeType > &) const
Determine type (inside/outside/mixed) for point. unknown if.
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:310
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit > > &) const
Get all intersections in order from start to end.
Various functions to operate on Lists.
static const Vector max
Definition: Vector.H:82
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
messageStream Info
virtual void getField(const List< pointIndexHit > &, labelList &) const
WIP. From a set of hits (points and.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
virtual void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Find first intersection on segment from start to end.
label size() const
Return number of elements in list.
Definition: DLListBaseI.H:77
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
label n
virtual void getRegion(const List< pointIndexHit > &, labelList &region) const
From a set of points and indices get the region.
void clear()
Clear the list, i.e. set size to zero.
Definition: List.C:379
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Return any intersection on segment from start to end.
dimensioned< Type > cmptDivide(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void distribute(const List< treeBoundBox > &, const bool keepNonLocal, autoPtr< mapDistribute > &faceMap, autoPtr< mapDistribute > &pointMap)
Set bounds of surface. Bounds currently set as list of.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
virtual void boundingSpheres(pointField &centres, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared), one per element.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & shrink()
Shrink the allocated space to the number of elements used.
Definition: DynamicListI.H:258
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
label findLower(const ListType &, typename ListType::const_reference, const label stary, const BinaryOp &bop)
Find last element < given value in sorted list and return index,.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:55
static const Vector min
Definition: Vector.H:83
bool isDict(const word &) const
Check if entry is a sub-dictionary.
Definition: dictionary.C:607
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
error FatalError
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual const wordList & regions() const
Names of regions.
virtual void setField(const labelList &values)
WIP. Store element-wise field.
A List obtained as a section of another List.
Definition: SubList.H:53
virtual label size() const
Range of local indices that can be returned.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
label nPoints
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
static autoPtr< coordinateSystem > New(const objectRegistry &obr, const dictionary &dict)
Select constructed from dictionary and objectRegistry.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
A class for managing temporary objects.
Definition: PtrList.H:118
defineTypeNameAndDebug(combustionModel, 0)