Skip to content

GitLab

Commit 2104cce0 authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

convergence implemented as test and torus geometry improved

parent be3f0e91
Hide whitespace changes
Inline Side-by-side
doc/grids/torus3.msh 0 → 100644
View file @ 2104cce0
$MeshFormat
2.2 0 8
$EndMeshFormat
$Nodes
128
0 -1.570572708 -1.570572708 -0.975245423
1 -1.526219626 -2.267431503 -0.6799716893
2 -2.0347697 -2.0347697 -0.4793956186
3 -0.8319823739 -2.453488486 -0.806880817
4 -1.565060182 -2.559411383 -4.052379965e-17
5 -0.7990153789 -1.617518811 -0.9806247566
6 -0.8536219861 -0.8536219861 -0.6094868082
7 -3.929774378e-17 -1.650087367 -0.9367823382
8 -2.267431503 -1.526219626 -0.6799716893
9 -2.453488486 -0.8319823739 -0.806880817
10 -2.559411383 -1.565060182 2.774770665e-17
11 -1.617518811 -0.7990153789 -0.9806247566
12 -1.650087367 -5.002661795e-17 -0.9367823382
13 -3.898116401e-17 -1 3.552713679e-17
14 -0.7071067812 -0.7071067812 3.885780586e-17
15 -0.8536219861 -0.8536219861 0.6094868082
16 -1 -3.343004889e-17 4.200632897e-17
17 -1.570572708 -1.570572708 0.975245423
18 -0.7990153789 -1.617518811 0.9806247566
19 -0.8319823739 -2.453488486 0.806880817
20 -5.002669454e-17 -1.650087367 0.9367823382
21 -0.8609991297 -2.873792007 -5.884184246e-17
22 -5.008339428e-17 -3 -2.664535264e-17
23 -2.121320344 -2.121320344 -5.551115123e-19
24 -2.0347697 -2.0347697 0.4793956186
25 -1.526219626 -2.267431503 0.6799716893
26 -1.617518811 -0.7990153789 0.9806247566
27 -2.453488486 -0.8319823739 0.806880817
28 -1.650087367 4.767296027e-17 0.9367823382
29 -2.267431503 -1.526219626 0.6799716893
30 -2.873792007 -0.8609991297 3.871904945e-17
31 -3 -5.026554023e-17 -5.107025914e-17
32 -4.864133429e-17 -2.545423873 -0.8381603658
33 0.8319823739 -2.453488486 -0.806880817
34 0.8609991297 -2.873792007 1.863089574e-17
35 0.8319823739 -2.453488486 0.806880817
36 1.565060182 -2.559411383 -4.496351188e-17
37 -5.859893473e-17 -2.545423873 0.8381603658
38 0.7990153789 -1.617518811 -0.9806247566
39 0.8536219861 -0.8536219861 -0.6094868082
40 1.570572708 -1.570572708 -0.975245423
41 -4.022603835e-17 -1.152144551 -0.5302274389
42 -3.450143156e-17 -1.152144551 0.5302274389
43 0.8536219861 -0.8536219861 0.6094868082
44 1 -3.475602814e-17 3.719247132e-17
45 0.7071067812 -0.7071067812 1.960237528e-17
46 0.7990153789 -1.617518811 0.9806247566
47 1.570572708 -1.570572708 0.975245423
48 1.526219626 -2.267431503 -0.6799716893
49 2.0347697 -2.0347697 -0.4793956186
50 2.267431503 -1.526219626 -0.6799716893
51 2.453488486 -0.8319823739 -0.806880817
52 2.559411383 -1.565060182 4.163287168e-17
53 1.617518811 -0.7990153789 -0.9806247566
54 1.650087367 -5.224744176e-17 -0.9367823382
55 2.873792007 -0.8609991297 4.884984495e-17
56 2.453488486 -0.8319823739 0.806880817
57 3 -4.786294821e-17 6.527764434e-17
58 2.267431503 -1.526219626 0.6799716893
59 2.0347697 -2.0347697 0.4793956186
60 2.121320344 -2.121320344 -4.857225733e-18
61 1.526219626 -2.267431503 0.6799716893
62 1.617518811 -0.7990153789 0.9806247566
63 1.650087367 2.324785001e-17 0.9367823382
64 -2.545423873 -5.641302345e-17 -0.8381603658
65 -2.453488486 0.8319823739 -0.806880817
66 -1.617518811 0.7990153789 -0.9806247566
67 -0.8536219861 0.8536219861 -0.6094868082
68 -1.570572708 1.570572708 -0.975245423
69 -1.152144551 -5.216010299e-17 -0.5302274389
70 -2.873792007 0.8609991297 -7.549513311e-17
71 -2.453488486 0.8319823739 0.806880817
72 -2.559411383 1.565060182 -4.383691933e-17
73 -2.545423873 5.571947967e-17 0.8381603658
74 -1.152144551 2.416846557e-17 0.5302274389
75 -0.8536219861 0.8536219861 0.6094868082
76 4.539578586e-17 1 3.552713679e-17
77 -0.7071067812 0.7071067812 4.662936703e-17
78 -1.617518811 0.7990153789 0.9806247566
79 -1.570572708 1.570572708 0.975245423
80 1.152144551 -4.910766149e-17 -0.5302274389
81 0.8536219861 0.8536219861 -0.6094868082
82 1.617518811 0.7990153789 -0.9806247566
83 2.453488486 0.8319823739 -0.806880817
84 1.570572708 1.570572708 -0.975245423
85 2.545423873 -5.502508096e-17 -0.8381603658
86 0.7071067812 0.7071067812 3.996802889e-17
87 0.8536219861 0.8536219861 0.6094868082
88 1.152144551 -2.024164674e-17 0.5302274389
89 2.545423873 -2.754712886e-17 0.8381603658
90 2.453488486 0.8319823739 0.806880817
91 2.873792007 0.8609991297 5.367241339e-17
92 2.559411383 1.565060182 1.054760109e-17
93 1.617518811 0.7990153789 0.9806247566
94 1.570572708 1.570572708 0.975245423
95 -0.7990153789 1.617518811 -0.9806247566
96 -0.8319823739 2.453488486 -0.806880817
97 -2.611354105e-17 1.650087367 -0.9367823382
98 -2.267431503 1.526219626 -0.6799716893
99 -2.0347697 2.0347697 -0.4793956186
100 -1.526219626 2.267431503 -0.6799716893
101 -1.565060182 2.559411383 1.674082988e-17
102 -0.8609991297 2.873792007 3.212709455e-17
103 -0.8319823739 2.453488486 0.806880817
104 5.427757004e-17 3 -5.329070512e-17
105 -0.7990153789 1.617518811 0.9806247566
106 5.655447571e-17 1.650087367 0.9367823382
107 -1.526219626 2.267431503 0.6799716893
108 -2.0347697 2.0347697 0.4793956186
109 -2.267431503 1.526219626 0.6799716893
110 -2.121320344 2.121320344 5.551115123e-19
111 2.860828768e-17 1.152144551 -0.5302274389
112 -1.410185133e-17 1.152144551 0.5302274389
113 0.7990153789 1.617518811 -0.9806247566
114 0.8319823739 2.453488486 -0.806880817
115 -4.30904661e-17 2.545423873 -0.8381603658
116 6.172168468e-17 2.545423873 0.8381603658
117 0.8319823739 2.453488486 0.806880817
118 0.8609991297 2.873792007 -4.996022519e-17
119 1.565060182 2.559411383 -4.441015529e-17
120 0.7990153789 1.617518811 0.9806247566
121 2.267431503 1.526219626 -0.6799716893
122 2.0347697 2.0347697 -0.4793956186
123 1.526219626 2.267431503 -0.6799716893
124 2.121320344 2.121320344 9.436895709e-18
125 2.0347697 2.0347697 0.4793956186
126 2.267431503 1.526219626 0.6799716893
127 1.526219626 2.267431503 0.6799716893
$EndNodes
$Elements
256
1 2 0 0 1 2
2 2 0 0 3 1
3 2 0 4 2 1
4 2 0 4 1 3
5 2 0 0 6 5
6 2 0 0 5 3
7 2 0 7 5 6
8 2 0 7 3 5
9 2 0 0 2 8
10 2 0 0 8 9
11 2 0 10 8 2
12 2 0 10 9 8
13 2 0 0 11 6
14 2 0 0 9 11
15 2 0 12 6 11
16 2 0 12 11 9
17 2 0 13 14 15
18 2 0 13 6 14
19 2 0 16 15 14
20 2 0 16 14 6
21 2 0 17 18 15
22 2 0 17 19 18
23 2 0 20 15 18
24 2 0 20 18 19
25 2 0 4 21 19
26 2 0 4 3 21
27 2 0 22 19 21
28 2 0 22 21 3
29 2 0 4 23 2
30 2 0 4 24 23
31 2 0 10 2 23
32 2 0 10 23 24
33 2 0 4 19 25
34 2 0 4 25 24
35 2 0 17 24 25
36 2 0 17 25 19
37 2 0 17 15 26
38 2 0 17 26 27
39 2 0 28 26 15
40 2 0 28 27 26
41 2 0 10 29 27
42 2 0 10 24 29
43 2 0 17 29 24
44 2 0 17 27 29
45 2 0 10 27 30
46 2 0 10 30 9
47 2 0 31 30 27
48 2 0 31 9 30
49 2 0 22 32 33
50 2 0 22 3 32
51 2 0 7 33 32
52 2 0 7 32 3
53 2 0 22 33 34
54 2 0 22 34 35
55 2 0 36 34 33
56 2 0 36 35 34
57 2 0 22 37 19
58 2 0 22 35 37
59 2 0 20 37 35
60 2 0 20 19 37
61 2 0 7 38 33
62 2 0 7 39 38
63 2 0 40 33 38
64 2 0 40 38 39
65 2 0 7 6 41
66 2 0 7 41 39
67 2 0 13 39 41
68 2 0 13 41 6
69 2 0 13 15 42
70 2 0 13 42 43
71 2 0 20 42 15
72 2 0 20 43 42
73 2 0 44 45 43
74 2 0 44 39 45
75 2 0 13 45 39
76 2 0 13 43 45
77 2 0 20 35 46
78 2 0 20 46 43
79 2 0 47 46 35
80 2 0 47 43 46
81 2 0 40 48 33
82 2 0 40 49 48
83 2 0 36 33 48
84 2 0 36 48 49
85 2 0 40 51 50
86 2 0 40 50 49
87 2 0 52 49 50
88 2 0 52 50 51
89 2 0 40 53 51
90 2 0 40 39 53
91 2 0 54 51 53
92 2 0 54 53 39
93 2 0 52 55 56
94 2 0 52 51 55
95 2 0 57 56 55
96 2 0 57 55 51
97 2 0 52 56 58
98 2 0 52 58 59
99 2 0 47 58 56
100 2 0 47 59 58
101 2 0 36 60 59
102 2 0 36 49 60
103 2 0 52 60 49
104 2 0 52 59 60
105 2 0 36 59 61
106 2 0 36 61 35
107 2 0 47 35 61
108 2 0 47 61 59
109 2 0 47 56 62
110 2 0 47 62 43
111 2 0 63 62 56
112 2 0 63 43 62
113 2 0 12 64 65
114 2 0 12 9 64
115 2 0 31 65 64
116 2 0 31 64 9
117 2 0 12 65 66
118 2 0 12 66 67
119 2 0 68 66 65
120 2 0 68 67 66
121 2 0 12 69 6
122 2 0 12 67 69
123 2 0 16 69 67
124 2 0 16 6 69
125 2 0 31 70 65
126 2 0 31 71 70
127 2 0 72 65 70
128 2 0 72 70 71
129 2 0 31 27 73
130 2 0 31 73 71
131 2 0 28 71 73
132 2 0 28 73 27
133 2 0 16 74 15
134 2 0 16 75 74
135 2 0 28 15 74
136 2 0 28 74 75
137 2 0 76 75 77
138 2 0 76 77 67
139 2 0 16 67 77
140 2 0 16 77 75
141 2 0 28 78 71
142 2 0 28 75 78
143 2 0 79 71 78
144 2 0 79 78 75
145 2 0 54 80 81
146 2 0 54 39 80
147 2 0 44 81 80
148 2 0 44 80 39
149 2 0 54 81 82
150 2 0 54 82 83
151 2 0 84 82 81
152 2 0 84 83 82
153 2 0 54 85 51
154 2 0 54 83 85
155 2 0 57 85 83
156 2 0 57 51 85
157 2 0 44 86 81
158 2 0 44 87 86
159 2 0 76 81 86
160 2 0 76 86 87
161 2 0 44 43 88
162 2 0 44 88 87
163 2 0 63 87 88
164 2 0 63 88 43
165 2 0 57 89 56
166 2 0 57 90 89
167 2 0 63 56 89
168 2 0 63 89 90
169 2 0 57 83 91
170 2 0 57 91 90
171 2 0 92 91 83
172 2 0 92 90 91
173 2 0 63 93 87
174 2 0 63 90 93
175 2 0 94 87 93
176 2 0 94 93 90
177 2 0 68 96 95
178 2 0 68 95 67
179 2 0 97 95 96
180 2 0 97 67 95
181 2 0 68 65 98
182 2 0 68 98 99
183 2 0 72 98 65
184 2 0 72 99 98
185 2 0 68 100 96
186 2 0 68 99 100
187 2 0 101 100 99
188 2 0 101 96 100
189 2 0 101 103 102
190 2 0 101 102 96
191 2 0 104 102 103
192 2 0 104 96 102
193 2 0 79 105 103
194 2 0 79 75 105
195 2 0 106 103 105
196 2 0 106 105 75
197 2 0 101 107 103
198 2 0 101 108 107
199 2 0 79 103 107
200 2 0 79 107 108
201 2 0 72 109 108
202 2 0 72 71 109
203 2 0 79 109 71
204 2 0 79 108 109
205 2 0 72 108 110
206 2 0 72 110 99
207 2 0 101 99 110
208 2 0 101 110 108
209 2 0 76 111 81
210 2 0 76 67 111
211 2 0 97 81 111
212 2 0 97 111 67
213 2 0 76 112 75
214 2 0 76 87 112
215 2 0 106 112 87
216 2 0 106 75 112
217 2 0 97 113 81
218 2 0 97 114 113
219 2 0 84 81 113
220 2 0 84 113 114
221 2 0 97 96 115
222 2 0 97 115 114
223 2 0 104 114 115
224 2 0 104 115 96
225 2 0 104 103 116
226 2 0 104 116 117
227 2 0 106 116 103
228 2 0 106 117 116
229 2 0 104 118 114
230 2 0 104 117 118
231 2 0 119 114 118
232 2 0 119 118 117
233 2 0 106 87 120
234 2 0 106 120 117
235 2 0 94 120 87
236 2 0 94 117 120
237 2 0 84 122 121
238 2 0 84 121 83
239 2 0 92 83 121
240 2 0 92 121 122
241 2 0 84 123 122
242 2 0 84 114 123
243 2 0 119 123 114
244 2 0 119 122 123
245 2 0 92 124 125
246 2 0 92 122 124
247 2 0 119 125 124
248 2 0 119 124 122
249 2 0 92 125 126
250 2 0 92 126 90
251 2 0 94 126 125
252 2 0 94 90 126
253 2 0 119 127 125
254 2 0 119 117 127
255 2 0 94 125 127
256 2 0 94 127 117
$EndElements
dune/curvedsurfacegrid/geometries/implicitsurface.hh
View file @ 2104cce0
......@@ -32,7 +32,7 @@ public:
/// \brief Evaluation of the closest-point projection
/**
* Simple iterative algorithm proposed by Demlow and Dziuk in
* Iterative algorithm proposed by Demlow and Dziuk in
* AN ADAPTIVE FINITE ELEMENT METHOD FOR THE LAPLACE-BELTRAMI OPERATOR ON IMPLICITLY DEFINED SURFACES
**/
template <class Domain>
......@@ -83,6 +83,69 @@ private:
DFunctor dphi_;
};
/// \brief Closest-point projection to surface given by implicit function using a simple projection algorithm
/**
* Surface S is given by zero-levelset of function F.
* We assume that F is differentiable in order to evaluate normals and to
* do an iterative projection.
**/
template <class F>
class SimpleImplicitSurfaceProjection
{
using Functor = F;
using DFunctor = std::decay_t<decltype(derivative(std::declval<F>()))>;
public:
/// \brief Constructor for a given differentiable function with surface = { x : phi(x) == 0 }
SimpleImplicitSurfaceProjection (const F& phi, int maxIter = 10)
: maxIter_(maxIter)
, phi_(phi)
, dphi_(derivative(phi_))
{}
/// \brief Evaluation of the closest-point projection
/**
* Simple iterative algorithm proposed by Ingo Nitschke in
* Diskretes Äußeres Kalkül (DEC) auf Oberflächen ohne Rand
**/
template <class Domain>
Domain operator() (Domain x) const
{
using std::sqrt;
using T = typename Domain::value_type;
T tol = sqrt(std::numeric_limits<T>::epsilon());
for (int i = 0; i < maxIter_; ++i) {
auto phi_x = phi_(x);
auto grad_phi_x = dphi_(x);
auto normalize = phi_x / grad_phi_x.two_norm2();
if (sqrt(phi_x * normalize) < tol)
break;
x -= grad_phi_x * normalize;
}
return x;
}
/// \brief Normal vector given by grad(F)/|grad(F)|
template <class Domain>
Domain normal (const Domain& x0) const
{
auto grad_phi = dphi_(x0);
return grad_phi / grad_phi.two_norm();
}
private:
int maxIter_;
Functor phi_;
DFunctor dphi_;
};
/// \brief Construct a grid function representing a projection to an implicit surface
template <class Grid, class F>
auto implicitSurfaceGridFunction (const F& phi, int maxIter = 10)
......@@ -90,6 +153,14 @@ auto implicitSurfaceGridFunction (const F& phi, int maxIter = 10)
return analyticGridFunction<Grid>(ImplicitSurfaceProjection<F>{phi, maxIter});
}
/// \brief Construct a grid function representing a projection to an implicit surface by
/// additionally giving the projection implementation as template
template <class Grid, template <class> Projection, class F>
auto implicitSurfaceGridFunction (const F& phi, int maxIter = 10)
{
return analyticGridFunction<Grid>(Projection<F>{phi, maxIter});
}
} // end namespace Dune
#endif // DUNE_CURVEDSURFACEGRID_IMPLICIT_SURFACE_PROJECTION_HH
dune/curvedsurfacegrid/geometries/torus.hh
View file @ 2104cce0
......@@ -9,8 +9,6 @@
#include <dune/common/math.hh>
#include <dune/curvedsurfacegrid/gridfunctions/analyticgridfunction.hh>
#include "implicitsurface.hh"
namespace Dune {
// torus functor
......@@ -54,30 +52,36 @@ public:
TorusProjection (T R, T r)
: R_(R)
, r_(r)
, implicit_(Phi{R,r},100)
{}
/// \brief Closest point projection
Domain operator() (const Domain& x) const
{
return implicit_(x);
using std::sqrt;
auto scale1 = R_ / sqrt(x[0]*x[0] + x[1]*x[1]);
Domain center{x[0] * scale1, x[1] * scale1, 0};
Domain out{x[0] - center[0], x[1] - center[1], x[2]};
out *= r_ / out.two_norm();
return out + center;
}
/// \brief Normal vector
Domain normal (const Domain& x) const
{
return implicit_.normal(x);
using std::sqrt;
auto X = (*this)(x);
T factor = (1 - R_/sqrt(X[0]*X[0] + X[1]*X[1]))/r_;
return { X[0]*factor, X[1]*factor, X[2]/r_ };
}
/// \brief Mean curvature
T mean_curvature (const FieldVector<T,3>& x) const
{
using std::sqrt;
auto X = implicit_(x);
T nrm = sqrt(X[0]*X[0] + X[1]*X[1]);
T zr = X[2]/r_;
T cos_v = (nrm < R_ ? -1 : 1) * sqrt(1 - zr*zr);
return (R_ + 2*r_*cos_v)/(2*r_*(R_ + r_*cos_v));
auto X = (*this)(x);
return (2 - R_/sqrt(X[0]*X[0] + X[1]*X[1]))/(2*r_);
}
/// \brief surface area of the torus = 4*pi^2*R*r
......@@ -88,7 +92,6 @@ public:
private:
T R_, r_;
ImplicitSurfaceProjection<Phi> implicit_;
};
/// \brief construct a grid function representing a torus parametrization
......
dune/curvedsurfacegrid/test/convergence.cc
View file @ 2104cce0
......@@ -30,7 +30,7 @@ const int quad_order = order+6;
#if GEOMETRY_TYPE < 3
const int num_levels = 4;
#else
const int num_levels = 3;
const int num_levels = 4;
#endif
// Test Hausdorf-distance
......@@ -186,7 +186,7 @@ int main (int argc, char** argv)
std::unique_ptr<HostGrid> hostGrid = GmshReader<HostGrid>::read( DUNE_GRID_PATH "ellipsoid3.msh");
auto gridFct = ellipsoidGridFunction<HostGrid>(1.0, 0.75, 1.25);
#elif GEOMETRY_TYPE == 3
std::unique_ptr<HostGrid> hostGrid = GmshReader<HostGrid>::read( DUNE_GRID_PATH "torus.msh");
std::unique_ptr<HostGrid> hostGrid = GmshReader<HostGrid>::read( DUNE_GRID_PATH "torus2.msh");
auto gridFct = torusGridFunction<HostGrid>(2.0, 1.0);
#endif
......@@ -197,9 +197,11 @@ int main (int argc, char** argv)
if (i > 0)
grid.globalRefine(1);
#if WRITE_FILES
using DC = LagrangeDataCollector<typename Grid::LeafGridView, order>;
VtkUnstructuredGridWriter<typename Grid::LeafGridView, DC> writer(grid.leafGridView());
writer.write("level_" + std::to_string(i) + ".vtu");
#endif
std::cout << "level = " << i << std::endl;
inf_errors.push_back( inf_error(grid) );
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment