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_sphere_projector_new
34 : #define __H__UG_sphere_projector_new
35 :
36 : #include "refinement_projector.h"
37 :
38 : namespace ug{
39 :
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 SphereProjector : public RefinementProjector {
50 : public:
51 0 : SphereProjector () :
52 : m_center (0, 0, 0),
53 0 : m_radius (-1),
54 0 : m_influenceRadius (-1)
55 : {}
56 :
57 0 : SphereProjector (const vector3& center) :
58 : m_center (center),
59 0 : m_radius (-1),
60 0 : m_influenceRadius (-1)
61 : {}
62 :
63 : SphereProjector (const vector3& center,
64 0 : number radius) :
65 : m_center (center),
66 0 : m_radius (radius),
67 0 : m_influenceRadius (-1)
68 : {}
69 :
70 : SphereProjector (const vector3& center,
71 : number radius,
72 0 : number influenceRadius) :
73 : m_center (center),
74 0 : m_radius (radius),
75 0 : m_influenceRadius (influenceRadius)
76 : {}
77 :
78 : /** \sa ug::RefinementProjector::RefinementProjector*/
79 0 : SphereProjector (SPIGeometry3d geometry,
80 0 : const vector3& center) :
81 : RefinementProjector (geometry),
82 : m_center (center),
83 0 : m_radius (-1),
84 0 : m_influenceRadius (-1)
85 0 : {}
86 :
87 : /** \sa ug::RefinementProjector::RefinementProjector*/
88 : SphereProjector (SPIGeometry3d geometry,
89 : const vector3& center,
90 : number radius) :
91 : RefinementProjector (geometry),
92 : m_center (center),
93 : m_radius (radius),
94 : m_influenceRadius (-1)
95 : {}
96 :
97 : /** \sa ug::RefinementProjector::RefinementProjector*/
98 : template <class TGeomProvider>
99 0 : SphereProjector (const TGeomProvider& geometry,
100 : const vector3& center,
101 : number radius,
102 : number influenceRadius) :
103 : RefinementProjector (geometry),
104 : m_center (center),
105 0 : m_radius (radius),
106 0 : m_influenceRadius (influenceRadius)
107 0 : {}
108 :
109 0 : virtual ~SphereProjector () {}
110 :
111 0 : void set_center (const vector3& center) {m_center = center;}
112 0 : const vector3& center () const {return m_center;}
113 :
114 0 : void set_radius (number radius) {m_radius = radius;}
115 0 : number radius () const {return m_radius;}
116 :
117 0 : void set_influence_radius (number influenceRadius) {m_influenceRadius = influenceRadius;}
118 0 : number influence_radius () const {return m_influenceRadius;}
119 :
120 : /// called when a new vertex was created from an old edge.
121 0 : virtual number new_vertex(Vertex* vrt, Edge* parent)
122 : {
123 0 : return perform_projection(vrt, parent);
124 : }
125 :
126 : /// called when a new vertex was created from an old face.
127 0 : virtual number new_vertex(Vertex* vrt, Face* parent)
128 : {
129 0 : return perform_projection(vrt, parent);
130 : }
131 :
132 : /// called when a new vertex was created from an old volume.
133 0 : virtual number new_vertex(Vertex* vrt, Volume* parent)
134 : {
135 0 : return perform_projection(vrt, parent);
136 : }
137 :
138 : private:
139 :
140 : template <class TElem>
141 0 : number perform_projection(Vertex* vrt, TElem* parent)
142 : {
143 : // first calculate the average distance of corners of the parent to the
144 : // parents center
145 0 : typename TElem::ConstVertexArray vrts = parent->vertices();
146 0 : size_t numVrts = parent->num_vertices();
147 :
148 0 : if(numVrts == 0){
149 0 : set_pos(vrt, vector3(0, 0, 0));
150 0 : return 1;
151 : }
152 :
153 : number avDist = 0;
154 : vector3 parentCenter(0, 0, 0);
155 :
156 0 : for(size_t i = 0; i < numVrts; ++i){
157 0 : vector3 p = pos(vrts[i]);
158 0 : avDist += VecDistance(p, m_center);
159 : parentCenter += p;
160 : }
161 :
162 0 : avDist /= (number)numVrts;
163 0 : VecScale(parentCenter, parentCenter, 1. / (number)numVrts);
164 :
165 : // calculate projection
166 : vector3 proj;
167 : VecSubtract(proj, parentCenter, m_center);
168 : number len = VecLength(proj);
169 0 : if(len > SMALL * avDist){ // if avDist is very small, len may be small, too
170 0 : VecScale(proj, proj, avDist / len);
171 : proj += m_center;
172 : set_pos(vrt, proj);
173 : }
174 : else
175 : set_pos(vrt, parentCenter);
176 :
177 0 : if(m_influenceRadius > 0) {
178 0 : if(m_radius > m_influenceRadius){
179 0 : const number dist = m_radius - m_influenceRadius;
180 0 : if(dist > 0)
181 0 : return clip<number>((len - m_influenceRadius) / dist, 0, 1);
182 0 : return len > m_radius ? 1 : 0;
183 : }
184 0 : else if(m_radius >= 0){
185 0 : const number dist = m_influenceRadius - m_radius;
186 0 : if(dist > 0)
187 0 : return clip<number>(1 - (len - m_radius) / dist, 0, 1);
188 0 : return len < m_radius ? 1 : 0;
189 : }
190 : else
191 0 : return clip<number>(1 - len / m_influenceRadius, 0, 1);
192 : }
193 :
194 : return 1;
195 : }
196 :
197 :
198 : friend class boost::serialization::access;
199 :
200 : template <class Archive>
201 0 : void serialize( Archive& ar, const unsigned int version)
202 : {
203 0 : ar & make_nvp("center", m_center);
204 0 : ar & make_nvp("radius", m_radius);
205 0 : ar & make_nvp("influence radius", m_influenceRadius);
206 : UG_EMPTY_BASE_CLASS_SERIALIZATION(SphereProjector, RefinementProjector);
207 0 : }
208 :
209 : vector3 m_center;
210 : number m_radius;
211 : number m_influenceRadius;
212 : };
213 :
214 :
215 : }// end of namespace
216 :
217 : #endif //__H__UG_sphere_projector_new
|
