face.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 "face.H"
27 #include "triFace.H"
28 #include "triPointRef.H"
29 #include "mathematicalConstants.H"
30 #include "Swap.H"
31 #include "ConstCirculator.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 const char* const Foam::face::typeName = "face";
36 
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
39 
41 Foam::face::calcEdges(const pointField& points) const
42 {
43  tmp<vectorField> tedges(new vectorField(size()));
44  vectorField& edges = tedges.ref();
45 
46  forAll(*this, i)
47  {
48  label ni = fcIndex(i);
49 
50  point thisPt = points[operator[](i)];
51  point nextPt = points[operator[](ni)];
52 
53  vector vec(nextPt - thisPt);
54  vec /= Foam::mag(vec) + vSmall;
55 
56  edges[i] = vec;
57  }
58 
59  return tedges;
60 }
61 
62 
63 Foam::scalar Foam::face::edgeCos
64 (
65  const vectorField& edges,
66  const label index
67 ) const
68 {
69  label leftEdgeI = left(index);
70  label rightEdgeI = right(index);
71 
72  // Note negate on left edge to get correct left-pointing edge.
73  return -(edges[leftEdgeI] & edges[rightEdgeI]);
74 }
75 
76 
77 Foam::label Foam::face::mostConcaveAngle
78 (
79  const pointField& points,
80  const vectorField& edges,
81  scalar& maxAngle
82 ) const
83 {
84  vector a(area(points));
85 
86  label index = 0;
87  maxAngle = -great;
88 
89  forAll(edges, i)
90  {
91  label leftEdgeI = left(i);
92  label rightEdgeI = right(i);
93 
94  vector edgeNormal = edges[rightEdgeI] ^ edges[leftEdgeI];
95 
96  scalar edgeCos = edges[leftEdgeI] & edges[rightEdgeI];
97  scalar edgeAngle = acos(max(-1.0, min(1.0, edgeCos)));
98 
99  scalar angle;
100 
101  if ((edgeNormal & a) > 0)
102  {
103  // Concave angle.
104  angle = constant::mathematical::pi + edgeAngle;
105  }
106  else
107  {
108  // Convex angle. Note '-' to take into account that rightEdge
109  // and leftEdge are head-to-tail connected.
110  angle = constant::mathematical::pi - edgeAngle;
111  }
112 
113  if (angle > maxAngle)
114  {
115  maxAngle = angle;
116  index = i;
117  }
118  }
119 
120  return index;
121 }
122 
123 
124 Foam::label Foam::face::split
125 (
126  const face::splitMode mode,
127  const pointField& points,
128  label& triI,
129  label& quadI,
130  faceList& triFaces,
131  faceList& quadFaces
132 ) const
133 {
134  label oldIndices = (triI + quadI);
135 
136  if (size() <= 2)
137  {
139  << "Serious problem: asked to split a face with < 3 vertices"
140  << abort(FatalError);
141  }
142  if (size() == 3)
143  {
144  // Triangle. Just copy.
145  if (mode == COUNTTRIANGLE || mode == COUNTQUAD)
146  {
147  triI++;
148  }
149  else
150  {
151  triFaces[triI++] = *this;
152  }
153  }
154  else if (size() == 4)
155  {
156  if (mode == COUNTTRIANGLE)
157  {
158  triI += 2;
159  }
160  if (mode == COUNTQUAD)
161  {
162  quadI++;
163  }
164  else if (mode == SPLITTRIANGLE)
165  {
166  // Start at point with largest internal angle.
167  const vectorField edges(calcEdges(points));
168 
169  scalar minAngle;
170  label startIndex = mostConcaveAngle(points, edges, minAngle);
171 
172  label nextIndex = fcIndex(startIndex);
173  label splitIndex = fcIndex(nextIndex);
174 
175  // Create triangles
176  face triFace(3);
177  triFace[0] = operator[](startIndex);
178  triFace[1] = operator[](nextIndex);
179  triFace[2] = operator[](splitIndex);
180 
181  triFaces[triI++] = triFace;
182 
183  triFace[0] = operator[](splitIndex);
184  triFace[1] = operator[](fcIndex(splitIndex));
185  triFace[2] = operator[](startIndex);
186 
187  triFaces[triI++] = triFace;
188  }
189  else
190  {
191  quadFaces[quadI++] = *this;
192  }
193  }
194  else
195  {
196  // General case. Like quad: search for largest internal angle.
197 
198  const vectorField edges(calcEdges(points));
199 
200  scalar minAngle = 1;
201  label startIndex = mostConcaveAngle(points, edges, minAngle);
202 
203  scalar bisectAngle = minAngle/2;
204  vector rightEdge = edges[right(startIndex)];
205 
206  //
207  // Look for opposite point which as close as possible bisects angle
208  //
209 
210  // split candidate starts two points away.
211  label index = fcIndex(fcIndex(startIndex));
212 
213  label minIndex = index;
214  scalar minDiff = constant::mathematical::pi;
215 
216  for (label i = 0; i < size() - 3; i++)
217  {
218  vector splitEdge
219  (
220  points[operator[](index)]
221  - points[operator[](startIndex)]
222  );
223  splitEdge /= Foam::mag(splitEdge) + vSmall;
224 
225  const scalar splitCos = splitEdge & rightEdge;
226  const scalar splitAngle = acos(max(-1.0, min(1.0, splitCos)));
227  const scalar angleDiff = fabs(splitAngle - bisectAngle);
228 
229  if (angleDiff < minDiff)
230  {
231  minDiff = angleDiff;
232  minIndex = index;
233  }
234 
235  // Go to next candidate
236  index = fcIndex(index);
237  }
238 
239 
240  // Split into two subshapes.
241  // face1: startIndex to minIndex
242  // face2: minIndex to startIndex
243 
244  // Get sizes of the two subshapes
245  label diff = 0;
246  if (minIndex > startIndex)
247  {
248  diff = minIndex - startIndex;
249  }
250  else
251  {
252  // Folded around
253  diff = minIndex + size() - startIndex;
254  }
255 
256  label nPoints1 = diff + 1;
257  label nPoints2 = size() - diff + 1;
258 
259  // Collect face1 points
260  face face1(nPoints1);
261 
262  index = startIndex;
263  for (label i = 0; i < nPoints1; i++)
264  {
265  face1[i] = operator[](index);
266  index = fcIndex(index);
267  }
268 
269  // Collect face2 points
270  face face2(nPoints2);
271 
272  index = minIndex;
273  for (label i = 0; i < nPoints2; i++)
274  {
275  face2[i] = operator[](index);
276  index = fcIndex(index);
277  }
278 
279  // Split faces
280  face1.split(mode, points, triI, quadI, triFaces, quadFaces);
281  face2.split(mode, points, triI, quadI, triFaces, quadFaces);
282  }
283 
284  return (triI + quadI - oldIndices);
285 }
286 
287 
288 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
289 
291 :
292  labelList(f)
293 {}
294 
295 
296 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
297 
298 int Foam::face::compare(const face& a, const face& b)
299 {
300  // Basic rule: we assume that the sequence of labels in each list
301  // will be circular in the same order (but not necessarily in the
302  // same direction or from the same starting point).
303 
304  // Trivial reject: faces are different size
305  label sizeA = a.size();
306  label sizeB = b.size();
307 
308  if (sizeA != sizeB || sizeA == 0)
309  {
310  return 0;
311  }
312  else if (sizeA == 1)
313  {
314  if (a[0] == b[0])
315  {
316  return 1;
317  }
318  else
319  {
320  return 0;
321  }
322  }
323 
324  ConstCirculator<face> aCirc(a);
325  ConstCirculator<face> bCirc(b);
326 
327  // Rotate face b until its element matches the starting element of face a.
328  do
329  {
330  if (aCirc() == bCirc())
331  {
332  // Set bCirc fulcrum to its iterator and increment the iterators
333  bCirc.setFulcrumToIterator();
334  ++aCirc;
335  ++bCirc;
336 
337  break;
338  }
340 
341  // If the circulator has stopped then faces a and b do not share a matching
342  // point. Doesn't work on matching, single element face.
343  if (!bCirc.circulate())
344  {
345  return 0;
346  }
347 
348  // Look forwards around the faces for a match
349  do
350  {
351  if (aCirc() != bCirc())
352  {
353  break;
354  }
355  }
356  while
357  (
360  );
361 
362  // If the circulator has stopped then faces a and b matched.
363  if (!aCirc.circulate())
364  {
365  return 1;
366  }
367  else
368  {
369  // Reset the circulators back to their fulcrum
370  aCirc.setIteratorToFulcrum();
371  bCirc.setIteratorToFulcrum();
372  ++aCirc;
373  --bCirc;
374  }
375 
376  // Look backwards around the faces for a match
377  do
378  {
379  if (aCirc() != bCirc())
380  {
381  break;
382  }
383  }
384  while
385  (
388  );
389 
390  // If the circulator has stopped then faces a and b matched.
391  if (!aCirc.circulate())
392  {
393  return -1;
394  }
395 
396  return 0;
397 }
398 
399 
400 bool Foam::face::sameVertices(const face& a, const face& b)
401 {
402  label sizeA = a.size();
403  label sizeB = b.size();
404 
405  // Trivial reject: faces are different size
406  if (sizeA != sizeB)
407  {
408  return false;
409  }
410  // Check faces with a single vertex
411  else if (sizeA == 1)
412  {
413  if (a[0] == b[0])
414  {
415  return true;
416  }
417  else
418  {
419  return false;
420  }
421  }
422 
423  forAll(a, i)
424  {
425  // Count occurrences of a[i] in a
426  label aOcc = 0;
427  forAll(a, j)
428  {
429  if (a[i] == a[j]) aOcc++;
430  }
431 
432  // Count occurrences of a[i] in b
433  label bOcc = 0;
434  forAll(b, j)
435  {
436  if (a[i] == b[j]) bOcc++;
437  }
438 
439  // Check if occurrences of a[i] in a and b are the same
440  if (aOcc != bOcc) return false;
441  }
442 
443  return true;
444 }
445 
446 
447 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
448 
450 {
451  if (size() > 1)
452  {
453  label ci = 0;
454  for (label i=1; i<size(); i++)
455  {
456  if (operator[](i) != operator[](ci))
457  {
458  operator[](++ci) = operator[](i);
459  }
460  }
461 
462  if (operator[](ci) != operator[](0))
463  {
464  ci++;
465  }
466 
467  setSize(ci);
468  }
469 
470  return size();
471 }
472 
473 
475 {
476  const label n = size();
477 
478  if (n > 2)
479  {
480  for (label i=1; i < (n+1)/2; ++i)
481  {
482  Swap(operator[](i), operator[](n-i));
483  }
484  }
485 }
486 
487 
489 {
490  // Calculate the centre by breaking the face into triangles and
491  // area-weighted averaging their centres
492 
493  const label nPoints = size();
494 
495  // If the face is a triangle, do a direct calculation
496  if (nPoints == 3)
497  {
498  return
499  (1.0/3.0)
500  *(
501  points[operator[](0)]
502  + points[operator[](1)]
503  + points[operator[](2)]
504  );
505  }
506 
507 
508  point centrePoint = Zero;
509  for (label pI=0; pI<nPoints; ++pI)
510  {
511  centrePoint += points[operator[](pI)];
512  }
513  centrePoint /= nPoints;
514 
515  scalar sumA = 0;
516  vector sumAc = Zero;
517 
518  for (label pI=0; pI<nPoints; ++pI)
519  {
520  const point& nextPoint = points[operator[]((pI + 1) % nPoints)];
521 
522  // Calculate 3*triangle centre
523  const vector ttc
524  (
525  points[operator[](pI)]
526  + nextPoint
527  + centrePoint
528  );
529 
530  // Calculate 2*triangle area
531  const scalar ta = Foam::mag
532  (
533  (points[operator[](pI)] - centrePoint)
534  ^ (nextPoint - centrePoint)
535  );
536 
537  sumA += ta;
538  sumAc += ta*ttc;
539  }
540 
541  if (sumA > vSmall)
542  {
543  return sumAc/(3.0*sumA);
544  }
545  else
546  {
547  return centrePoint;
548  }
549 }
550 
551 
553 {
554  const label nPoints = size();
555 
556  // Calculate the area by summing the face triangle areas.
557  // Changed to deal with small concavity by using a central decomposition
558  //
559 
560  // If the face is a triangle, do a direct calculation to avoid round-off
561  // error-related problems
562  //
563  if (nPoints == 3)
564  {
565  return triPointRef
566  (
567  p[operator[](0)],
568  p[operator[](1)],
569  p[operator[](2)]
570  ).area();
571  }
572 
573  label pI;
574 
575  point centrePoint = Zero;
576  for (pI = 0; pI < nPoints; ++pI)
577  {
578  centrePoint += p[operator[](pI)];
579  }
580  centrePoint /= nPoints;
581 
582  vector a = Zero;
583 
584  point nextPoint = centrePoint;
585 
586  for (pI = 0; pI < nPoints; ++pI)
587  {
588  if (pI < nPoints - 1)
589  {
590  nextPoint = p[operator[](pI + 1)];
591  }
592  else
593  {
594  nextPoint = p[operator[](0)];
595  }
596 
597  // Note: for best accuracy, centre point always comes last
598  //
599  a += triPointRef
600  (
601  p[operator[](pI)],
602  nextPoint,
603  centrePoint
604  ).area();
605  }
606 
607  return a;
608 }
609 
610 
612 {
613  const vector a = area(points);
614  const scalar maga = Foam::mag(a);
615  return maga > 0 ? a/maga : Zero;
616 }
617 
618 
620 {
621  // Reverse the label list and return
622  // The starting points of the original and reverse face are identical.
623 
624  const labelList& f = *this;
625  labelList newList(size());
626 
627  newList[0] = f[0];
628 
629  for (label pointi = 1; pointi < newList.size(); pointi++)
630  {
631  newList[pointi] = f[size() - pointi];
632  }
633 
634  return face(move(newList));
635 }
636 
637 
639 {
640  const labelList& f = *this;
641 
642  forAll(f, localIdx)
643  {
644  if (f[localIdx] == globalIndex)
645  {
646  return localIdx;
647  }
648  }
649 
650  return -1;
651 }
652 
653 
654 Foam::scalar Foam::face::sweptVol
655 (
656  const pointField& oldPoints,
657  const pointField& newPoints
658 ) const
659 {
660  // This Optimization causes a small discrepancy between the swept-volume of
661  // opposite faces of complex cells with triangular faces opposing polygons.
662  // It could be used without problem for tetrahedral cells
663  // if (size() == 3)
664  // {
665  // return
666  // (
667  // triPointRef
668  // (
669  // oldPoints[operator[](0)],
670  // oldPoints[operator[](1)],
671  // oldPoints[operator[](2)]
672  // ).sweptVol
673  // (
674  // triPointRef
675  // (
676  // newPoints[operator[](0)],
677  // newPoints[operator[](1)],
678  // newPoints[operator[](2)]
679  // )
680  // )
681  // );
682  // }
683 
684  scalar sv = 0;
685 
686  // Calculate the swept volume by breaking the face into triangles and
687  // summing their swept volumes.
688  // Changed to deal with small concavity by using a central decomposition
689 
690  point centreOldPoint = centre(oldPoints);
691  point centreNewPoint = centre(newPoints);
692 
693  label nPoints = size();
694 
695  for (label pi=0; pi<nPoints-1; ++pi)
696  {
697  // Note: for best accuracy, centre point always comes last
698  sv += triPointRef
699  (
700  centreOldPoint,
701  oldPoints[operator[](pi)],
702  oldPoints[operator[](pi + 1)]
703  ).sweptVol
704  (
706  (
707  centreNewPoint,
708  newPoints[operator[](pi)],
709  newPoints[operator[](pi + 1)]
710  )
711  );
712  }
713 
714  sv += triPointRef
715  (
716  centreOldPoint,
717  oldPoints[operator[](nPoints-1)],
718  oldPoints[operator[](0)]
719  ).sweptVol
720  (
722  (
723  centreNewPoint,
724  newPoints[operator[](nPoints-1)],
725  newPoints[operator[](0)]
726  )
727  );
728 
729  return sv;
730 }
731 
732 
734 (
735  const pointField& p,
736  const point& refPt,
737  scalar density
738 ) const
739 {
740  // If the face is a triangle, do a direct calculation
741  if (size() == 3)
742  {
743  return triPointRef
744  (
745  p[operator[](0)],
746  p[operator[](1)],
747  p[operator[](2)]
748  ).inertia(refPt, density);
749  }
750 
751  const point ctr = centre(p);
752 
753  tensor J = Zero;
754 
755  forAll(*this, i)
756  {
757  J += triPointRef
758  (
759  p[operator[](i)],
760  p[operator[](fcIndex(i))],
761  ctr
762  ).inertia(refPt, density);
763  }
764 
765  return J;
766 }
767 
768 
770 {
771  const labelList& points = *this;
772 
773  edgeList e(points.size());
774 
775  for (label pointi = 0; pointi < points.size() - 1; ++pointi)
776  {
777  e[pointi] = edge(points[pointi], points[pointi + 1]);
778  }
779 
780  // Add last edge
781  e.last() = edge(points.last(), points[0]);
782 
783  return e;
784 }
785 
786 
788 {
789  forAll(*this, i)
790  {
791  if (operator[](i) == e.start())
792  {
793  if (operator[](rcIndex(i)) == e.end())
794  {
795  // Reverse direction
796  return -1;
797  }
798  else if (operator[](fcIndex(i)) == e.end())
799  {
800  // Forward direction
801  return 1;
802  }
803 
804  // No match
805  return 0;
806  }
807  else if (operator[](i) == e.end())
808  {
809  if (operator[](rcIndex(i)) == e.start())
810  {
811  // Forward direction
812  return 1;
813  }
814  else if (operator[](fcIndex(i)) == e.start())
815  {
816  // Reverse direction
817  return -1;
818  }
819 
820  // No match
821  return 0;
822  }
823  }
824 
825  // Not found
826  return 0;
827 }
828 
829 
831 {
832  return nTriangles();
833 }
834 
835 
837 (
838  const pointField& points,
839  label& triI,
840  faceList& triFaces
841 ) const
842 {
843  label quadI = 0;
844  faceList quadFaces;
845 
846  return split(SPLITTRIANGLE, points, triI, quadI, triFaces, quadFaces);
847 }
848 
849 
851 (
852  const pointField& points,
853  label& triI,
854  label& quadI
855 ) const
856 {
857  faceList triFaces;
858  faceList quadFaces;
859 
860  return split(COUNTQUAD, points, triI, quadI, triFaces, quadFaces);
861 }
862 
863 
865 (
866  const pointField& points,
867  label& triI,
868  label& quadI,
869  faceList& triFaces,
870  faceList& quadFaces
871 ) const
872 {
873  return split(SPLITQUAD, points, triI, quadI, triFaces, quadFaces);
874 }
875 
876 
878 {
879  const edgeList& eds = f.edges();
880 
881  label longestEdgeI = -1;
882  scalar longestEdgeLength = -small;
883 
884  forAll(eds, edI)
885  {
886  scalar edgeLength = eds[edI].mag(pts);
887 
888  if (edgeLength > longestEdgeLength)
889  {
890  longestEdgeI = edI;
891  longestEdgeLength = edgeLength;
892  }
893  }
894 
895  return longestEdgeI;
896 }
897 
898 
899 // ************************************************************************* //
label nTrianglesQuads(const pointField &points, label &nTris, label &nQuads) const
Number of triangles and quads after splitting.
Definition: face.C:851
dimensionedScalar acos(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A triangle primitive used to calculate face areas and swept volumes.
Definition: triangle.H:57
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
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
edgeList edges() const
Return edges in face point ordering,.
Definition: face.C:769
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
T & operator[](const label)
Return element of UList.
Definition: UListI.H:167
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
face()
Construct null.
Definition: faceI.H:44
label rcIndex(const label i) const
Return the reverse circular index, i.e. the previous index.
Definition: UListI.H:65
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
Definition: UListI.H:58
List< face > faceList
Definition: faceListFwd.H:43
tensor inertia(const pointField &, const point &refPt=vector::zero, scalar density=1.0) const
Return the inertia tensor, with optional reference.
Definition: face.C:734
point centre(const pointField &) const
Centre point of face.
Definition: face.C:488
label which(const label globalIndex) const
Navigation through face vertices.
Definition: face.C:638
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
void setIteratorToFulcrum()
Set the iterator to the current position of the fulcrum.
label triangles(const pointField &points, label &triI, faceList &triFaces) const
Split into triangles using existing points.
Definition: face.C:837
vector normal(const pointField &) const
Return unit normal.
Definition: face.C:611
label collapse()
Collapse face by removing duplicate point labels.
Definition: face.C:449
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
face reverseFace() const
Return face with reverse direction.
Definition: face.C:619
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
void Swap(T &a, T &b)
Definition: Swap.H:43
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
int edgeDirection(const edge &) const
Return the edge direction on the face.
Definition: face.C:787
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
label nPoints
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:68
face triFace(3)
static int compare(const face &, const face &)
Compare faces.
Definition: face.C:298
label trianglesQuads(const pointField &points, label &triI, label &quadI, faceList &triFaces, faceList &quadFaces) const
Split into triangles and quads.
Definition: face.C:865
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
vector area(const pointField &) const
Return vector area.
Definition: face.C:552
static const char *const typeName
Definition: face.H:143
labelList f(nPoints)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
scalar sweptVol(const pointField &oldPoints, const pointField &newPoints) const
Return the volume swept out by the face when its points move.
Definition: face.C:655
label longestEdge(const face &f, const pointField &pts)
Find the longest edge on a face. Face point labels index into pts.
Definition: face.C:877
static bool sameVertices(const face &, const face &)
Return true if the faces have the same vertices.
Definition: face.C:400
void setSize(const label)
Reset size of List.
Definition: List.C:281
bool circulate(const CirculatorBase::direction dir=CirculatorBase::direction::none)
Circulate around the list in the given direction.
vector point
Point is a vector.
Definition: point.H:41
label end() const
Return end vertex label.
Definition: edgeI.H:92
Swap its arguments.
void flip()
Flip the face in-place.
Definition: face.C:474
triangle< point, const point & > triPointRef
Definition: triPointRef.H:44
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
Field< vector > vectorField
Specialisation of Field<T> for vector.
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
Walks over a container as if it were circular. The container must have the following members defined:...
label nTriangles() const
Number of triangles after splitting.
Definition: faceI.H:130
A class for managing temporary objects.
Definition: PtrList.H:53
T & last()
Return the last element of the list.
Definition: UListI.H:128
void setFulcrumToIterator()
Set the fulcrum to the current position of the iterator.
label size() const
Return the number of elements in the UList.
Definition: ListI.H:170
label start() const
Return start vertex label.
Definition: edgeI.H:81