Line data Source code
1 : /*
2 : * Copyright (c) 2010-2015: 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__SUBDIVISION_RULES_PIECEWISE_LOOP__
34 : #define __H__UG__SUBDIVISION_RULES_PIECEWISE_LOOP__
35 :
36 : #include <vector>
37 : #include <cassert>
38 : #include "lib_grid/lg_base.h"
39 :
40 : namespace ug
41 : {
42 :
43 : /// \addtogroup lib_grid_algorithms_refinement_subdivision
44 : /// @{
45 :
46 : /// A singleton that stores all rules for a piecewise-loop subdivision surface.
47 : class SubdivRules_PLoop
48 : {
49 : public:
50 : struct NeighborInfo{
51 0 : NeighborInfo() {}
52 : NeighborInfo(Vertex* n, size_t cval) :
53 : nbr(n), creaseValence(cval) {}
54 :
55 : Vertex* nbr;
56 : /** 0 means that the neighbor is not a crease vertex.
57 : * > 0: The valence of the crease regarding only the
58 : * part on the side of the center-vertex.*/
59 : size_t creaseValence;
60 : };
61 :
62 : public:
63 : /// returns the only instance to this singleton.
64 0 : static SubdivRules_PLoop& inst()
65 : {
66 0 : static SubdivRules_PLoop subdivRules;
67 0 : return subdivRules;
68 : }
69 :
70 : ////////////////////////////////
71 : // WEIGHTS
72 : number ref_even_inner_center_weight(size_t valence) const
73 0 : {return 1. - (number)valence * get_beta(valence);}
74 :
75 : number ref_even_inner_nbr_weight(size_t valence) const
76 0 : {return get_beta(valence);}
77 :
78 : /// returns weights for center, nbr1 and nbr2.
79 : vector3 ref_even_crease_weights() const
80 : {return vector3(0.75, 0.125, 0.125);}
81 :
82 : vector4 ref_odd_inner_weights() const
83 0 : {return vector4(0.375, 0.375, 0.125, 0.125);}
84 :
85 : /// weights of an odd vertex on an inner edge that is connected to a crease.
86 : /** The weight for the vertex on the crease is in v.x(), the inner edge vertex
87 : * in v.y() and the two indirectly connected vertex weights are in v.z() and v.w.
88 : * creaseValence specifies the number of associated edges of the crease vertex.
89 : *
90 : * Rules are taken from:
91 : * "Piecewise Smooth Subdivision Surfaces with Normal Control",
92 : * H. Biermann, Adi Levin, Denis Zorin.*/
93 0 : vector4 ref_odd_inner_weights(size_t creaseValence) const
94 : {
95 : assert(creaseValence > 2 && "Bad crease valence. Underlying grid is not a surface grid.");
96 0 : if(creaseValence == 4)
97 : return ref_odd_inner_weights();
98 0 : number gamma = 0.5 - 0.25 * cos(PI / (number)(creaseValence - 1));
99 0 : return vector4(0.75 - gamma, gamma, 0.125, 0.125);
100 :
101 : //number c = 0.25 * cos((2.*PI) / (number)(creaseValence - 1));
102 : //return vector4(0.5 - c, 0.25 + c, 0.125, 0.125);
103 : }
104 :
105 : vector2 ref_odd_crease_weights() const
106 : {return vector2(0.5, 0.5);}
107 :
108 : number proj_inner_center_weight(size_t valence) const
109 0 : {return 1.0 - (number)valence / (0.375 / get_beta(valence) + valence);}
110 :
111 : number proj_inner_nbr_weight(size_t valence) const
112 0 : {return 1.0 / (0.375 / get_beta(valence) + valence);}
113 :
114 : vector3 proj_crease_weights() const
115 : {return vector3(2./3., 1./6., 1./6.);}
116 :
117 : /** nbrInfos have to be specified in clockwise or counter-clockwise order.*/
118 0 : void proj_inner_crease_nbr_weights(number& centerWgtOut, number* nbrWgtsOut,
119 : NeighborInfo* nbrInfos, size_t numNbrs) const
120 : {
121 : number wcntrProj = proj_inner_center_weight(numNbrs);
122 : number wnbrProj = proj_inner_nbr_weight(numNbrs);
123 :
124 : // initialize all weights with 0
125 0 : centerWgtOut = 0;
126 0 : for(size_t i = 0; i < numNbrs; ++i)
127 0 : nbrWgtsOut[i] = 0;
128 :
129 0 : for(size_t i = 0; i < numNbrs; ++i){
130 0 : NeighborInfo& nbrInfo = nbrInfos[i];
131 : vector4 oddWeights(0.5, 0.5, 0, 0);
132 :
133 0 : if(nbrInfo.creaseValence == 0){
134 : // Calculate the weights for odd-refinement on the given edge
135 : oddWeights = ref_odd_inner_weights();
136 : }
137 : else{
138 0 : oddWeights = ref_odd_inner_weights(nbrInfo.creaseValence);
139 : }
140 :
141 0 : nbrWgtsOut[i] += oddWeights.x() * wnbrProj;
142 0 : centerWgtOut += oddWeights.y() * wnbrProj;
143 0 : nbrWgtsOut[next_ind(i, numNbrs)] += oddWeights.z() * wnbrProj;
144 0 : nbrWgtsOut[prev_ind(i, numNbrs)] += oddWeights.w() * wnbrProj;
145 : }
146 :
147 :
148 : // add scaled weights for evem refinement mask
149 0 : centerWgtOut += wcntrProj * ref_even_inner_center_weight(numNbrs);
150 0 : for(size_t i = 0; i < numNbrs; ++i)
151 0 : nbrWgtsOut[i] += wcntrProj * ref_even_inner_nbr_weight(numNbrs);
152 0 : }
153 : /*
154 : number proj_inner_crease_nbr_center_weight(size_t valence);
155 : number proj_inner_crease_nbr_nbr_weight(size_t valence);
156 : number proj_inner_crease_nbr_nbr_weight(size_t valence, size_t cValence);
157 : */
158 : /// returns beta as it is used in the subdivision masks.
159 : /** performs a lookup if the valency is small enough.
160 : * calculates a fresh beta else.*/
161 : number get_beta(size_t valency) const;
162 : /*
163 : ////////////////////////////////
164 : // EVEN MASKS
165 : template <class TAAPos>
166 : typename TAAPos::ValueType
167 : apply_even_mask(Grid& grid, Vertex* center,
168 : TAAPos& aaPos);
169 :
170 : template <class TAAPos>
171 : typename TAAPos::ValueType
172 : apply_even_crease_mask(Vertex* center, Vertex* nbr1,
173 : Vertex* nbr2, TAAPos& aaPos);
174 :
175 : ////////////////////////////////
176 : // ODD MASKS
177 : template <class TAAPos>
178 : typename TAAPos::ValueType
179 : apply_odd_mask(Vertex* vrt, Edge* parent,
180 : TAAPos& aaPos);
181 :
182 : template <class TAAPos>
183 : typename TAAPos::ValueType
184 : apply_odd_crease_mask(Vertex* vrt, Edge* parent,
185 : TAAPos& aaPos);
186 :
187 : template <class TAAPos>
188 : typename TAAPos::ValueType
189 : apply_odd_crease_nbr_mask(Grid& grid, Vertex* vrt,
190 : Edge* parent, TAAPos& aaPos);
191 :
192 : ////////////////////////////////
193 : // PROJECTION
194 : template <class TAAPos>
195 : typename TAAPos::ValueType
196 : project_inner_vertex(Grid& grid, Vertex* vrt,
197 : TAAPos& aaPos);
198 :
199 : template <class TAAPos>
200 : typename TAAPos::ValueType
201 : project_inner_vertex(Vertex* vrt, Vertex* nbrs,
202 : int* nbrCreaseValencies, int numNbrs,
203 : TAAPos& aaPos);
204 :
205 : template <class TAAPos>
206 : typename TAAPos::ValueType
207 : project_crease_vertex(Vertex* vrt, Vertex* nbr1,
208 : Vertex* nbr2, TAAPos& aaPos);
209 : */
210 : private:
211 : /// calculates beta as it is used in the subdivision masks.
212 : number calculate_beta(size_t valency) const;
213 :
214 : private:
215 : /// private constructor prohibits multiple instantiation.
216 : SubdivRules_PLoop();
217 :
218 : /// private copy constructor prohibits copy-construction.
219 : SubdivRules_PLoop(const SubdivRules_PLoop& src);
220 :
221 : /// private assignment operator prohibits assignment.
222 : SubdivRules_PLoop& operator=(const SubdivRules_PLoop& src);
223 :
224 : /// returns the next index in a cyclic index set
225 0 : inline size_t next_ind(size_t ind, size_t numInds) const {return (ind + 1) % numInds;}
226 : /// returns the previous index in a cyclic index set
227 0 : inline size_t prev_ind(size_t ind, size_t numInds) const {return (ind + numInds - 1) % numInds;}
228 :
229 :
230 : private:
231 : std::vector<number> m_betas;//< precalculated betas.
232 : };
233 :
234 : /// @} // end of add_to_group
235 :
236 : }// end of namespace
237 :
238 : #endif
|