pointEdgeCollapseI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2012-2016 OpenFOAM Foundation
9  Copyright (C) 2019-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "polyMesh.H"
30 #include "transform.H"
31 
32 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33 
34 template<class TrackingData>
35 inline bool Foam::pointEdgeCollapse::update
36 (
37  const pointEdgeCollapse& w2,
38  const scalar tol,
39  TrackingData& td
40 )
41 {
42  if (!w2.valid(td))
43  {
45  << "problem." << abort(FatalError);
46  }
47 
48  if (!valid(td))
49  {
50  operator=(w2);
51  return true;
52  }
53 
54  if (w2.collapseIndex_ == -1 || collapseIndex_ == -1)
55  {
56  // Not marked for collapse; only happens on edges.
57  return false;
58  }
59 
60  if (w2.collapsePriority_ < collapsePriority_)
61  {
62  return false;
63  }
64  else if (w2.collapsePriority_ > collapsePriority_)
65  {
66  operator=(w2);
67  return true;
68  }
69 
70  // Get overwritten by w2 if it has a higher priority
71  if (w2.collapseIndex_ < collapseIndex_)
72  {
73  operator=(w2);
74  return true;
75  }
76  else if (w2.collapseIndex_ == collapseIndex_)
77  {
78  bool identicalPoint = samePoint(w2.collapsePoint_);
79 
80  bool nearer = (magSqr(w2.collapsePoint_) < magSqr(collapsePoint_));
81 
82  if (nearer)
83  {
84  operator=(w2);
85  }
86 
87  if (identicalPoint)
88  {
89  return false;
90  }
91  else
92  {
93  return nearer;
94  }
95  }
96 
97  return false;
98 }
99 
100 
101 inline bool Foam::pointEdgeCollapse::samePoint(const point& pt) const
102 {
103  bool isLegal1 = (cmptMin(collapsePoint_) < 0.5*GREAT);
104  bool isLegal2 = (cmptMin(pt) < 0.5*GREAT);
105 
106  if (isLegal1 && isLegal2)
107  {
108  return mag(collapsePoint_ - pt) < 1e-9;//SMALL;
109  }
110  else
111  {
112  return isLegal1 == isLegal2;
113  }
114 }
115 
116 
117 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
118 
120 :
121  collapsePoint_(GREAT, GREAT, GREAT),
122  collapseIndex_(-2),
123  collapsePriority_(-2)
124 {}
125 
126 
128 (
129  const point& collapsePoint,
130  const label collapseIndex,
131  const label collapsePriority
132 )
133 :
134  collapsePoint_(collapsePoint),
135  collapseIndex_(collapseIndex),
136  collapsePriority_(collapsePriority)
137 {}
138 
139 
140 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
141 
142 template<class TrackingData>
143 inline bool Foam::pointEdgeCollapse::valid(TrackingData& td) const
144 {
145  return collapseIndex_ != -2;
146 }
147 
148 
149 template<class TrackingData>
151 (
152  const polyPatch& patch,
153  const label patchPointi,
154  const point& coord,
155  TrackingData& td
156 )
157 {
158  collapsePoint_ -= coord;
159 }
160 
161 
162 template<class TrackingData>
164 (
165  const tensor& rotTensor,
166  TrackingData& td
167 )
168 {
169  collapsePoint_ = Foam::transform(rotTensor, collapsePoint_);
170 }
171 
172 
173 // Update absolute geometric quantities. Note that distance (dist_)
174 // is not affected by leaving/entering domain.
175 template<class TrackingData>
177 (
178  const polyPatch& patch,
179  const label patchPointi,
180  const point& coord,
181  TrackingData& td
182 )
183 {
184  // back to absolute form
185  collapsePoint_ += coord;
186 }
187 
188 
189 // Update this with information from connected edge
190 template<class TrackingData>
192 (
193  const polyMesh& mesh,
194  const label pointi,
195  const label edgeI,
196  const pointEdgeCollapse& edgeInfo,
197  const scalar tol,
198  TrackingData& td
199 )
200 {
201  return update(edgeInfo, tol, td);
202 }
203 
204 
205 // Update this with new information on same point
206 template<class TrackingData>
208 (
209  const polyMesh& mesh,
210  const label pointi,
211  const pointEdgeCollapse& newPointInfo,
212  const scalar tol,
213  TrackingData& td
214 )
215 {
216  return update(newPointInfo, tol, td);
217 }
218 
219 
220 // Update this with new information on same point. No extra information.
221 template<class TrackingData>
223 (
224  const pointEdgeCollapse& newPointInfo,
225  const scalar tol,
226  TrackingData& td
227 )
228 {
229  return update(newPointInfo, tol, td);
230 }
231 
232 
233 // Update this with information from connected point
234 template<class TrackingData>
236 (
237  const polyMesh& mesh,
238  const label edgeI,
239  const label pointi,
240  const pointEdgeCollapse& pointInfo,
241  const scalar tol,
242  TrackingData& td
243 )
244 {
245  return update(pointInfo, tol, td);
246 }
247 
248 
249 template<class TrackingData>
251 (
252  const pointEdgeCollapse& rhs,
253  TrackingData& td
254 ) const
255 {
256  return operator==(rhs);
257 }
258 
259 
260 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
261 
262 inline bool Foam::pointEdgeCollapse::operator==
263 (
264  const pointEdgeCollapse& rhs
265 ) const
266 {
267  return
268  collapseIndex_ == rhs.collapseIndex_
269  && collapsePriority_ == rhs.collapsePriority_
270  && samePoint(rhs.collapsePoint_);
271 }
272 
273 
274 inline bool Foam::pointEdgeCollapse::operator!=
275 (
276  const pointEdgeCollapse& rhs
277 ) const
278 {
279  return !(*this == rhs);
280 }
281 
282 
283 // ************************************************************************* //
Foam::pointEdgeCollapse::updateEdge
bool updateEdge(const polyMesh &mesh, const label edgeI, const label pointi, const pointEdgeCollapse &pointInfo, const scalar tol, TrackingData &td)
Influence of point on edge.
Definition: pointEdgeCollapseI.H:236
Foam::Tensor< scalar >
Foam::pointEdgeCollapse::equal
bool equal(const pointEdgeCollapse &, TrackingData &) const
Test for equality, with TrackingData.
Definition: pointEdgeCollapseI.H:251
Foam::pointEdgeCollapse::transform
void transform(const tensor &rotTensor, TrackingData &td)
Apply rotation matrix to origin.
Definition: pointEdgeCollapseI.H:164
update
mesh update()
Foam::pointEdgeCollapse::enterDomain
void enterDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert relative origin to absolute by adding entering point.
Definition: pointEdgeCollapseI.H:177
Foam::transform
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:521
polyMesh.H
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::pointEdgeCollapse
Determines length of string of edges walked to point.
Definition: pointEdgeCollapse.H:61
Foam::magSqr
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Foam::cmptMin
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:302
Foam::operator==
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::FatalError
error FatalError
w2
#define w2
Definition: blockCreate.C:35
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::pointEdgeCollapse::pointEdgeCollapse
pointEdgeCollapse()
Default construct.
Definition: pointEdgeCollapseI.H:119
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::pointEdgeCollapse::valid
bool valid(TrackingData &td) const
Changed or contains original (invalid) value.
Definition: pointEdgeCollapseI.H:143
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::pointEdgeCollapse::leaveDomain
void leaveDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert origin to relative vector to leaving point.
Definition: pointEdgeCollapseI.H:151
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::Vector< scalar >
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
transform.H
3D tensor transformation operations.
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::pointEdgeCollapse::updatePoint
bool updatePoint(const polyMesh &mesh, const label pointi, const label edgeI, const pointEdgeCollapse &edgeInfo, const scalar tol, TrackingData &td)
Influence of edge on point.
Definition: pointEdgeCollapseI.H:192