CV2DI.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) 2013-2016 OpenFOAM Foundation
9 Copyright (C) 2019 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30
31inline Foam::label Foam::CV2D::insertPoint
32(
33 const point2D& p,
34 const label type
35)
36{
37 unsigned int nVert = number_of_vertices();
38
39 return insertPoint(toPoint(p), nVert, type);
40}
41
42
43inline Foam::label Foam::CV2D::insertPoint
44(
45 const point2D& p,
46 const label index,
47 const label type
48)
49{
50 return insertPoint(toPoint(p), index, type);
51}
52
53
54inline Foam::label Foam::CV2D::insertPoint
55(
56 const Point& p,
57 const label index,
58 const label type
59)
60{
61 unsigned int nVert = number_of_vertices();
62
63 Vertex_handle vh = insert(p);
64
65 if (nVert == number_of_vertices())
66 {
68 << "Failed to insert point " << toPoint2D(p) << endl;
69 }
70 else
71 {
72 vh->index() = index;
73 vh->type() = type;
74 }
75
76 return vh->index();
77}
78
79
80inline bool Foam::CV2D::insertMirrorPoint
81(
82 const point2D& nearSurfPt,
83 const point2D& surfPt
84)
85{
86 point2D mirrorPoint(2*surfPt - nearSurfPt);
87
88 if (qSurf_.outside(toPoint3D(mirrorPoint)))
89 {
90 insertPoint(mirrorPoint, Vb::MIRROR_POINT);
91 return true;
92 }
93
94 return false;
95}
96
97
98inline void Foam::CV2D::insertPointPair
99(
100 const scalar ppDist,
101 const point2D& surfPt,
102 const vector2D& n
103)
104{
105 vector2D ppDistn = ppDist*n;
106
107 label master = insertPoint
108 (
109 surfPt - ppDistn,
110 number_of_vertices() + 1
111 );
112
113 insertPoint(surfPt + ppDistn, master);
114}
115
116
117// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118
120{
121 return controls_;
122}
123
124
126{
127 return reinterpret_cast<const point2D&>(p);
128}
129
130
132{
133 point2DField temp(p.size());
134 forAll(temp, pointi)
135 {
136 temp[pointi] = point2D(p[pointi].x(), p[pointi].y());
137 }
138 return temp;
139}
140
141
143{
144 return Foam::point(p.x(), p.y(), z_);
145}
146
147
148#ifdef CGAL_INEXACT
149
151{
152 return reinterpret_cast<point2DFromPoint>(P);
153}
154
155
157{
158 return reinterpret_cast<PointFromPoint2D>(p);
159}
160
161#else
162
164{
165 return point2D(CGAL::to_double(P.x()), CGAL::to_double(P.y()));
166}
167
168
170{
171 return Point(p.x(), p.y());
172}
173
174#endif
175
176
178{
179 return Foam::point(CGAL::to_double(P.x()), CGAL::to_double(P.y()), z_);
180}
181
182
183inline void Foam::CV2D::movePoint(const Vertex_handle& vh, const Point& P)
184{
185 int i = vh->index();
186 int t = vh->type();
187
188 remove(vh);
189
190 Vertex_handle newVh = insert(P);
191
192 newVh->index() = i;
193 newVh->type() = t;
194
195 // label i = vh->index();
196 // move(vh, P);
197 // vh->index() = i;
198
199 //vh->set_point(P);
200 //fast_restore_Delaunay(vh);
201}
202
203
204// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
205
207{
208 return boundaryTriangle
209 (
210 *fc->vertex(0),
211 *fc->vertex(1),
212 *fc->vertex(2)
213 );
214}
215
216
218{
219 return outsideTriangle
220 (
221 *fc->vertex(0),
222 *fc->vertex(1),
223 *fc->vertex(2)
224 );
225}
226
227
228// ************************************************************************* //
CGAL::Point_3< K > Point
scalar y
label n
Fb::Face_handle Face_handle
Definition: indexedFace.H:71
const Point & PointFromPoint2D
Definition: CV2D.H:357
void movePoint(const Vertex_handle &vh, const Point &P)
Definition: CV2DI.H:183
const cv2DControls & meshControls() const
Definition: CV2DI.H:119
PointFromPoint2D toPoint(const point2D &) const
Definition: CV2DI.H:169
Foam::point toPoint3D(const point2D &) const
Definition: CV2DI.H:142
const point2D & toPoint2D(const Foam::point &) const
Definition: CV2DI.H:125
Generic templated field type.
Definition: Field.H:82
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Controls for the 2D CV mesh generator.
Definition: cv2DControls.H:61
volScalarField & p
#define WarningInFunction
Report a warning using Foam::Warning.
Vector2D< scalar > vector2D
A 2D vector of scalars obtained from the generic Vector2D.
Definition: vector2D.H:51
PointFrompoint toPoint(const Foam::point &p)
bool boundaryTriangle(const CV2D::Face_handle fc)
Definition: CV2DI.H:206
vector point
Point is a vector.
Definition: point.H:43
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
vector2D point2D
Point2D is a vector.
Definition: point2D.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
bool outsideTriangle(const CV2D::Face_handle fc)
Definition: CV2DI.H:217
srcOptions insert("case", fileName(rootDirSource/caseDirSource))
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333