Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
reference_element_aux.icc
Go to the documentation of this file.
1#ifndef _RHEOLEF_REFERENCE_ELEMENT_AUX_ICC
2#define _RHEOLEF_REFERENCE_ELEMENT_AUX_ICC
23// minimal auto-contained rheolef helpers for msh2geo
24// used in:
25// - reference_element.cc
26// - msh2geo.cc
27// in order to avoid code duplication
28//
29namespace rheolef {
30
31// used by msh2geo.cc that do not include reference_element.h
32#ifndef _RHEOLEF_REFERENCE_ELEMENT_H
33constexpr size_t reference_element__p = 0;
34constexpr size_t reference_element__e = 1;
35constexpr size_t reference_element__t = 2;
36constexpr size_t reference_element__q = 3;
37constexpr size_t reference_element__T = 4;
38constexpr size_t reference_element__P = 5;
39constexpr size_t reference_element__H = 6;
40constexpr size_t reference_element__max_variant = 7;
41#else // _RHEOLEF_REFERENCE_ELEMENT_H
42constexpr size_t reference_element__p = reference_element::p;
43constexpr size_t reference_element__e = reference_element::e;
44constexpr size_t reference_element__t = reference_element::t;
45constexpr size_t reference_element__q = reference_element::q;
46constexpr size_t reference_element__T = reference_element::T;
47constexpr size_t reference_element__P = reference_element::P;
48constexpr size_t reference_element__H = reference_element::H;
49constexpr size_t reference_element__max_variant = reference_element::max_variant;
50#endif // _RHEOLEF_REFERENCE_ELEMENT_H
51
52static size_t _reference_element_dimension_by_variant[reference_element__max_variant] = {0, 1, 2, 2, 3, 3, 3};
53static size_t _reference_element_n_edge_by_variant [reference_element__max_variant] = {0, 1, 3, 4, 6, 9, 12};
54static size_t _reference_element_n_face_by_variant [reference_element__max_variant] = {0, 0, 1, 1, 4, 5, 6};
55
56static
57size_t
58reference_element_variant (const char* name_table, char name)
59{
60 for (size_t i = 0; i < reference_element__max_variant; i++) {
61 if (name_table [i] == name) return i;
62 }
63 return reference_element__max_variant;
64}
65static
66size_t
67reference_element_dimension_by_name(char element_name) {
68 switch (element_name) {
69 case 'p': return 0;
70 case 'e': return 1;
71 case 't':
72 case 'q': return 2;
73 case 'T':
74 case 'P':
75 case 'H': return 3;
76 default: error_macro("unexpected element name `"<<element_name<<"' (ascii="<<int(element_name)<<")");
77 return 0;
78 }
79}
80static
81size_t
82reference_element_dimension_by_variant (size_t variant)
83{
84 return _reference_element_dimension_by_variant[variant];
85}
86static
87size_t
88reference_element_variant (char element_name)
89{
90 switch (element_name) {
91 case 'p': return reference_element__p;
92 case 'e': return reference_element__e;
93 case 't': return reference_element__t;
94 case 'q': return reference_element__q;
95 case 'T': return reference_element__T;
96 case 'P': return reference_element__P;
97 case 'H': return reference_element__H;
98 default: error_macro("unexpected element name `"<<element_name<<"'\"=char("<<int(element_name)<<")");
99 return 0;
100 }
101}
102static
103size_t
104reference_element_n_vertex (size_t variant)
105{
106 switch (variant) {
107 case reference_element__p: return 1;
108 case reference_element__e: return 2;
109 case reference_element__t: return 3;
110 case reference_element__q: return 4;
111 case reference_element__T: return 4;
112 case reference_element__P: return 6;
113 case reference_element__H: return 8;
114 default: error_macro ("unexpected element variant `"<<variant<<"'"); return 0;
115 }
116}
117static
118size_t
119reference_element_n_node (size_t variant, size_t order)
120{
121 switch (variant) {
122 case reference_element__p: return 1;
123 case reference_element__e: return order+1;
124 case reference_element__t: return (order+1)*(order+2)/2;
125 case reference_element__q: return (order+1)*(order+1);
126 case reference_element__T: return (order+1)*(order+2)*(order+3)/6;
127 case reference_element__P: return (order+1)*(order+1)*(order+2)/2;
128 case reference_element__H: return (order+1)*(order+1)*(order+1);
129 default: error_macro ("unexpected element variant `"<<variant<<"'"); return 0;
130 }
131}
132static
133void
134reference_element_init_local_nnode_by_variant (
135 size_t order,
136 std::array<size_t,reference_element__max_variant>& sz)
137{
138 check_macro (order > 0, "invalid order="<<order<<": expect order > 0");
139 sz [reference_element__p] = 1;
140 sz [reference_element__e] = order-1;
141 sz [reference_element__t] = (order-1)*(order-2)/2;
142 sz [reference_element__q] = (order-1)*(order-1);
143 sz [reference_element__T] = (order-1)*(order-2)*(order-3)/6;
144 sz [reference_element__P] = (order-1)*(order-1)*(order-2)/2;
145 sz [reference_element__H] = (order-1)*(order-1)*(order-1);
146}
147// -----------------------------------------------------------------------------
148// reference_element::first_inod_by_dimension
149// -----------------------------------------------------------------------------
150static const size_t
152 reference_element__p,
153 reference_element__e,
154 reference_element__t,
155 reference_element__T,
156 reference_element__max_variant
157};
158static
159size_t
160reference_element_first_variant_by_dimension (size_t dim)
161{
162 return _first_variant_by_dimension[dim];
163}
164static
165size_t
166reference_element_last_variant_by_dimension (size_t dim)
167{
169}
170// -----------------------------------------------------------------------------
171// reference_element::first_inod_by_variant
172// -----------------------------------------------------------------------------
173static
174size_t
175reference_element_p_first_inod_by_variant (
176 size_t order,
177 size_t subgeo_variant)
178{
179 if (subgeo_variant == reference_element__p) return 0;
180 return 1;
181}
182static
183size_t
184reference_element_e_first_inod_by_variant (
185 size_t order,
186 size_t subgeo_variant)
187{
188 if (subgeo_variant == reference_element__p) return 0;
189 if (subgeo_variant == reference_element__e) return 2;
190 return order+1;
191}
192static
193size_t
194reference_element_t_first_inod_by_variant (
195 size_t order,
196 size_t subgeo_variant)
197{
198 if (subgeo_variant == reference_element__p) return 0;
199 if (subgeo_variant == reference_element__e) return 3;
200 if (subgeo_variant == reference_element__t) return 3 + 3*(order-1);
201 return (order+1)*(order+2)/2;
202}
203static
204size_t
205reference_element_q_first_inod_by_variant (
206 size_t order,
207 size_t subgeo_variant)
208{
209 if (subgeo_variant == reference_element__p) return 0;
210 if (subgeo_variant == reference_element__e) return 4;
211 if (subgeo_variant == reference_element__t) return 4 + 4*(order-1);
212 if (subgeo_variant == reference_element__q) return 4 + 4*(order-1);
213 return (order+1)*(order+1);
214}
215static
216size_t
217reference_element_T_first_inod_by_variant (
218 size_t order,
219 size_t subgeo_variant)
220{
221 if (subgeo_variant == reference_element__p) return 0;
222 if (subgeo_variant == reference_element__e) return 4;
223 if (subgeo_variant == reference_element__t) return 4 + 6*(order-1);
224 if (subgeo_variant == reference_element__q) return 4 + 6*(order-1) + 4*((order-1)*(order-2)/2);
225 if (subgeo_variant == reference_element__T) return 4 + 6*(order-1) + 4*((order-1)*(order-2)/2);
226 return (order+1)*(order+2)*(order+3)/6;
227}
228static
229size_t
230reference_element_P_first_inod_by_variant (
231 size_t order,
232 size_t subgeo_variant)
233{
234 if (subgeo_variant == reference_element__p) return 0;
235 if (subgeo_variant == reference_element__e) return 6;
236 if (subgeo_variant == reference_element__t) return 6 + 9*(order-1);
237 if (subgeo_variant == reference_element__q) return 6 + 9*(order-1) + 2*((order-1)*(order-2)/2);
238 if (subgeo_variant == reference_element__T) return 6 + 9*(order-1) + 2*((order-1)*(order-2)/2) + 3*(order-1)*(order-1);
239 if (subgeo_variant == reference_element__P) return 6 + 9*(order-1) + 2*((order-1)*(order-2)/2) + 3*(order-1)*(order-1);
240 return ((order+1)*(order+1)*(order+2))/2;
241}
242static
243size_t
244reference_element_H_first_inod_by_variant (
245 size_t order,
246 size_t subgeo_variant)
247{
248 if (subgeo_variant == reference_element__p) return 0;
249 if (subgeo_variant == reference_element__e) return 8;
250 if (subgeo_variant == reference_element__t) return 8 + 12*(order-1);
251 if (subgeo_variant == reference_element__q) return 8 + 12*(order-1);
252 if (subgeo_variant == reference_element__T) return 8 + 12*(order-1) + 6*(order-1)*(order-1);
253 if (subgeo_variant == reference_element__P) return 8 + 12*(order-1) + 6*(order-1)*(order-1);
254 if (subgeo_variant == reference_element__H) return 8 + 12*(order-1) + 6*(order-1)*(order-1);
255 return (order+1)*(order+1)*(order+1);
256}
257static
258size_t
259reference_element_first_inod_by_variant (
260 size_t variant,
261 size_t order,
262 size_t subgeo_variant)
263{
264#define _RHEOLEF_geo_element_auto_case(NAME) \
265 case reference_element__##NAME: \
266 return reference_element_##NAME##_first_inod_by_variant (order, subgeo_variant);
267
268 switch (variant) {
276 default: error_macro ("unexpected element variant `"<<variant<<"'"); return 0;
277 }
278#undef _RHEOLEF_geo_element_auto_case
279}
280static
281size_t
282reference_element_last_inod_by_variant (
283 size_t variant,
284 size_t order,
285 size_t subgeo_variant)
286{
287 return reference_element_first_inod_by_variant (variant, order, subgeo_variant+1);
288}
289// -----------------------------------------------------------------------------
290// reference_element::ilat2loc_inod
291// -----------------------------------------------------------------------------
292static
293size_t
294reference_element_p_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
295{
296 return 0;
297}
298static
299size_t
300reference_element_e_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
301{
302 size_t i = ilat[0];
303 if (i == 0) return 0;
304 if (i == order) return 1;
305 return 2 + (i - 1);
306}
307static
308size_t
309reference_element_t_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
310{
311 size_t i = ilat[0];
312 size_t j = ilat[1];
313 if (j == 0) { // horizontal edge [0,1]
314 if (i == 0) return 0;
315 if (i == order) return 1;
316 return 3 + 0*(order-1) + (i - 1);
317 }
318 if (i + j == order) { // oblique edge ]1,2]
319 if (j == order) return 2;
320 return 3 + 1*(order-1) + (j - 1);
321 }
322 size_t j1 = order - j;
323 if (i == 0) { // vertical edge ]2,0[
324 return 3 + 2*(order-1) + (j1 - 1);
325 }
326 // internal face ]0,1[x]0,2[: i, then j
327 size_t n_face_node = (order-1)*(order-2)/2;
328 size_t ij = (n_face_node - (j1-1)*j1/2) + (i - 1);
329 return 3 + 3*(order-1) + ij;
330}
331static
332size_t
333reference_element_q_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
334{
335 size_t i = ilat[0];
336 size_t j = ilat[1];
337 if (j == 0) { // horizontal edge [0,1]
338 if (i == 0) return 0;
339 if (i == order) return 1;
340 return 4 + 0*(order-1) + (i - 1);
341 }
342 if (i == order) { // vertical edge ]1,2]
343 if (j == order) return 2;
344 return 4 + 1*(order-1) + (j - 1);
345 }
346 size_t i1 = order - i;
347 if (j == order) { // horizontal edge ]2,3]
348 if (i == 0) return 3;
349 return 4 + 2*(order-1) + (i1 - 1);
350 }
351 size_t j1 = order - j;
352 if (i == 0) { // vertical edge ]3,0[
353 return 4 + 3*(order-1) + (j1 - 1);
354 }
355 // internal face ]0,1[x]0,3[: i, then j
356 size_t ij = (order-1)*(j-1) + (i-1);
357 size_t n_bdry_node = 4 + 4*(order-1);
358 return n_bdry_node + ij;
359}
360static
361size_t
362reference_element_T_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
363{
364 size_t i = ilat[0];
365 size_t j = ilat[1];
366 size_t k = ilat[2];
367 size_t n_tot_edge_node = 6*(order-1);
368 size_t n_face_node = (order-1)*(order-2)/2;
369 if (k == 0) {
370 // horizontal face 021 = [0,2]x[0,1]
371 if (j == 0) { // left edge [0,1]
372 if (i == 0) return 0;
373 if (i == order) return 1;
374 return 4 + 0*(order-1) + (i - 1);
375 }
376 if (i + j == order) { // oblique edge ]1,2[
377 if (j == order) return 2;
378 return 4 + 1*(order-1) + (j - 1);
379 }
380 size_t j1 = order - j;
381 if (i == 0) { // right edge ]0,2[
382 return 4 + 2*(order-1) + (j1 - 1);
383 }
384 // internal face 021 = ]0,2[x]0,1[: j, then i
385 size_t i1 = order - i;
386 size_t ji = (n_face_node - i1*(i1-1)/2) + (j - 1);
387 return 4 + n_tot_edge_node + 0*n_face_node + ji;
388 }
389 if (i == 0) {
390 // back face 032 =]0,3]x]0,2[
391 if (j == 0) { // vertical edge ]0,3]
392 if (k == order) return 3;
393 return 4 + 3*(order-1) + (k - 1);
394 }
395 if (j + k == order) { // oblique edge ]2,3[
396 return 4 + 5*(order-1) + (k - 1);
397 }
398 // internal face 032 =]0,3]x]0,2[, k, then j
399 size_t j1 = order - j;
400 size_t kj = (n_face_node - j1*(j1-1)/2) + (k - 1);
401 return 4 + n_tot_edge_node + 1*n_face_node + kj;
402 }
403 if (j == 0) {
404 // left face 013 = ]0,1[x]0,3[
405 if (i + k == order) { // oblique edge ]1,3[
406 return 4 + 4*(order-1) + (k - 1);
407 }
408 // internal face 013 = ]0,1[x]0,3[: i, then k
409 size_t k1 = order - k;
410 size_t ik = (n_face_node - k1*(k1-1)/2) + (i - 1);
411 return 4 + n_tot_edge_node + 2*n_face_node + ik;
412 }
413 if (i + j + k == order) {
414 // internal to oblique face 123 = ]1,2[x]1,3], j, then k
415 size_t k1 = order - k;
416 size_t jk = (n_face_node - k1*(k1-1)/2) + (j - 1);
417 return 4 + n_tot_edge_node + 3*n_face_node + jk;
418 }
419 // internal to volume: i, then j, then k
420 size_t n_tot_face_node = 4*n_face_node;
421 size_t n_tot_node = (order+1)*(order+2)*(order+3)/6;
422 size_t n_volume_node = (order-1)*(order-2)*(order-3)/6;
423 size_t k1 = order - k;
424 size_t j1 = order - j - k;
425 size_t n_level_k_node = (k1-1)*(k1-2)/2;
426 size_t ijk = (n_volume_node - k1*(k1-1)*(k1-2)/6)
427 + (n_level_k_node - j1*(j1-1)/2)
428 + (i - 1);
429 return 4 + n_tot_edge_node + n_tot_face_node + ijk;
430}
431static
432size_t
433reference_element_P_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
434{
435 size_t i = ilat[0];
436 size_t j = ilat[1];
437 size_t k = ilat[2];
438 size_t i1 = order - i;
439 size_t j1 = order - j;
440 size_t k1 = order - k;
441 size_t n_node_edge = order-1;
442 size_t n_tot_edge_node = 9*(order-1);
443 size_t n_node_face_tri = (order-1)*(order-2)/2;
444 size_t n_node_face_qua = (order-1)*(order-1);
445 if (k == 0) {
446 // horizontal face 021 = [0,1]x[0,2]
447 if (j == 0) { // left edge [0,1]
448 if (i == 0) return 0;
449 if (i == order) return 1;
450 return 6 + 0*(order-1) + (i - 1);
451 }
452 if (i + j == order) { // oblique edge ]1,2[
453 if (j == order) return 2;
454 return 6 + 1*(order-1) + (j - 1);
455 }
456 if (i == 0) { // right edge ]0,2[
457 return 6 + 2*(order-1) + (j1 - 1);
458 }
459 // internal face 021 = ]0,2[x]0,1[: j, then i
460 size_t ji = (n_node_face_tri - i1*(i1-1)/2) + (j - 1);
461 return 6 + n_tot_edge_node + 0*n_node_face_tri + ji;
462 }
463 if (k == order) {
464 // horizontal face 354 = [3,4]x[3,5]
465 if (j == 0) { // left edge [3,4]
466 if (i == 0) return 3;
467 if (i == order) return 4;
468 return 6 + 6*(order-1) + (i - 1);
469 }
470 if (i + j == order) { // oblique edge ]4,5[
471 if (j == order) return 5;
472 return 6 + 7*(order-1) + (j - 1);
473 }
474 if (i == 0) { // right edge ]3,5[
475 return 6 + 8*(order-1) + (j1 - 1);
476 }
477 // internal face 354 = ]3,4[x]3,5[: i, then j
478 size_t ij = (n_node_face_tri - j1*(j1-1)/2) + (i - 1);
479 return 6 + n_tot_edge_node + 1*n_node_face_tri + ij;
480 }
481 if (j == 0) {
482 // left face 0143 = [0,1]x]0,3[
483 if (i == 0) {
484 // vertical edge ]0,3[
485 return 6 + 3*(order-1) + (k - 1);
486 }
487 if (i == order) {
488 // vertical edge ]1,4[
489 return 6 + 4*(order-1) + (k - 1);
490 }
491 // internal face 01243 = ]0,1[x]0,3[: i then k
492 size_t ik = n_node_edge*(k-1) + (i-1);
493 return 6 + n_tot_edge_node + 2*n_node_face_tri + 0*n_node_face_qua + ik;
494 }
495 if (i == 0) {
496 // back face 0352 = ]0,2]x]0,3[
497 if (j == order) {
498 // vertical edge ]2,5[
499 return 6 + 5*(order-1) + (k - 1);
500 }
501 // internal face 0352 = ]0,2[x]0,3[: k then j
502 size_t kj = n_node_edge*(j-1) + (k-1);
503 return 6 + n_tot_edge_node + 2*n_node_face_tri + 2*n_node_face_qua + kj;
504 }
505 if (i + j == order) {
506 // internal face 1254 = ]1,2[x]1,4[: j then k
507 size_t jk = n_node_edge*(k-1) + (j-1);
508 return 6 + n_tot_edge_node + 2*n_node_face_tri + 1*n_node_face_qua + jk;
509 }
510 // internal volume ]0,1[x]0,2[x][0,3[: i, then j, then k
511 size_t n_tot_face_node = 2*n_node_face_tri + 3*n_node_face_qua;
512 size_t ij = (n_node_face_tri - (j1-1)*j1/2) + (i - 1);
513 size_t ijk = n_node_face_tri*(k-1) + ij;
514 return 6 + n_tot_edge_node + n_tot_face_node + ijk;
515}
516static
517size_t
518reference_element_H_ilat2loc_inod (size_t order, const point_basic<size_t>& ilat)
519{
520 size_t i = ilat[0];
521 size_t j = ilat[1];
522 size_t k = ilat[2];
523 size_t i1 = order - i;
524 size_t j1 = order - j;
525 size_t k1 = order - k;
526 size_t n_node_edge = order-1;
527 size_t n_node_face = (order-1)*(order-1);
528 if (k == 0) { // bottom face [0,3]x[0,1]
529 if (j == 0) { // bottom-left edge [0,1]
530 if (i == 0) return 0;
531 if (i == order) return 1;
532 return 8 + 0*n_node_edge + (i - 1);
533 }
534 if (j == order) { // bottom-right edge [3,2]
535 if (i == 0) return 3;
536 if (i == order) return 2;
537 return 8 + 2*n_node_edge + (i1 - 1);
538 }
539 if (i == 0) { // bottom-back edge ]3,0[
540 return 8 + 3*n_node_edge + (j1 - 1);
541 }
542 if (i == order) { // bottom-front edge ]1,2[
543 return 8 + 1*n_node_edge + (j - 1);
544 }
545 // internal face ]0,3[x]0,1[: j then i (TODO: not checked since Pk(H) only with k <= 2 yet
546 size_t ji = n_node_edge*(i-1) + (j-1);
547 return 8 + 12*n_node_edge + 0*n_node_face + ji;
548 }
549 if (k == order) { // top face [4,5]x[4,7]
550 if (j == 0) { // top-left edge [4,5]
551 if (i == 0) return 4;
552 if (i == order) return 5;
553 return 8 + 8*n_node_edge + (i - 1);
554 }
555 if (j == order) { // top-right edge [7,6]
556 if (i == 0) return 7;
557 if (i == order) return 6;
558 return 8 + 10*n_node_edge + (i1 - 1);
559 }
560 if (i == 0) { // top-back edge ]7,4[
561 return 8 + 11*n_node_edge + (j1 - 1);
562 }
563 if (i == order) { // top-front edge ]5,6[
564 return 8 + 9*n_node_edge + (j - 1);
565 }
566 // internal face ]4,5[x]4,7[: i then j
567 size_t ij = n_node_edge*(j-1) + (i-1);
568 return 8 + 12*n_node_edge + 3*n_node_face + ij;
569 }
570 if (j == 0) { // left face ]0,1[x[0,4]
571 if (i == 0) { // left-back edge [0,4]
572 return 8 + 4*n_node_edge + (k - 1);
573 }
574 if (i == order) { // left-front edge [1,5]
575 return 8 + 5*n_node_edge + (k - 1);
576 }
577 // internal face ]0,1[x]0,4[: i then k
578 size_t ik = n_node_edge*(k-1) + (i-1);
579 return 8 + 12*n_node_edge + 2*n_node_face + ik;
580 }
581 if (j == order) { // right face ]2,3[x[2,6]
582 if (i == 0) { // right-back edge [3,7]
583 return 8 + 7*n_node_edge + (k - 1);
584 }
585 if (i == order) { // right-front edge [2,6]
586 return 8 + 6*n_node_edge + (k - 1);
587 }
588 // internal face ]2,3[x]2,6[: i1 then k
589 size_t i1k = n_node_edge*(k-1) + (i1-1);
590 return 8 + 12*n_node_edge + 5*n_node_face + i1k;
591 }
592 if (i == 0) { // internal back face ]0,4[x[0,3]: k then j
593 size_t kj = n_node_edge*(j-1) + (k-1);
594 return 8 + 12*n_node_edge + 1*n_node_face + kj;
595 }
596 if (i == order) { // internal front face ]1,2[x[1,5]: j then k
597 size_t jk = n_node_edge*(k-1) + (j-1);
598 return 8 + 12*n_node_edge + 4*n_node_face + jk;
599 }
600 // internal volume: i then j then k
601 size_t n_node_bdry = 8 + 12*n_node_edge + 6*n_node_face;
602 size_t ijk = (order-1)*((order-1)*(k-1) + (j-1)) + (i-1);
603 return n_node_bdry + ijk;
604}
605
606} // namespace rheolef
607#endif // _RHEOLEF_REFERENCE_ELEMENT_AUX_ICC
static const variant_type H
static const variant_type q
static const variant_type e
static const variant_type max_variant
static const variant_type p
static const variant_type T
static const variant_type P
static const variant_type t
constexpr size_t reference_element__H
constexpr size_t reference_element__e
constexpr size_t reference_element__P
static const size_t _first_variant_by_dimension[5]
static size_t _reference_element_n_face_by_variant[reference_element__max_variant]
constexpr size_t reference_element__T
constexpr size_t reference_element__p
constexpr size_t reference_element__t
constexpr size_t reference_element__q
constexpr size_t reference_element__max_variant
static size_t _reference_element_n_edge_by_variant[reference_element__max_variant]
static size_t _reference_element_dimension_by_variant[reference_element__max_variant]
#define error_macro(message)
Definition dis_macros.h:49
Expr1::float_type T
Definition field_expr.h:230
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
This file is part of Rheolef.
#define _RHEOLEF_geo_element_auto_case(VARIANT)
Definition sphere.icc:25