controlMeshRefinement.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) 2013-2015 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 "controlMeshRefinement.H"
28 #include "OFstream.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(controlMeshRefinement, 0);
35 }
36 
37 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
38 
39 
40 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41 
42 Foam::scalar Foam::controlMeshRefinement::calcFirstDerivative
43 (
44  const Foam::point& a,
45  const scalar& cellSizeA,
46  const Foam::point& b,
47  const scalar& cellSizeB
48 ) const
49 {
50  return (cellSizeA - cellSizeB)/mag(a - b);
51 }
52 
53 
54 //Foam::scalar Foam::controlMeshRefinement::calcSecondDerivative
55 //(
56 // const Foam::point& a,
57 // const scalar& cellSizeA,
58 // const Foam::point& midPoint,
59 // const scalar& cellSizeMid,
60 // const Foam::point& b,
61 // const scalar& cellSizeB
62 //) const
63 //{
64 // return (cellSizeA - 2*cellSizeMid + cellSizeB)/magSqr((a - b)/2);
65 //}
66 
67 
68 bool Foam::controlMeshRefinement::detectEdge
69 (
70  const Foam::point& startPt,
71  const Foam::point& endPt,
72  pointHit& pointFound,
73  const scalar tolSqr,
74  const scalar secondDerivTolSqr
75 ) const
76 {
77  Foam::point a(startPt);
78  Foam::point b(endPt);
79 
80  Foam::point midPoint = (a + b)/2.0;
81 
82  label nIterations = 0;
83 
84  while (true)
85  {
86  nIterations++;
87 
88  if
89  (
90  magSqr(a - b) < tolSqr
91  )
92  {
93  pointFound.setPoint(midPoint);
94  pointFound.setHit();
95 
96  return true;
97  }
98 
99  // Split into two regions
100 
101  scalar cellSizeA = sizeControls_.cellSize(a);
102  scalar cellSizeB = sizeControls_.cellSize(b);
103 
104 // if (magSqr(cellSizeA - cellSizeB) < 1e-6)
105 // {
106 // return false;
107 // }
108 
109  scalar cellSizeMid = sizeControls_.cellSize(midPoint);
110 
111  // Region 1
112  Foam::point midPoint1 = (a + midPoint)/2.0;
113  const scalar cellSizeMid1 = sizeControls_.cellSize(midPoint1);
114 
115 // scalar firstDerivative1 =
116 // calcFirstDerivative(cellSizeA, cellSizeMid);
117 
118  scalar secondDerivative1 =
119  calcSecondDerivative
120  (
121  a,
122  cellSizeA,
123  midPoint1,
124  cellSizeMid1,
125  midPoint,
126  cellSizeMid
127  );
128 
129  // Region 2
130  Foam::point midPoint2 = (midPoint + b)/2.0;
131  const scalar cellSizeMid2 = sizeControls_.cellSize(midPoint2);
132 
133 // scalar firstDerivative2 =
134 // calcFirstDerivative(f, cellSizeMid, cellSizeB);
135 
136  scalar secondDerivative2 =
137  calcSecondDerivative
138  (
139  midPoint,
140  cellSizeMid,
141  midPoint2,
142  cellSizeMid2,
143  b,
144  cellSizeB
145  );
146 
147  // Neither region appears to have an inflection
148  // To be sure should use higher order derivatives
149  if
150  (
151  magSqr(secondDerivative1) < secondDerivTolSqr
152  && magSqr(secondDerivative2) < secondDerivTolSqr
153  )
154  {
155  return false;
156  }
157 
158  // Pick region with greatest second derivative
159  if (magSqr(secondDerivative1) > magSqr(secondDerivative2))
160  {
161  b = midPoint;
162  midPoint = midPoint1;
163  }
164  else
165  {
166  a = midPoint;
167  midPoint = midPoint2;
168  }
169  }
170 }
171 
172 
173 Foam::pointHit Foam::controlMeshRefinement::findDiscontinuities
174 (
175  const linePointRef& l
176 ) const
177 {
179 
180  const scalar tolSqr = sqr(1e-3);
181  const scalar secondDerivTolSqr = sqr(1e-3);
182 
183  detectEdge
184  (
185  l.start(),
186  l.end(),
187  p,
188  tolSqr,
189  secondDerivTolSqr
190  );
191 
192  return p;
193 }
194 
195 
196 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
197 
198 
199 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
200 
201 Foam::controlMeshRefinement::controlMeshRefinement
202 (
203  cellShapeControl& shapeController
204 )
205 :
206  shapeController_(shapeController),
207  mesh_(shapeController.shapeControlMesh()),
208  sizeControls_(shapeController.sizeAndAlignment()),
209  geometryToConformTo_(sizeControls_.geometryToConformTo())
210 {}
211 
212 
213 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
214 
216 {}
217 
218 
219 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
220 
222 (
223  const autoPtr<backgroundMeshDecomposition>& decomposition
224 )
225 {
226  if (shapeController_.shapeControlMesh().vertexCount() > 0)
227  {
228  // Mesh already populated.
229  Info<< "Cell size and alignment mesh already populated." << endl;
230  return;
231  }
232 
233  autoPtr<boundBox> overallBoundBox;
234 
235  // Need to pass in the background mesh decomposition so that can test if
236  // a point to insert is on the processor.
237  if (Pstream::parRun())
238  {
239 // overallBoundBox.set(new boundBox(decomposition().procBounds()));
240  }
241  else
242  {
243 // overallBoundBox.set
244 // (
245 // new boundBox(geometryToConformTo_.geometry().bounds())
246 // );
247 //
248 // mesh_.insertBoundingPoints
249 // (
250 // overallBoundBox(),
251 // sizeControls_
252 // );
253  }
254 
255  Map<label> priorityMap;
256 
257  const PtrList<cellSizeAndAlignmentControl>& controlFunctions =
258  sizeControls_.controlFunctions();
259 
260  forAll(controlFunctions, fI)
261  {
262  const cellSizeAndAlignmentControl& controlFunction =
263  controlFunctions[fI];
264 
265  const Switch& forceInsertion =
266  controlFunction.forceInitialPointInsertion();
267 
268  Info<< "Inserting points from " << controlFunction.name()
269  << " (" << controlFunction.type() << ")" << endl;
270  Info<< " Force insertion is " << forceInsertion.asText() << endl;
271 
272  pointField pts;
273  scalarField sizes;
274  triadField alignments;
275 
276  controlFunction.initialVertices(pts, sizes, alignments);
277 
278  Info<< " Got initial vertices list of size " << pts.size() << endl;
279 
280  List<Vb> vertices(pts.size());
281 
282  // Clip the minimum size
283  for (label vI = 0; vI < pts.size(); ++vI)
284  {
285  vertices[vI] = Vb(pts[vI], Vb::vtInternalNearBoundary);
286 
287  label maxPriority = -1;
288  scalar size = sizeControls_.cellSize(pts[vI], maxPriority);
289 
290  if (maxPriority > controlFunction.maxPriority())
291  {
292  vertices[vI].targetCellSize() = max
293  (
294  size,
295  shapeController_.minimumCellSize()
296  );
297  }
298 // else if (maxPriority == controlFunction.maxPriority())
299 // {
300 // vertices[vI].targetCellSize() = max
301 // (
302 // min(sizes[vI], size),
303 // shapeController_.minimumCellSize()
304 // );
305 // }
306  else
307  {
308  vertices[vI].targetCellSize() = max
309  (
310  sizes[vI],
311  shapeController_.minimumCellSize()
312  );
313  }
314 
315  vertices[vI].alignment() = alignments[vI];
316  }
317 
318  Info<< " Clipped minimum size" << endl;
319 
320  pts.clear();
321  sizes.clear();
322  alignments.clear();
323 
324  PackedBoolList keepVertex(vertices.size(), true);
325 
326  forAll(vertices, vI)
327  {
328  bool keep = true;
329 
330  pointFromPoint pt = topoint(vertices[vI].point());
331 
332  if (Pstream::parRun())
333  {
334  keep = decomposition().positionOnThisProcessor(pt);
335  }
336 
337  if (keep && geometryToConformTo_.wellOutside(pt, SMALL))
338  {
339  keep = false;
340  }
341 
342  if (!keep)
343  {
344  keepVertex[vI] = false;
345  }
346  }
347 
348  inplaceSubset(keepVertex, vertices);
349 
350  const label preInsertedSize = mesh_.number_of_vertices();
351 
352  Info<< " Check sizes" << endl;
353 
354  forAll(vertices, vI)
355  {
356  bool insertPoint = false;
357 
358  pointFromPoint pt(topoint(vertices[vI].point()));
359 
360  if
361  (
362  mesh_.dimension() < 3
363  || mesh_.is_infinite
364  (
365  mesh_.locate(vertices[vI].point())
366  )
367  )
368  {
369  insertPoint = true;
370  }
371 
372  const scalar interpolatedCellSize = shapeController_.cellSize(pt);
373  const triad interpolatedAlignment =
374  shapeController_.cellAlignment(pt);
375  const scalar calculatedCellSize = vertices[vI].targetCellSize();
376  const triad calculatedAlignment = vertices[vI].alignment();
377 
378  if (debug)
379  {
380  Info<< "Point = " << pt << nl
381  << " Size(interp) = " << interpolatedCellSize << nl
382  << " Size(calc) = " << calculatedCellSize << nl
383  << " Align(interp) = " << interpolatedAlignment << nl
384  << " Align(calc) = " << calculatedAlignment << nl
385  << endl;
386  }
387 
388  const scalar sizeDiff =
389  mag(interpolatedCellSize - calculatedCellSize);
390  const scalar alignmentDiff =
391  diff(interpolatedAlignment, calculatedAlignment);
392 
393  if (debug)
394  {
395  Info<< " size difference = " << sizeDiff << nl
396  << ", alignment difference = " << alignmentDiff << endl;
397  }
398 
399  // TODO: Also need to base it on the alignments
400  if
401  (
402  sizeDiff/interpolatedCellSize > 0.1
403  || alignmentDiff > 0.15
404  )
405  {
406  insertPoint = true;
407  }
408 
409  if (forceInsertion || insertPoint)
410  {
411  const label oldSize = mesh_.vertexCount();
412 
413  cellShapeControlMesh::Vertex_handle insertedVert = mesh_.insert
414  (
415  pt,
416  calculatedCellSize,
417  vertices[vI].alignment(),
419  );
420 
421  if (oldSize == mesh_.vertexCount() - 1)
422  {
423  priorityMap.insert
424  (
425  insertedVert->index(),
426  controlFunction.maxPriority()
427  );
428  }
429  }
430  }
431 
432  //mesh_.rangeInsertWithInfo(vertices.begin(), vertices.end());
433 
434  Info<< " Inserted "
435  << returnReduce
436  (
437  label(mesh_.number_of_vertices()) - preInsertedSize,
438  sumOp<label>()
439  )
440  << "/" << returnReduce(vertices.size(), sumOp<label>())
441  << endl;
442  }
443 
444 
445 
446  forAll(controlFunctions, fI)
447  {
448  const cellSizeAndAlignmentControl& controlFunction =
449  controlFunctions[fI];
450 
451  const Switch& forceInsertion =
452  controlFunction.forceInitialPointInsertion();
453 
454  Info<< "Inserting points from " << controlFunction.name()
455  << " (" << controlFunction.type() << ")" << endl;
456  Info<< " Force insertion is " << forceInsertion.asText() << endl;
457 
458  DynamicList<Foam::point> extraPts;
459  DynamicList<scalar> extraSizes;
460 
461  controlFunction.cellSizeFunctionVertices(extraPts, extraSizes);
462 
463  List<Vb> vertices(extraPts.size());
464 
465  // Clip the minimum size
466  for (label vI = 0; vI < extraPts.size(); ++vI)
467  {
468  vertices[vI] = Vb(extraPts[vI], Vb::vtUnassigned);
469 
470  label maxPriority = -1;
471  scalar size = sizeControls_.cellSize(extraPts[vI], maxPriority);
472 
473  if (maxPriority > controlFunction.maxPriority())
474  {
475  vertices[vI].targetCellSize() = max
476  (
477  size,
478  shapeController_.minimumCellSize()
479  );
480  }
481  else if (maxPriority == controlFunction.maxPriority())
482  {
483  vertices[vI].targetCellSize() = max
484  (
485  min(extraSizes[vI], size),
486  shapeController_.minimumCellSize()
487  );
488  }
489  else
490  {
491  vertices[vI].targetCellSize() = max
492  (
493  extraSizes[vI],
494  shapeController_.minimumCellSize()
495  );
496  }
497  }
498 
499  PackedBoolList keepVertex(vertices.size(), true);
500 
501  forAll(vertices, vI)
502  {
503  bool keep = true;
504 
505  pointFromPoint pt = topoint(vertices[vI].point());
506 
507  if (Pstream::parRun())
508  {
509  keep = decomposition().positionOnThisProcessor(pt);
510  }
511 
512  if (keep && geometryToConformTo_.wellOutside(pt, SMALL))
513  {
514  keep = false;
515  }
516 
517  if (!keep)
518  {
519  keepVertex[vI] = false;
520  }
521  }
522 
523  inplaceSubset(keepVertex, vertices);
524 
525  const label preInsertedSize = mesh_.number_of_vertices();
526 
527  forAll(vertices, vI)
528  {
529  bool insertPoint = false;
530 
531  pointFromPoint pt(topoint(vertices[vI].point()));
532 
533  if
534  (
535  mesh_.dimension() < 3
536  || mesh_.is_infinite
537  (
538  mesh_.locate(vertices[vI].point())
539  )
540  )
541  {
542  insertPoint = true;
543  }
544 
545  const scalar interpolatedCellSize = shapeController_.cellSize(pt);
546  const scalar calculatedCellSize = vertices[vI].targetCellSize();
547 
548  if (debug)
549  {
550  Info<< "Point = " << pt << nl
551  << " Size(interp) = " << interpolatedCellSize << nl
552  << " Size(calc) = " << calculatedCellSize << nl
553  << endl;
554  }
555 
556  const scalar sizeDiff =
557  mag(interpolatedCellSize - calculatedCellSize);
558 
559  if (debug)
560  {
561  Info<< " size difference = " << sizeDiff << endl;
562  }
563 
564  // TODO: Also need to base it on the alignments
565  if (sizeDiff/interpolatedCellSize > 0.1)
566  {
567  insertPoint = true;
568  }
569 
570  if (forceInsertion || insertPoint)
571  {
572  // Check the priority
573 
574 // cellShapeControlMesh::Cell_handle ch =
575 // mesh_.locate(toPoint<cellShapeControlMesh::Point>(pt));
576 
577 // if (mesh_.is_infinite(ch))
578 // {
579 // continue;
580 // }
581 
582 // const label newPtPriority = controlFunction.maxPriority();
583 
584 // label highestPriority = -1;
585 // for (label cI = 0; cI < 4; ++cI)
586 // {
587 // if (mesh_.is_infinite(ch->vertex(cI)))
588 // {
589 // continue;
590 // }
591 
592 // const label vertPriority =
593 // priorityMap[ch->vertex(cI)->index()];
594 
595 // if (vertPriority > highestPriority)
596 // {
597 // highestPriority = vertPriority;
598 // }
599 // }
600 
601 // if (newPtPriority >= highestPriority)
602 // {
603 // const label oldSize = mesh_.vertexCount();
604 //
605 // cellShapeControlMesh::Vertex_handle insertedVert =
606  mesh_.insert
607  (
608  pt,
609  calculatedCellSize,
610  vertices[vI].alignment(),
612  );
613 
614 // if (oldSize == mesh_.vertexCount() - 1)
615 // {
616 // priorityMap.insert
617 // (
618 // insertedVert->index(),
619 // newPtPriority
620 // );
621 // }
622 // }
623  }
624  }
625 
626  //mesh_.rangeInsertWithInfo(vertices.begin(), vertices.end());
627 
628  Info<< " Inserted extra points "
629  << returnReduce
630  (
631  label(mesh_.number_of_vertices()) - preInsertedSize,
632  sumOp<label>()
633  )
634  << "/" << returnReduce(vertices.size(), sumOp<label>())
635  << endl;
636  }
637 
638  // Change cell size function of bounding points to be consistent
639  // with their nearest neighbours
640 // for
641 // (
642 // CellSizeDelaunay::Finite_vertices_iterator vit =
643 // mesh_.finite_vertices_begin();
644 // vit != mesh_.finite_vertices_end();
645 // ++vit
646 // )
647 // {
648 // if (vit->uninitialised())
649 // {
650 // // Get its adjacent vertices
651 // std::list<CellSizeDelaunay::Vertex_handle> adjacentVertices;
652 //
653 // mesh_.adjacent_vertices
654 // (
655 // vit,
656 // std::back_inserter(adjacentVertices)
657 // );
658 //
659 // scalar totalCellSize = 0;
660 // label nVerts = 0;
661 //
662 // for
663 // (
664 // std::list<CellSizeDelaunay::Vertex_handle>::iterator avit =
665 // adjacentVertices.begin();
666 // avit != adjacentVertices.end();
667 // ++avit
668 // )
669 // {
670 // if (!(*avit)->uninitialised())
671 // {
672 // totalCellSize += (*avit)->targetCellSize();
673 // nVerts++;
674 // }
675 // }
676 //
677 // Pout<< "Changing " << vit->info();
678 //
679 // vit->targetCellSize() = totalCellSize/nVerts;
680 // vit->type() = Vb::vtInternalNearBoundary;
681 //
682 // Pout<< "to " << vit->info() << endl;
683 // }
684 // }
685 }
686 
687 
689 (
690  const autoPtr<backgroundMeshDecomposition>& decomposition
691 )
692 {
693  Info<< "Iterate over "
694  << returnReduce(label(mesh_.number_of_finite_edges()), sumOp<label>())
695  << " cell size mesh edges" << endl;
696 
697  DynamicList<Vb> verts(mesh_.number_of_vertices());
698 
699  label count = 0;
700 
701  for
702  (
703  CellSizeDelaunay::Finite_edges_iterator eit =
704  mesh_.finite_edges_begin();
705  eit != mesh_.finite_edges_end();
706  ++eit
707  )
708  {
709  if (count % 10000 == 0)
710  {
711  Info<< count << " edges, inserted " << verts.size()
712  << " Time = " << mesh_.time().elapsedCpuTime()
713  << endl;
714  }
715  count++;
716 
717  CellSizeDelaunay::Cell_handle c = eit->first;
718  CellSizeDelaunay::Vertex_handle vA = c->vertex(eit->second);
719  CellSizeDelaunay::Vertex_handle vB = c->vertex(eit->third);
720 
721  if
722  (
723  mesh_.is_infinite(vA)
724  || mesh_.is_infinite(vB)
725  || (vA->referred() && vB->referred())
726  || (vA->referred() && (vA->procIndex() > vB->procIndex()))
727  || (vB->referred() && (vB->procIndex() > vA->procIndex()))
728  )
729  {
730  continue;
731  }
732 
733  pointFromPoint ptA(topoint(vA->point()));
734  pointFromPoint ptB(topoint(vB->point()));
735 
736  linePointRef l(ptA, ptB);
737 
738  const pointHit hitPt = findDiscontinuities(l);
739 
740  if (hitPt.hit())
741  {
742  const Foam::point& pt = hitPt.hitPoint();
743 
744  if (!geometryToConformTo_.inside(pt))
745  {
746  continue;
747  }
748 
749  if (Pstream::parRun())
750  {
751  if (!decomposition().positionOnThisProcessor(pt))
752  {
753  continue;
754  }
755  }
756 
757  verts.append
758  (
759  Vb
760  (
761  toPoint(pt),
763  )
764  );
765 
766  verts.last().targetCellSize() = sizeControls_.cellSize(pt);
767  verts.last().alignment() = triad::unset;
768  }
769  }
770 
771  mesh_.insertPoints(verts, false);
772 
773  return verts.size();
774 }
775 
776 
777 // * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
778 
779 
780 // * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * * //
781 
782 
783 // * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
784 
785 
786 // ************************************************************************* //
void inplaceSubset(const UList< T > &select, const T &value, ListType &)
Inplace extract elements of List when select is a certain value.
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:394
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
Field< bool > wellOutside(const pointField &samplePts, const scalarField &testDistSqr) const
Check if point is outside surfaces to conform to by at least.
pointFromPoint topoint(const Point &P)
const double e
Elementary charge.
Definition: doubleFloat.H:78
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Field< bool > inside(const pointField &samplePts) const
Check if points are inside surfaces to conform to.
cellSizeAndAlignmentControl(const Time &runTime, const word &name, const dictionary &controlFunctionDict, const conformationSurfaces &geometryToConformTo, const scalar &defaultCellSize)
Construct from dictionary and references to conformalVoronoiMesh and.
An indexed form of CGAL::Triangulation_vertex_base_3<K> used to keep track of the Delaunay vertices i...
Definition: indexedVertex.H:51
void initialMeshPopulation(const autoPtr< backgroundMeshDecomposition > &decomposition)
CGAL::indexedVertex< K > Vb
label refineMesh(const autoPtr< backgroundMeshDecomposition > &decomposition)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
line< point, const point & > linePointRef
Line using referred points.
Definition: linePointRef.H:45
static const triad unset
Definition: triad.H:99
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
PointFrompoint toPoint(const Foam::point &p)
~controlMeshRefinement()
Destructor.
Field< triad > triadField
Specialisation of Field<T> for triad.
Definition: triadField.H:49
dimensioned< scalar > magSqr(const dimensioned< Type > &)
CellSizeDelaunay::Vertex_handle Vertex_handle
static direction size()
Return the number of elements in the VectorSpace = Ncmpts.
Definition: VectorSpaceI.H:93
static const char nl
Definition: Ostream.H:262
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
vector point
Point is a vector.
Definition: point.H:41
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
PointHit< point > pointHit
Definition: pointHit.H:41
volScalarField & p
Namespace for OpenFOAM.