Line data Source code
1 : /*
2 : * Copyright (c) 2016: G-CSC, Goethe University Frankfurt
3 : * Author: Sebastian Reiter
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 : #ifndef __H__UG_cylinder_projector_new
34 : #define __H__UG_cylinder_projector_new
35 :
36 : #include "common/math/misc/math_util.h"
37 : #include "refinement_projector.h"
38 :
39 : namespace ug{
40 : /// Projects new vertices onto a sphere during refinement
41 : /** For projection during refinement the radius property is ignored. Instead
42 : * the distance to the center of a newly inserted vertex is calculated
43 : * as the average distance of the vertices of the parent element to the center.
44 : * The radius property thus defaults to -1.
45 : *
46 : * You may still specify a radius. This radius can be used for auto-fitting of
47 : * the center and for reprojecting a set of vertices onto the sphere.
48 : */
49 : class CylinderProjector : public RefinementProjector {
50 : public:
51 0 : CylinderProjector () :
52 : m_center (0, 0, 0),
53 : m_axis (0, 0, 1),
54 0 : m_radius (-1),
55 0 : m_influenceRadius (-1)
56 : {}
57 :
58 : CylinderProjector (const vector3& center,
59 0 : const vector3& axis) :
60 : m_center (center),
61 : m_axis (axis),
62 0 : m_radius (-1),
63 0 : m_influenceRadius (-1)
64 : {}
65 :
66 : CylinderProjector (const vector3& center,
67 : const vector3& axis,
68 0 : number radius) :
69 : m_center (center),
70 : m_axis (axis),
71 0 : m_radius (radius),
72 0 : m_influenceRadius (-1)
73 : {}
74 :
75 : CylinderProjector (const vector3& center,
76 : const vector3& axis,
77 : number radius,
78 0 : number influenceRadius) :
79 : m_center (center),
80 : m_axis (axis),
81 0 : m_radius (radius),
82 0 : m_influenceRadius (influenceRadius)
83 : {}
84 :
85 : /** \sa ug::RefinementProjector::RefinementProjector*/
86 0 : CylinderProjector (SPIGeometry3d geometry,
87 : const vector3& center,
88 : const vector3& axis,
89 : number radius,
90 0 : number influenceRadius) :
91 : RefinementProjector (geometry),
92 : m_center (center),
93 : m_axis (axis),
94 0 : m_radius (radius),
95 0 : m_influenceRadius (influenceRadius)
96 0 : {}
97 :
98 0 : virtual ~CylinderProjector () {}
99 :
100 0 : void set_center (const vector3& center) {m_center = center;}
101 0 : const vector3& center () const {return m_center;}
102 :
103 0 : void set_axis (const vector3& axis) {m_axis = axis;}
104 0 : const vector3& axis () const {return m_axis;}
105 :
106 0 : void set_radius (number radius) {m_radius = radius;}
107 0 : number radius () const {return m_radius;}
108 :
109 0 : void set_influence_radius (number influenceRadius) {m_influenceRadius = influenceRadius;}
110 0 : number influence_radius () const {return m_influenceRadius;}
111 :
112 : /// called when a new vertex was created from an old edge.
113 0 : virtual number new_vertex(Vertex* vrt, Edge* parent)
114 : {
115 0 : return perform_projection(vrt, parent);
116 : }
117 :
118 : /// called when a new vertex was created from an old face.
119 0 : virtual number new_vertex(Vertex* vrt, Face* parent)
120 : {
121 0 : return perform_projection(vrt, parent);
122 : }
123 :
124 : /// called when a new vertex was created from an old volume.
125 0 : virtual number new_vertex(Vertex* vrt, Volume* parent)
126 : {
127 0 : return perform_projection(vrt, parent);
128 : }
129 :
130 : private:
131 :
132 : template <class TElem>
133 0 : number perform_projection(Vertex* vrt, TElem* parent)
134 : {
135 : // calculate the new position by linear interpolation and project that point
136 : // onto the cylinder.
137 0 : typename TElem::ConstVertexArray vrts = parent->vertices();
138 0 : size_t numVrts = parent->num_vertices();
139 :
140 0 : if(numVrts == 0){
141 0 : set_pos(vrt, vector3(0, 0, 0));
142 0 : return 1;
143 : }
144 :
145 : number avDist = 0;
146 : vector3 parentCenter (0, 0, 0);
147 :
148 0 : for(size_t i = 0; i < numVrts; ++i){
149 0 : vector3 p = pos(vrts[i]);
150 0 : avDist += DistancePointToRay(p, m_center, m_axis);
151 : parentCenter += p;
152 : }
153 :
154 0 : avDist /= (number)numVrts;
155 0 : VecScale(parentCenter, parentCenter, 1. / (number)numVrts);
156 :
157 : vector3 proj, v;
158 0 : ProjectPointToRay(proj, parentCenter, m_center, m_axis);
159 : VecSubtract(v, parentCenter, proj);
160 : number len = VecLength(v);
161 0 : if(len > SMALL * avDist){ // if avDist is very small, len may be small, too
162 0 : VecScale(v, v, avDist/len);
163 : proj += v;
164 : set_pos(vrt, proj);
165 : }
166 : else
167 : set_pos(vrt, parentCenter);
168 :
169 0 : if(m_influenceRadius > 0) {
170 0 : if(m_radius > m_influenceRadius){
171 0 : const number dist = m_radius - m_influenceRadius;
172 0 : if(dist > 0)
173 0 : return clip<number>((len - m_influenceRadius) / dist, 0, 1);
174 0 : return len > m_radius ? 1 : 0;
175 : }
176 0 : else if(m_radius >= 0){
177 0 : const number dist = m_influenceRadius - m_radius;
178 0 : if(dist > 0)
179 0 : return clip<number>(1 - (len - m_radius) / dist, 0, 1);
180 0 : return len < m_radius ? 1 : 0;
181 : }
182 : else
183 0 : return clip<number>(1 - len / m_influenceRadius, 0, 1);
184 : }
185 : return 1;
186 : }
187 :
188 :
189 : friend class boost::serialization::access;
190 :
191 : template <class Archive>
192 0 : void serialize( Archive& ar, const unsigned int version)
193 : {
194 0 : ar & make_nvp("center", m_center);
195 0 : ar & make_nvp("axis", m_axis);
196 0 : ar & make_nvp("radius", m_radius);
197 0 : ar & make_nvp("influence radius", m_influenceRadius);
198 : UG_EMPTY_BASE_CLASS_SERIALIZATION(CylinderProjector, RefinementProjector);
199 0 : }
200 :
201 : vector3 m_center;
202 : vector3 m_axis;
203 : number m_radius;
204 : number m_influenceRadius;
205 : };
206 :
207 : }// end of namespace
208 :
209 : #endif //__H__UG_cylinder_projector_new
|
