searchableSurfaces.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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 
26 #include "searchableSurfaces.H"
28 #include "ListOps.H"
29 #include "Time.H"
30 #include "DynamicField.H"
31 #include "PatchTools.H"
32 #include "triSurfaceMesh.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(searchableSurfaces, 0);
39 }
40 
41 
42 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 bool Foam::searchableSurfaces::connected
45 (
46  const triSurface& s,
47  const label edgeI,
48  const pointIndexHit& hit
49 )
50 {
51  const triFace& localFace = s.localFaces()[hit.index()];
52  const edge& e = s.edges()[edgeI];
53 
54  forAll(localFace, i)
55  {
56  if (e.otherVertex(localFace[i]) != -1)
57  {
58  return true;
59  }
60  }
61 
62  return false;
63 }
64 
65 
66 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
67 
68 // Construct with length.
70 :
72  regionNames_(size),
73  allSurfaces_(identity(size))
74 {}
75 
76 
77 //Foam::searchableSurfaces::searchableSurfaces
78 //(
79 // const IOobject& io,
80 // const PtrList<dictionary>& dicts
81 //)
82 //:
83 // PtrList<searchableSurface>(dicts.size()),
84 // regionNames_(dicts.size()),
85 // allSurfaces_(identity(dicts.size()))
86 //{
87 // forAll(dicts, surfI)
88 // {
89 // const dictionary& dict = dicts[surfI];
90 //
91 // // Make IOobject with correct name
92 // autoPtr<IOobject> namedIO(io.clone());
93 // namedIO().rename(dict.lookup("name"));
94 //
95 // // Create and hook surface
96 // set
97 // (
98 // surfI,
99 // searchableSurface::New
100 // (
101 // dict.lookup("type"),
102 // namedIO(),
103 // dict
104 // )
105 // );
106 // const searchableSurface& s = operator[](surfI);
107 //
108 // // Construct default region names by prepending surface name
109 // // to region name.
110 // const wordList& localNames = s.regions();
111 //
112 // wordList globalNames(localNames.size());
113 // forAll(localNames, regionI)
114 // {
115 // globalNames[regionI] = s.name() + '_' + localNames[regionI];
116 // }
117 //
118 // // See if dictionary provides any global region names.
119 // if (dict.found("regions"))
120 // {
121 // const dictionary& regionsDict = dict.subDict("regions");
122 //
123 // forAllConstIter(dictionary, regionsDict, iter)
124 // {
125 // const word& key = iter().keyword();
126 //
127 // if (regionsDict.isDict(key))
128 // {
129 // // Get the dictionary for region iter.key()
130 // const dictionary& regionDict = regionsDict.subDict(key);
131 //
132 // label index = findIndex(localNames, key);
133 //
134 // if (index == -1)
135 // {
136 // FatalErrorInFunction
137 // << "Unknown region name " << key
138 // << " for surface " << s.name() << endl
139 // << "Valid region names are " << localNames
140 // << exit(FatalError);
141 // }
142 //
143 // globalNames[index] = word(regionDict.lookup("name"));
144 // }
145 // }
146 // }
147 //
148 // // Now globalNames contains the names of the regions.
149 // Info<< "Surface:" << s.name() << " has regions:"
150 // << endl;
151 // forAll(globalNames, regionI)
152 // {
153 // Info<< " " << globalNames[regionI] << endl;
154 // }
155 //
156 // // Create reverse lookup
157 // forAll(globalNames, regionI)
158 // {
159 // regionNames_.insert
160 // (
161 // globalNames[regionI],
162 // labelPair(surfI, regionI)
163 // );
164 // }
165 // }
166 //}
167 
168 
170 (
171  const IOobject& io,
172  const dictionary& topDict,
173  const bool singleRegionName
174 )
175 :
176  PtrList<searchableSurface>(topDict.size()),
177  names_(topDict.size()),
178  regionNames_(topDict.size()),
179  allSurfaces_(identity(topDict.size()))
180 {
181  label surfI = 0;
182  forAllConstIter(dictionary, topDict, iter)
183  {
184  const word& key = iter().keyword();
185 
186  if (!topDict.isDict(key))
187  {
189  << "Found non-dictionary entry " << iter()
190  << " in top-level dictionary " << topDict
191  << exit(FatalError);
192  }
193 
194  const dictionary& dict = topDict.subDict(key);
195 
196  names_[surfI] = key;
197  dict.readIfPresent("name", names_[surfI]);
198 
199  // Make IOobject with correct name
200  autoPtr<IOobject> namedIO(io.clone());
201  // Note: we would like to e.g. register triSurface 'sphere.stl' as
202  // 'sphere'. Unfortunately
203  // no support for having object read from different location than
204  // their object name. Maybe have stlTriSurfaceMesh which appends .stl
205  // when reading/writing?
206  namedIO().rename(key); // names_[surfI]
207 
208  // Create and hook surface
209  set
210  (
211  surfI,
213  (
214  dict.lookup("type"),
215  namedIO(),
216  dict
217  )
218  );
219  const searchableSurface& s = operator[](surfI);
220 
221  // Construct default region names by prepending surface name
222  // to region name.
223  const wordList& localNames = s.regions();
224 
225  wordList& rNames = regionNames_[surfI];
226  rNames.setSize(localNames.size());
227 
228  if (singleRegionName && localNames.size() == 1)
229  {
230  rNames[0] = names_[surfI];
231  }
232  else
233  {
234  forAll(localNames, regionI)
235  {
236  rNames[regionI] = names_[surfI] + '_' + localNames[regionI];
237  }
238  }
239 
240  // See if dictionary provides any global region names.
241  if (dict.found("regions"))
242  {
243  const dictionary& regionsDict = dict.subDict("regions");
244 
245  forAllConstIter(dictionary, regionsDict, iter)
246  {
247  const word& key = iter().keyword();
248 
249  if (regionsDict.isDict(key))
250  {
251  // Get the dictionary for region iter.keyword()
252  const dictionary& regionDict = regionsDict.subDict(key);
253 
254  label index = findIndex(localNames, key);
255 
256  if (index == -1)
257  {
259  << "Unknown region name " << key
260  << " for surface " << s.name() << endl
261  << "Valid region names are " << localNames
262  << exit(FatalError);
263  }
264 
265  rNames[index] = word(regionDict.lookup("name"));
266  }
267  }
268  }
269 
270  surfI++;
271  }
272 
273  // Trim (not really necessary since we don't allow non-dictionary entries)
275  names_.setSize(surfI);
276  regionNames_.setSize(surfI);
277  allSurfaces_.setSize(surfI);
278 }
279 
280 
281 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
282 
284 (
285  const word& wantedName
286 ) const
287 {
288  return findIndex(names_, wantedName);
289 }
290 
291 
293 (
294  const word& surfaceName,
295  const word& regionName
296 ) const
297 {
298  label surfaceIndex = findSurfaceID(surfaceName);
299 
300  return findIndex(this->operator[](surfaceIndex).regions(), regionName);
301 }
302 
303 
304 // Find any intersection
306 (
307  const pointField& start,
308  const pointField& end,
309  labelList& hitSurfaces,
310  List<pointIndexHit>& hitInfo
311 ) const
312 {
314  (
315  *this,
316  allSurfaces_,
317  start,
318  end,
319  hitSurfaces,
320  hitInfo
321  );
322 }
323 
324 
326 (
327  const pointField& start,
328  const pointField& end,
329  labelListList& hitSurfaces,
330  List<List<pointIndexHit>>& hitInfo
331 ) const
332 {
334  (
335  *this,
336  allSurfaces_,
337  start,
338  end,
339  hitSurfaces,
340  hitInfo
341  );
342 }
343 
344 
345 //Find intersections of edge nearest to both endpoints.
347 (
348  const pointField& start,
349  const pointField& end,
350  labelList& surface1,
351  List<pointIndexHit>& hit1,
352  labelList& surface2,
353  List<pointIndexHit>& hit2
354 ) const
355 {
357  (
358  *this,
359  allSurfaces_,
360  start,
361  end,
362  surface1,
363  hit1,
364  surface2,
365  hit2
366  );
367 }
368 
369 
371 (
372  const pointField& samples,
373  const scalarField& nearestDistSqr,
374  labelList& nearestSurfaces,
375  List<pointIndexHit>& nearestInfo
376 ) const
377 {
379  (
380  *this,
381  allSurfaces_,
382  samples,
383  nearestDistSqr,
384  nearestSurfaces,
385  nearestInfo
386  );
387 }
388 
389 
391 (
392  const pointField& samples,
393  const scalarField& nearestDistSqr,
394  const labelList& regionIndices,
395  labelList& nearestSurfaces,
396  List<pointIndexHit>& nearestInfo
397 ) const
398 {
400  (
401  *this,
402  allSurfaces_,
403  samples,
404  nearestDistSqr,
405  regionIndices,
406  nearestSurfaces,
407  nearestInfo
408  );
409 }
410 
411 
413 {
415  (
416  *this,
417  allSurfaces_
418  );
419 }
420 
421 
422 bool Foam::searchableSurfaces::checkClosed(const bool report) const
423 {
424  if (report)
425  {
426  Info<< "Checking for closedness." << endl;
427  }
428 
429  bool hasError = false;
430 
431  forAll(*this, surfI)
432  {
433  if (!operator[](surfI).hasVolumeType())
434  {
435  hasError = true;
436 
437  if (report)
438  {
439  Info<< " " << names()[surfI]
440  << " : not closed" << endl;
441  }
442 
443  if (isA<triSurface>(operator[](surfI)))
444  {
445  const triSurface& s = dynamic_cast<const triSurface&>
446  (
447  operator[](surfI)
448  );
449  const labelListList& edgeFaces = s.edgeFaces();
450 
451  label nSingleEdges = 0;
452  forAll(edgeFaces, edgeI)
453  {
454  if (edgeFaces[edgeI].size() == 1)
455  {
456  nSingleEdges++;
457  }
458  }
459 
460  label nMultEdges = 0;
461  forAll(edgeFaces, edgeI)
462  {
463  if (edgeFaces[edgeI].size() > 2)
464  {
465  nMultEdges++;
466  }
467  }
468 
469  if (report && (nSingleEdges != 0 || nMultEdges != 0))
470  {
471  Info<< " connected to one face : "
472  << nSingleEdges << nl
473  << " connected to >2 faces : "
474  << nMultEdges << endl;
475  }
476  }
477  }
478  }
479 
480  if (report)
481  {
482  Info<< endl;
483  }
484 
485  return returnReduce(hasError, orOp<bool>());
486 }
487 
488 
490 {
491  if (report)
492  {
493  Info<< "Checking for normal orientation." << endl;
494  }
495 
496  bool hasError = false;
497 
498  forAll(*this, surfI)
499  {
500  if (isA<triSurface>(operator[](surfI)))
501  {
502  const triSurface& s = dynamic_cast<const triSurface&>
503  (
504  operator[](surfI)
505  );
506 
507  labelHashSet borderEdge(s.size()/1000);
508  PatchTools::checkOrientation(s, false, &borderEdge);
509  // Colour all faces into zones using borderEdge
510  labelList normalZone;
511  label nZones = PatchTools::markZones(s, borderEdge, normalZone);
512  if (nZones > 1)
513  {
514  hasError = true;
515 
516  if (report)
517  {
518  Info<< " " << names()[surfI]
519  << " : has multiple orientation zones ("
520  << nZones << ")" << endl;
521  }
522  }
523  }
524  }
525 
526  if (report)
527  {
528  Info<< endl;
529  }
530 
531  return returnReduce(hasError, orOp<bool>());
532 }
533 
534 
536 (
537  const scalar maxRatio,
538  const bool report
539 ) const
540 {
541  if (report)
542  {
543  Info<< "Checking for size." << endl;
544  }
545 
546  bool hasError = false;
547 
548  forAll(*this, i)
549  {
550  const boundBox& bb = operator[](i).bounds();
551 
552  for (label j = i+1; j < size(); j++)
553  {
554  scalar ratio = bb.mag()/operator[](j).bounds().mag();
555 
556  if (ratio > maxRatio || ratio < 1.0/maxRatio)
557  {
558  hasError = true;
559 
560  if (report)
561  {
562  Info<< " " << names()[i]
563  << " bounds differ from " << names()[j]
564  << " by more than a factor 100:" << nl
565  << " bounding box : " << bb << nl
566  << " bounding box : " << operator[](j).bounds()
567  << endl;
568  }
569  break;
570  }
571  }
572  }
573 
574  if (report)
575  {
576  Info<< endl;
577  }
578 
579  return returnReduce(hasError, orOp<bool>());
580 }
581 
582 
584 (
585  const scalar tolerance,
586  const autoPtr<writer<scalar>>& setWriter,
587  const bool report
588 ) const
589 {
590  if (report)
591  {
592  Info<< "Checking for intersection." << endl;
593  }
594 
595  // cpuTime timer;
596 
597  bool hasError = false;
598 
599  forAll(*this, i)
600  {
601  if (isA<triSurfaceMesh>(operator[](i)))
602  {
603  const triSurfaceMesh& s0 = dynamic_cast<const triSurfaceMesh&>
604  (
605  operator[](i)
606  );
607  const edgeList& edges0 = s0.edges();
608  const pointField& localPoints0 = s0.localPoints();
609 
610  // Collect all the edge vectors
611  pointField start(edges0.size());
612  pointField end(edges0.size());
613  forAll(edges0, edgeI)
614  {
615  const edge& e = edges0[edgeI];
616  start[edgeI] = localPoints0[e[0]];
617  end[edgeI] = localPoints0[e[1]];
618  }
619 
620  // Test all other surfaces for intersection
621  forAll(*this, j)
622  {
623  List<pointIndexHit> hits;
624 
625  if (i == j)
626  {
627  // Slightly shorten the edges to prevent finding lots of
628  // intersections. The fast triangle intersection routine
629  // has problems with rays passing through a point of the
630  // triangle.
631  vectorField d
632  (
633  max(tolerance, 10*s0.tolerance())
634  *(end-start)
635  );
636  start += d;
637  end -= d;
638  }
639 
640  operator[](j).findLineAny(start, end, hits);
641 
642  label nHits = 0;
643  DynamicField<point> intersections(edges0.size()/100);
644  DynamicField<scalar> intersectionEdge(intersections.capacity());
645 
646  forAll(hits, edgeI)
647  {
648  if
649  (
650  hits[edgeI].hit()
651  && (i != j || !connected(s0, edgeI, hits[edgeI]))
652  )
653  {
654  intersections.append(hits[edgeI].hitPoint());
655  intersectionEdge.append(1.0*edgeI);
656  nHits++;
657  }
658  }
659 
660  if (nHits > 0)
661  {
662  if (report)
663  {
664  Info<< " " << names()[i]
665  << " intersects " << names()[j]
666  << " at " << nHits
667  << " locations."
668  << endl;
669 
670  // vtkSetWriter<scalar> setWriter;
671  if (setWriter.valid())
672  {
673  scalarField dist(mag(intersections));
674  coordSet track
675  (
676  names()[i] + '_' + names()[j],
677  "xyz",
678  move(intersections),
679  dist
680  );
681  wordList valueSetNames(1, "edgeIndex");
682  List<const scalarField*> valueSets
683  (
684  1,
685  &intersectionEdge
686  );
687 
688  fileName fName
689  (
690  setWriter().getFileName(track, valueSetNames)
691  );
692  Info<< " Writing intersection locations to "
693  << fName << endl;
694  OFstream os
695  (
696  s0.searchableSurface::time().path()
697  /fName
698  );
699  setWriter().write
700  (
701  track,
702  valueSetNames,
703  valueSets,
704  os
705  );
706  }
707  }
708 
709  hasError = true;
710  break;
711  }
712  }
713  }
714  }
715 
716  if (report)
717  {
718  Info<< endl;
719  }
720 
721  return returnReduce(hasError, orOp<bool>());
722 }
723 
724 
726 (
727  const scalar minQuality,
728  const bool report
729 ) const
730 {
731  if (report)
732  {
733  Info<< "Checking for triangle quality." << endl;
734  }
735 
736  bool hasError = false;
737 
738  forAll(*this, surfI)
739  {
740  if (isA<triSurface>(operator[](surfI)))
741  {
742  const triSurface& s = dynamic_cast<const triSurface&>
743  (
744  operator[](surfI)
745  );
746 
747  label nBadTris = 0;
748  forAll(s, facei)
749  {
750  const labelledTri& f = s[facei];
751 
752  scalar q = triPointRef
753  (
754  s.points()[f[0]],
755  s.points()[f[1]],
756  s.points()[f[2]]
757  ).quality();
758 
759  if (q < minQuality)
760  {
761  nBadTris++;
762  }
763  }
764 
765  if (nBadTris > 0)
766  {
767  hasError = true;
768 
769  if (report)
770  {
771  Info<< " " << names()[surfI]
772  << " : has " << nBadTris << " bad quality triangles "
773  << " (quality < " << minQuality << ")" << endl;
774  }
775  }
776  }
777  }
778 
779  if (report)
780  {
781  Info<< endl;
782  }
783 
784  return returnReduce(hasError, orOp<bool>());
785 
786 }
787 
788 
789 // Check all surfaces. Return number of failures.
791 (
792  const bool report
793 ) const
794 {
795  label noFailedChecks = 0;
796 
797  if (checkClosed(report))
798  {
799  noFailedChecks++;
800  }
801 
802  if (checkNormalOrientation(report))
803  {
804  noFailedChecks++;
805  }
806  return noFailedChecks;
807 }
808 
809 
811 (
812  const scalar maxRatio,
813  const scalar tol,
814  const autoPtr<writer<scalar>>& setWriter,
815  const scalar minQuality,
816  const bool report
817 ) const
818 {
819  label noFailedChecks = 0;
820 
821  if (maxRatio > 0 && checkSizes(maxRatio, report))
822  {
823  noFailedChecks++;
824  }
825 
826  if (checkIntersection(tol, setWriter, report))
827  {
828  noFailedChecks++;
829  }
830 
831  if (checkQuality(minQuality, report))
832  {
833  noFailedChecks++;
834  }
835 
836  return noFailedChecks;
837 }
838 
839 
841 (
842  const List<wordList>& patchTypes,
843  Ostream& os
844 ) const
845 {
846  Info<< "Statistics." << endl;
847  forAll(*this, surfI)
848  {
849  Info<< " " << names()[surfI] << ':' << endl;
850 
851  const searchableSurface& s = operator[](surfI);
852 
853  Info<< " type : " << s.type() << nl
854  << " size : " << s.globalSize() << nl;
855  if (isA<triSurfaceMesh>(s))
856  {
857  const triSurfaceMesh& ts = dynamic_cast<const triSurfaceMesh&>(s);
858  Info<< " edges : " << ts.nEdges() << nl
859  << " points : " << ts.points()().size() << nl;
860  }
861  Info<< " bounds : " << s.bounds() << nl
862  << " closed : " << Switch(s.hasVolumeType()) << endl;
863 
864  if (patchTypes.size() && patchTypes[surfI].size() >= 1)
865  {
866  wordList unique(HashSet<word>(patchTypes[surfI]).sortedToc());
867  Info<< " patches : ";
868  forAll(unique, i)
869  {
870  Info<< unique[i];
871  if (i < unique.size()-1)
872  {
873  Info<< ',';
874  }
875  }
876  Info<< endl;
877  }
878  }
879  Info<< endl;
880 }
881 
882 
883 // * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
884 
886 (
887  const word& surfName
888 ) const
889 {
890  const label surfI = findSurfaceID(surfName);
891 
892  if (surfI < 0)
893  {
895  << "Surface named " << surfName << " not found." << nl
896  << "Available surface names: " << names_ << endl
897  << abort(FatalError);
898  }
899 
900  return operator[](surfI);
901 }
902 
903 
905 (
906  const word& surfName
907 )
908 {
909  const label surfI = findSurfaceID(surfName);
910 
911  if (surfI < 0)
912  {
914  << "Surface named " << surfName << " not found." << nl
915  << "Available surface names: " << names_ << endl
916  << abort(FatalError);
917  }
918 
919  return operator[](surfI);
920 }
921 
922 
923 // ************************************************************************* //
static bool checkOrientation(const PrimitivePatch< FaceList, PointField > &, const bool report=false, labelHashSet *marked=0)
Check for orientation issues.
dictionary dict
autoPtr< IOobject > clone() const
Clone.
Definition: IOobject.H:268
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
const word & name() const
Return name.
Definition: IOobject.H:295
A class for handling file names.
Definition: fileName.H:79
virtual const wordList & regions() const =0
Names of regions.
label checkTopology(const bool report) const
All topological checks. Return number of failed checks.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
edgeList edges() const
Return edges in face point ordering,.
Definition: face.C:769
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static void findNearest(const PtrList< searchableSurface > &, const labelList &surfacesToTest, const pointField &, const scalarField &nearestDistSqr, labelList &surfaces, List< pointIndexHit > &)
Find nearest. Return -1 (and a miss()) or surface and nearest.
void findAnyIntersection(const pointField &start, const pointField &end, labelList &surfaces, List< pointIndexHit > &) const
Find any intersection. Return hit point information and.
virtual tmp< pointField > points() const
Get the points that define the surface.
Output to file stream.
Definition: OFstream.H:82
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
label findSurfaceID(const word &name) const
Find index of surface. Return -1 if not found.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
iterator end()
Return an iterator to end traversing the UPtrList.
Definition: UPtrListI.H:300
Base class of (analytical or triangulated) surface. Encapsulates all the search routines. WIP.
static void findAllIntersections(const PtrList< searchableSurface > &, const labelList &surfacesToTest, const pointField &start, const pointField &end, labelListList &surfaces, List< List< pointIndexHit >> &surfaceHits)
Find all intersections in order from start to end. Returns for.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const =0
Return any intersection on segment from start to end.
boundBox bounds() const
Calculate bounding box.
void findNearestIntersection(const pointField &start, const pointField &end, labelList &surface1, List< pointIndexHit > &hit1, labelList &surface2, List< pointIndexHit > &hit2) const
label findSurfaceRegionID(const word &surfaceName, const word &regionName) const
bool isDict(const word &) const
Check if entry is a sub-dictionary.
Definition: dictionary.C:653
static label markZones(const PrimitivePatch< FaceList, PointField > &, const BoolListType &borderEdge, labelList &faceZone)
Size and fills faceZone with zone of face.
fileName path() const
Return complete path.
Definition: IOobject.C:390
bool checkNormalOrientation(const bool report) const
Are all (triangulated) surfaces consistent normal orientation.
Various functions to operate on Lists.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:699
IOoject and searching on triSurface.
friend Ostream & operator(Ostream &, const UPtrList< T > &)
Write UPtrList to Ostream.
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))
Holds list of sampling positions.
Definition: coordSet.H:49
void findNearest(const pointField &, const scalarField &nearestDistSqr, labelList &surfaces, List< pointIndexHit > &) const
Find nearest. Return -1 (and a miss()) or surface and nearest.
static void findNearestIntersection(const PtrList< searchableSurface > &allSurfaces, const labelList &surfacesToTest, const pointField &start, const pointField &end, labelList &surface1, List< pointIndexHit > &hit1, labelList &surface2, List< pointIndexHit > &hit2)
Find intersections of edge nearest to both endpoints.
bool checkIntersection(const scalar tol, const autoPtr< writer< scalar >> &, const bool report) const
Do surfaces self-intersect or intersect others.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
face triFace(3)
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:90
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
const Field< PointType > & points() const
Return reference to global points.
const labelListList & edgeFaces() const
Return edge-face addressing.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
Dynamically sized Field.
Definition: DynamicField.H:49
searchableSurfaces(const label)
Construct with length specified. Fill later.
virtual label globalSize() const
Range of global indices that can be returned.
Triangle with additional region number.
Definition: labelledTri.H:57
static void findAnyIntersection(const PtrList< searchableSurface > &, const labelList &surfacesToTest, const pointField &start, const pointField &end, labelList &surfaces, List< pointIndexHit > &)
Find any intersection. Return hit point information and.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
label nEdges() const
Return number of edges in patch.
static const char nl
Definition: Ostream.H:265
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
void findAllIntersections(const pointField &start, const pointField &end, labelListList &surfaces, List< List< pointIndexHit >> &) const
Find all intersections in order from start to end. Returns for.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label checkGeometry(const scalar maxRatio, const scalar tolerance, const autoPtr< writer< scalar >> &setWriter, const scalar minQuality, const bool report) const
All geometric checks. Return number of failed checks.
void setSize(const label)
Reset size of List.
Definition: List.C:281
const searchableSurface & operator[](const word &) const
Return const reference to searchableSurface by name.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
const wordList & names() const
static autoPtr< searchableSurface > New(const word &surfaceType, const IOobject &io, const dictionary &dict)
Return a reference to the selected searchableSurface.
scalar tolerance() const
Return tolerance to use in searches.
triangle< point, const point & > triPointRef
Definition: triPointRef.H:44
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
bool checkQuality(const scalar minQuality, const bool report) const
Check triangle quality.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
const boundBox & bounds() const
Return const reference to boundBox.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
virtual Ostream & write(const token &)=0
Write next token to stream.
bool checkClosed(const bool report) const
Are all surfaces closed and manifold.
Triangulated surface description with patch information.
Definition: triSurface.H:66
void writeStats(const List< wordList > &, Ostream &) const
Write some stats.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
bool checkSizes(const scalar maxRatio, const bool report) const
Are all bounding boxes of similar size.
virtual bool hasVolumeType() const =0
Whether supports volume type below.
static boundBox bounds(const PtrList< searchableSurface > &allSurfaces, const labelList &surfacesToTest)
Find the boundBox of the selected surfaces.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583