GetFEM  5.4.2
bgeot_poly_composite.cc
1 /*===========================================================================
2 
3  Copyright (C) 2002-2020 Yves Renard
4 
5  This file is a part of GetFEM
6 
7  GetFEM is free software; you can redistribute it and/or modify it
8  under the terms of the GNU Lesser General Public License as published
9  by the Free Software Foundation; either version 3 of the License, or
10  (at your option) any later version along with the GCC Runtime Library
11  Exception either version 3.1 or (at your option) any later version.
12  This program is distributed in the hope that it will be useful, but
13  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15  License and GCC Runtime Library Exception for more details.
16  You should have received a copy of the GNU Lesser General Public License
17  along with this program; if not, write to the Free Software Foundation,
18  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
19 
20 ===========================================================================*/
21 
22 #include "getfem/dal_singleton.h"
25 
26 namespace bgeot {
27 
28  inline scalar_type sfloor(scalar_type x)
29  { return (x >= 0) ? floor(x) : -floor(-x); }
30 
31 
33  const base_node &y) const {
34  size_type s = x.size();
35  scalar_type c1 = c_max, c2 = c_max * scalar_type(base);
36  GMM_ASSERT2(y.size() == s, "dimension error");
37 
38  base_node::const_iterator itx=x.begin(), itex=x.end(), ity=y.begin();
39  int ret = 0;
40  for (; itx != itex; ++itx, ++ity) {
41  long a = long(sfloor((*itx) * c1)), b = long(sfloor((*ity) * c1));
42  if ((scalar_type(gmm::abs(a)) > scalar_type(base))
43  || (scalar_type(gmm::abs(b)) > scalar_type(base))) {
44  exp_max++; c_max /= scalar_type(base);
45  return (*this)(x,y);
46  }
47  if (ret == 0) { if (a < b) ret = -1; else if (a > b) ret = 1; }
48  }
49  if (ret) return ret;
50 
51  for (int e = exp_max; e >= exp_min; --e, c1 *= scalar_type(base),
52  c2 *= scalar_type(base)) {
53  itx = x.begin(), itex = x.end(), ity = y.begin();
54  for (; itx != itex; ++itx, ++ity) {
55  int a = int(sfloor(((*itx) * c2) - sfloor((*itx) * c1)
56  * scalar_type(base)));
57  int b = int(sfloor(((*ity) * c2) - sfloor((*ity) * c1)
58  * scalar_type(base)));
59  if (a < b) return -1; else if (a > b) return 1;
60  }
61  }
62  return 0;
63  }
64 
65  mesh_precomposite::mesh_precomposite(const basic_mesh &m) : box_tree{1e-13} {
66  initialise(m);
67  }
68 
69  void mesh_precomposite::initialise(const basic_mesh &m) {
70  msh = &m;
71  det.resize(m.nb_convex()); orgs.resize(m.nb_convex());
72  gtrans.resize(m.nb_convex()); gtransinv.resize(m.nb_convex());
73  for (size_type i = 0; i <= m.points().index().last_true(); ++i)
74  vertices.add(m.points()[i]);
75 
76  base_node min, max;
77  for (dal::bv_visitor cv(m.convex_index()); !cv.finished(); ++cv) {
78 
79  pgeometric_trans pgt = m.trans_of_convex(cv);
80  size_type N = pgt->structure()->dim();
81  size_type P = m.dim();
82  GMM_ASSERT1(pgt->is_linear(), "Bad geometric transformation");
83 
84  base_matrix G(P, pgt->nb_points());
85  base_node X(N);
86 
87  m.points_of_convex(cv, G);
89  gtrans[cv] = ctx.K();
90  gtransinv[cv] = ctx.B();
91  det[cv] = ctx.J();
92 
93  auto points_of_convex = m.points_of_convex(cv);
94  orgs[cv] = pgt->transform(X, G);
95  bounding_box(min, max, points_of_convex);
96  box_to_convexes_map[box_tree.add_box(min, max)].push_back(cv);
97  }
98 
99  box_tree.build_tree();
100  }
101 
102  DAL_TRIPLE_KEY(base_poly_key, short_type, short_type,
103  std::vector<opt_long_scalar_type>);
104 
105  polynomial_composite::polynomial_composite(const mesh_precomposite &m,
106  bool lc, bool ff)
107  : mp(&m), local_coordinate(lc), faces_first(ff),
108  default_polys(mp->dim()+1) {
109  for (dim_type i = 0; i <= mp->dim(); ++i)
110  default_polys[i] = base_poly(i, 0.);
111  }
112 
113  static void mult_diff_transposed(const base_matrix &M, const base_node &p,
114  const base_node &p1, base_node &p2) {
115  for (dim_type d = 0; d < p2.size(); ++d) {
116  p2[d] = 0;
117  auto col = mat_col(M, d);
118  for (dim_type i = 0; i < p1.size(); ++i)
119  p2[d] += col[i] * (p[i] - p1[i]);
120  }
121  }
122 
123  scalar_type polynomial_composite::eval(const base_node &p,
124  size_type l) const {
125  if (l != size_type(-1)) {
126  if (!local_coordinate) return poly_of_subelt(l).eval(p.begin());
127  base_node p1(gmm::mat_ncols(mp->gtransinv[l]));
128  mult_diff_transposed(mp->gtransinv[l], p, mp->orgs[l], p1);
129  return poly_of_subelt(l).eval(p1.begin());
130  }
131 
132  auto dim = mp->dim();
133  base_node p1_stored, p1, p2(dim);
134  size_type cv_stored(-1);
135 
136  auto &box_tree = mp->box_tree;
137  rtree::pbox_set box_list;
138  box_tree.find_boxes_at_point(p, box_list);
139 
140  while (box_list.size()) {
141  auto pmax_box = *box_list.begin();
142 
143  if (box_list.size() > 1) {
144  auto max_rate = -1.0;
145  for (auto &&pbox : box_list) {
146  auto box_rate = 1.0;
147  for (size_type i = 0; i < dim; ++i) {
148  auto h = pbox->max->at(i) - pbox->min->at(i);
149  if (h > 0) {
150  auto rate = std::min(pbox->max->at(i) - p[i],
151  p[i] - pbox->min->at(i)) / h;
152  box_rate = std::min(rate, box_rate);
153  }
154  }
155  if (box_rate > max_rate) { pmax_box = pbox; max_rate = box_rate; }
156  }
157  }
158 
159  for (auto cv : mp->box_to_convexes_map.at(pmax_box->id)) {
160  dim_type elt_dim = dim_type(gmm::mat_ncols(mp->gtransinv[cv]));
161  p1.resize(elt_dim);
162  mult_diff_transposed(mp->gtransinv[cv], p, mp->orgs[cv], p1);
163  scalar_type ddist(0);
164 
165  if (elt_dim < dim) {
166  p2 = mp->orgs[cv]; gmm::mult_add(mp->gtrans[cv], p1, p2);
167  ddist = gmm::vect_dist2(p2, p);
168  }
169  if (mp->trans_of_convex(cv)->convex_ref()->is_in(p1) < 1E-9
170  && ddist < 1E-9) {
171  if (!faces_first || elt_dim < dim) {
172  scalar_type res = to_scalar(poly_of_subelt(cv).eval(local_coordinate
173  ? p1.begin() : p.begin()));
174  return res;
175  }
176  p1_stored = p1; cv_stored = cv;
177  }
178  }
179 
180  if (box_list.size() == 1) break;
181  box_list.erase(pmax_box);
182  }
183  if (cv_stored != size_type(-1)) {
184  scalar_type res =
185  to_scalar(poly_of_subelt(cv_stored).eval(local_coordinate
186  ? p1_stored.begin(): p.begin()));
187  return res;
188  }
189  GMM_ASSERT1(false, "Element not found in composite polynomial: " << p);
190  }
191 
192 
193  void polynomial_composite::derivative(short_type k) {
194  if (local_coordinate) {
195  dim_type P = mp->dim();
196  base_vector e(P), f(P); e[k] = 1.0;
197  for (size_type ic = 0; ic < mp->nb_convex(); ++ic) {
198  dim_type N = dim_type(gmm::mat_ncols(mp->gtransinv[ic]));
199  f.resize(N);
200  gmm::mult(gmm::transposed(mp->gtransinv[ic]), e, f);
201  base_poly Der(N, 0), Q;
202  if (polytab.find(ic) != polytab.end()) {
203  auto &poly = poly_of_subelt(ic);
204  for (dim_type n = 0; n < N; ++n) {
205  Q = poly;
206  Q.derivative(n);
207  Der += Q * f[n];
208  }
209  if (Der.is_zero()) polytab.erase(ic); else set_poly_of_subelt(ic,Der);
210  }
211  }
212  }
213  else
214  for (size_type ic = 0; ic < mp->nb_convex(); ++ic) {
215  auto poly = poly_of_subelt(ic);
216  poly.derivative(k);
217  if (polytab.find(ic) != polytab.end()) set_poly_of_subelt(ic, poly);
218  }
219  }
220 
221  void polynomial_composite::set_poly_of_subelt(size_type l,
222  const base_poly &poly) {
223  auto poly_key =
224  std::make_shared<base_poly_key>(poly.degree(), poly.dim(), poly);
225  pstored_base_poly o(std::dynamic_pointer_cast<const stored_base_poly>
226  (dal::search_stored_object(poly_key)));
227  if (!o) {
228  o = std::make_shared<stored_base_poly>(poly);
229  dal::add_stored_object(poly_key, o);
230  }
231  polytab[l] = o;
232  }
233 
234  const base_poly &polynomial_composite::poly_of_subelt(size_type l) const {
235  auto it = polytab.find(l);
236  if (it == polytab.end()) {
237  if (local_coordinate)
238  return default_polys[gmm::mat_ncols(mp->gtransinv[l])];
239  else
240  return default_polys[mp->dim()];
241  }
242  return *(it->second);
243  }
244 
245 
246  /* build a regularly refined mesh for a simplex of dimension <= 3.
247  All simplexes have the same orientation as the original simplex.
248  */
249  static void
250  structured_mesh_for_simplex_(pconvex_structure cvs,
251  pgeometric_trans opt_gt,
252  const std::vector<base_node> *opt_gt_pts,
253  short_type k, pbasic_mesh pm) {
254  scalar_type h = 1./k;
255  switch (cvs->dim()) {
256  case 1 :
257  {
258  base_node a(1), b(1);
259  for (short_type i = 0; i < k; ++i) {
260  a[0] = i * h; b[0] = (i+1) * h;
261  if (opt_gt) a = opt_gt->transform(a, *opt_gt_pts);
262  if (opt_gt) b = opt_gt->transform(b, *opt_gt_pts);
263  size_type na = pm->add_point(a);
264  size_type nb = pm->add_point(b);
265  pm->add_segment(na, nb);
266  }
267  }
268  break;
269  case 2 :
270  {
271  base_node A(2),B(2),C(2),D(2);
272  for (short_type i = 0; i < k; ++i) {
273  scalar_type a = i * h, b = (i+1) * h;
274  for (short_type l = 0; l+i < k; ++l) {
275  scalar_type c = l * h, d = (l+1) * h;
276  A[0] = a; A[1] = c;
277  B[0] = b; B[1] = c;
278  C[0] = a; C[1] = d;
279  D[0] = b; D[1] = d;
280  if (opt_gt) {
281  A = opt_gt->transform(A, *opt_gt_pts);
282  B = opt_gt->transform(B, *opt_gt_pts);
283  C = opt_gt->transform(C, *opt_gt_pts);
284  D = opt_gt->transform(D, *opt_gt_pts);
285  }
286  /* add triangle with correct orientation (det [B-A;C-A] > 0) */
287  size_type nA = pm->add_point(A);
288  size_type nB = pm->add_point(B);
289  size_type nC = pm->add_point(C);
290  size_type nD = pm->add_point(D);
291  pm->add_triangle(nA,nB,nC);
292  if (l+i+1 < k)
293  pm->add_triangle(nC,nB,nD);
294  }
295  }
296  }
297  break;
298  case 3 :
299  {
300  /* based on decompositions of small cubes
301  the number of tetrahedrons is k^3
302  */
303  base_node A,B,C,D,E,F,G,H;
304  for (short_type ci = 0; ci < k; ci++) {
305  scalar_type x = ci*h;
306  for (short_type cj = 0; cj < k-ci; cj++) {
307  scalar_type y=cj*h;
308  for (short_type ck = 0; ck < k-ci-cj; ck++) {
309  scalar_type z=ck*h;
310 
311  A = {x, y, z};
312  B = {x+h, y, z};
313  C = {x, y+h, z};
314  D = {x+h, y+h, z};
315  E = {x, y, z+h};
316  F = {x+h, y, z+h};
317  G = {x, y+h, z+h};
318  H = {x+h, y+h, z+h};
319  if (opt_gt) {
320  A = opt_gt->transform(A, *opt_gt_pts);
321  B = opt_gt->transform(B, *opt_gt_pts);
322  C = opt_gt->transform(C, *opt_gt_pts);
323  D = opt_gt->transform(D, *opt_gt_pts);
324  E = opt_gt->transform(E, *opt_gt_pts);
325  F = opt_gt->transform(F, *opt_gt_pts);
326  G = opt_gt->transform(G, *opt_gt_pts);
327  H = opt_gt->transform(H, *opt_gt_pts);
328  }
329  size_type t[8];
330  t[0] = pm->add_point(A);
331  t[1] = pm->add_point(B);
332  t[2] = pm->add_point(C);
333  t[4] = pm->add_point(E);
334  t[3] = t[5] = t[6] = t[7] = size_type(-1);
335  if (k > 1 && ci+cj+ck < k-1) {
336  t[3] = pm->add_point(D);
337  t[5] = pm->add_point(F);
338  t[6] = pm->add_point(G);
339  }
340  if (k > 2 && ci+cj+ck < k-2) {
341  t[7] = pm->add_point(H);
342  }
343  /**
344  Note that the orientation of each tetrahedron is the same
345  (det > 0)
346  */
347  pm->add_tetrahedron(t[0], t[1], t[2], t[4]);
348  if (k > 1 && ci+cj+ck < k-1) {
349  pm->add_tetrahedron(t[1], t[2], t[4], t[5]);
350  pm->add_tetrahedron(t[6], t[4], t[2], t[5]);
351  pm->add_tetrahedron(t[2], t[3], t[5], t[1]);
352  pm->add_tetrahedron(t[2], t[5], t[3], t[6]);
353  }
354  if (k > 2 && ci+cj+ck < k-2) {
355  pm->add_tetrahedron(t[3], t[5], t[7], t[6]);
356  }
357  }
358  }
359  }
360  }
361  break;
362  default :
363  GMM_ASSERT1(false, "Sorry, not implemented for simplices of "
364  "dimension " << int(cvs->dim()));
365  }
366  }
367 
368  static void structured_mesh_for_parallelepiped_
369  (pconvex_structure cvs, pgeometric_trans opt_gt,
370  const std::vector<base_node> *opt_gt_pts, short_type k, pbasic_mesh pm) {
371  scalar_type h = 1./k;
372  size_type n = cvs->dim();
373  size_type pow2n = (size_type(1) << n);
374  size_type powkn = 1; for (size_type i=0; i < n; ++i) powkn *= k;
375  std::vector<size_type> strides(n);
376  size_type nbpts = 1;
377  for (size_type i=0; i < n; ++i) { strides[i] = nbpts; nbpts *= (k+1); }
378  std::vector<short_type> kcnt; kcnt.resize(n+1,0);
379  std::vector<size_type> pids; pids.reserve(nbpts);
380  base_node pt(n);
381  size_type kk;
382  /* insert nodes and fill pids with their numbers */
383  while (kcnt[n] == 0) {
384  for (size_type z = 0; z < n; ++z)
385  pt[z] = h*kcnt[z];
386  if (opt_gt) pt = opt_gt->transform(pt, *opt_gt_pts);
387  pids.push_back(pm->add_point(pt));
388  kk=0;
389  while (kk <= n)
390  { kcnt[kk]++; if (kcnt[kk] == k+1) { kcnt[kk] = 0; kk++; } else break; }
391  }
392 
393  /*
394  insert convexes using node ids stored in 'pids'
395  */
396  std::vector<size_type> ppts(pow2n);
397  kcnt[n] = 0;
398  while (kcnt[n] == 0) {
399  for (kk = 0; kk < pow2n; ++kk) {
400  size_type pos = 0;
401  for (size_type z = 0; z < n; ++z) {
402  pos += kcnt[z]*strides[z];
403  if ((kk & (size_type(1) << z))) pos += strides[z];
404  }
405  ppts[kk] = pids.at(pos);
406  }
407  pm->add_convex(parallelepiped_linear_geotrans(n), ppts.begin());
408  kcnt[(kk = 0)]++;
409  while (kk < n && kcnt[kk] == k) { kcnt[kk] = 0; kcnt[++kk]++; }
410  }
411  }
412 
413  static void
414  structured_mesh_for_convex_(pconvex_structure cvs,
415  pgeometric_trans opt_gt,
416  const std::vector<base_node> *opt_gt_pts,
417  short_type k, pbasic_mesh pm) {
418  size_type nbp = basic_structure(cvs)->nb_points();
419  dim_type n = cvs->dim();
420  /* Identifying simplexes. */
421  if (nbp == size_type(n+1) &&
422  *key_of_stored_object(basic_structure(cvs))
423  ==*key_of_stored_object(simplex_structure(n))) {
424  // smc.pm->write_to_file(cout);
425  structured_mesh_for_simplex_(cvs,opt_gt,opt_gt_pts,k,pm);
426  /* Identifying parallelepipeds. */
427  } else if (nbp == (size_type(1) << n) &&
428  *key_of_stored_object(basic_structure(cvs))
429  == *key_of_stored_object(parallelepiped_structure(n))) {
430  structured_mesh_for_parallelepiped_(cvs,opt_gt,opt_gt_pts,k,pm);
431  } else if (nbp == size_type(2 * n) &&
432  *key_of_stored_object(basic_structure(cvs))
433  == *key_of_stored_object(prism_P1_structure(n))) {
434  GMM_ASSERT1(false, "Sorry, structured_mesh not implemented for prisms.");
435  } else {
436  GMM_ASSERT1(false, "This element is not taken into account. Contact us");
437  }
438  }
439 
440  /* extract the mesh_structure on faces */
441  static void structured_mesh_of_faces_(pconvex_ref cvr, dim_type f,
442  const basic_mesh &m,
443  mesh_structure &facem) {
444  facem.clear();
445  dal::bit_vector on_face;
446  for (size_type i = 0; i < m.nb_max_points(); ++i) {
447  if (m.is_point_valid(i)) {
448  if (gmm::abs(cvr->is_in_face(f, m.points()[i])) < 1e-12)
449  on_face.add(i);
450  }
451  }
452  //cerr << "on_face=" << on_face << endl;
453  for (dal::bv_visitor cv(m.convex_index()); !cv.finished(); ++cv) {
454  for (short_type ff=0; ff < m.structure_of_convex(cv)->nb_faces(); ++ff) {
455  mesh_structure::ind_pt_face_ct
456  ipts=m.ind_points_of_face_of_convex(cv,ff);
457  bool allin = true;
458  for (size_type i=0; i < ipts.size(); ++i) if (!on_face[ipts[i]])
459  { allin = false; break; }
460  if (allin) {
461  /* cerr << "ajout de la face " << ff << " du convexe " << cv << ":";
462  for (size_type i=0; i < ipts.size(); ++i)
463  cerr << on_face[ipts[i]] << "/" << ipts[i] << " ";
464  cerr << endl;
465  */
466  facem.add_convex(m.structure_of_convex(cv)->faces_structure()[ff],
467  ipts.begin());
468  }
469  }
470  }
471  }
472 
473  DAL_TRIPLE_KEY(str_mesh_key, pconvex_structure, short_type, bool);
474 
475  struct str_mesh_cv__ : virtual public dal::static_stored_object {
476  pconvex_structure cvs;
477  short_type n;
478  bool simplex_mesh; /* true if the convex has been splited into simplexes,
479  which were refined */
480  std::unique_ptr<basic_mesh> pm;
481  std::vector<std::unique_ptr<mesh_structure>> pfacem; /* array of mesh_structures for faces */
482  dal::bit_vector nodes_on_edges;
483  std::unique_ptr<mesh_precomposite> pmp;
484  str_mesh_cv__(pconvex_structure c, short_type k, bool smesh_) :
485  cvs(c), n(k), simplex_mesh(smesh_)
486  { DAL_STORED_OBJECT_DEBUG_CREATED(this, "cv mesh"); }
487  ~str_mesh_cv__() { DAL_STORED_OBJECT_DEBUG_DESTROYED(this,"cv mesh"); }
488  };
489 
490  typedef std::shared_ptr<const str_mesh_cv__> pstr_mesh_cv__;
491 
492  /**
493  * This function returns a mesh in pm which contains a refinement of the convex
494  * cvr if force_simplexification is false, refined convexes have the same
495  * basic_structure than cvr, if it is set to true, the cvr is decomposed
496  * into simplexes which are then refined.
497  * TODO: move it into another file and separate the pmesh_precomposite part ?
498  **/
499  void structured_mesh_for_convex(pconvex_ref cvr, short_type k,
500  pbasic_mesh &pm, pmesh_precomposite &pmp,
501  bool force_simplexification) {
502  size_type n = cvr->structure()->dim();
503  size_type nbp = basic_structure(cvr->structure())->nb_points();
504 
505  force_simplexification = force_simplexification || (nbp == n+1);
506  dal::pstatic_stored_object_key
507  pk = std::make_shared<str_mesh_key>(basic_structure(cvr->structure()),
508  k, force_simplexification);
509 
510  dal::pstatic_stored_object o = dal::search_stored_object_on_all_threads(pk);
511  pstr_mesh_cv__ psmc;
512  if (o)
513  psmc = std::dynamic_pointer_cast<const str_mesh_cv__>(o);
514  else {
515 
516  auto ppsmc = std::make_shared<str_mesh_cv__>
517  (basic_structure(cvr->structure()), k, force_simplexification);
518  str_mesh_cv__ &smc(*ppsmc);
519  psmc = ppsmc;
520 
521  smc.pm = std::make_unique<basic_mesh>();
522 
523  if (force_simplexification) {
524  // cout << "cvr = " << int(cvr->structure()->dim()) << " : "
525  // << cvr->structure()->nb_points() << endl;
526  const mesh_structure* splx_mesh
527  = basic_convex_ref(cvr)->simplexified_convex();
528  /* splx_mesh is correctly oriented */
529  for (size_type ic=0; ic < splx_mesh->nb_convex(); ++ic) {
530  std::vector<base_node> cvpts(splx_mesh->nb_points_of_convex(ic));
531  pgeometric_trans sgt
532  = simplex_geotrans(cvr->structure()->dim(), 1);
533  for (size_type j=0; j < cvpts.size(); ++j) {
534  size_type ip = splx_mesh->ind_points_of_convex(ic)[j];
535  cvpts[j] = basic_convex_ref(cvr)->points()[ip];
536  //cerr << "cvpts[" << j << "]=" << cvpts[j] << endl;
537  }
538  structured_mesh_for_convex_(splx_mesh->structure_of_convex(ic),
539  sgt, &cvpts, k, smc.pm.get());
540  }
541  } else {
542  structured_mesh_for_convex_(cvr->structure(), 0, 0, k, smc.pm.get());
543  }
544  smc.pfacem.resize(cvr->structure()->nb_faces());
545  for (dim_type f=0; f < cvr->structure()->nb_faces(); ++f) {
546  smc.pfacem[f] = std::make_unique<mesh_structure>();
547  structured_mesh_of_faces_(cvr, f, *(smc.pm), *(smc.pfacem[f]));
548  }
549 
550  smc.pmp = std::make_unique<mesh_precomposite>(*(smc.pm));
551  dal::add_stored_object(pk, psmc, basic_structure(cvr->structure()));
552  }
553  pm = psmc->pm.get();
554  pmp = psmc->pmp.get();
555  }
556 
557  const basic_mesh *
559  pbasic_mesh pm; pmesh_precomposite pmp;
560  structured_mesh_for_convex(cvr, k, pm, pmp, true);
561  return pm;
562  }
563 
564  const std::vector<std::unique_ptr<mesh_structure>>&
566  dal::pstatic_stored_object_key
567  pk = std::make_shared<str_mesh_key>(basic_structure(cvr->structure()),
568  k, true);
569  dal::pstatic_stored_object o = dal::search_stored_object(pk);
570  if (o) {
571  pstr_mesh_cv__ psmc = std::dynamic_pointer_cast<const str_mesh_cv__>(o);
572  return psmc->pfacem;
573  }
574  else GMM_ASSERT1(false,
575  "call refined_simplex_mesh_for_convex first (or fix me)");
576  }
577 
578 } /* end of namespace getfem. */
bgeot::basic_convex_ref
pconvex_ref basic_convex_ref(pconvex_ref cvr)
return the associated order 1 reference convex.
Definition: bgeot_convex_ref.h:138
bgeot::prism_P1_structure
pconvex_structure prism_P1_structure(dim_type nc)
Give a pointer on the structures of a prism of dimension d.
Definition: bgeot_convex_structure.h:206
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
dal_singleton.h
A simple singleton implementation.
bgeot::geotrans_interpolation_context
the geotrans_interpolation_context structure is passed as the argument of geometric transformation in...
Definition: bgeot_geometric_trans.h:411
gmm::vect_dist2
number_traits< typename linalg_traits< V1 >::value_type >::magnitude_type vect_dist2(const V1 &v1, const V2 &v2)
Euclidean distance between two vectors.
Definition: gmm_blas.h:597
bgeot::mesh_structure::ind_points_of_convex
const ind_set & ind_points_of_convex(size_type ic) const
Return a container to the list of points attached to convex ic.
Definition: bgeot_mesh_structure.h:107
bgeot::basic_structure
pconvex_structure basic_structure(pconvex_structure cv)
Original structure (if concerned)
Definition: bgeot_convex_structure.h:172
bgeot::short_type
gmm::uint16_type short_type
used as the common short type integer in the library
Definition: bgeot_config.h:72
bgeot::simplex_structure
pconvex_structure simplex_structure(dim_type nc)
Give a pointer on the structures of a simplex of dimension d.
Definition: bgeot_convex_structure.cc:148
dal::search_stored_object
pstatic_stored_object search_stored_object(pstatic_stored_object_key k)
Gives a pointer to an object from a key pointer.
Definition: dal_static_stored_objects.cc:177
bgeot::small_vector< scalar_type >
bgeot::mesh_structure::nb_convex
size_type nb_convex() const
The total number of convexes in the mesh.
Definition: bgeot_mesh_structure.h:96
bgeot::refined_simplex_mesh_for_convex_faces
const std::vector< std::unique_ptr< mesh_structure > > & refined_simplex_mesh_for_convex_faces(pconvex_ref cvr, short_type k)
simplexify the faces of a convex_ref
Definition: bgeot_poly_composite.cc:565
bgeot::refined_simplex_mesh_for_convex
const basic_mesh * refined_simplex_mesh_for_convex(pconvex_ref cvr, short_type k)
simplexify a convex_ref.
Definition: bgeot_poly_composite.cc:558
dal::add_stored_object
void add_stored_object(pstatic_stored_object_key k, pstatic_stored_object o, permanence perm)
Add an object with two optional dependencies.
Definition: dal_static_stored_objects.cc:284
bgeot::imbricated_box_less::operator()
int operator()(const base_node &x, const base_node &y) const
comparaison function
Definition: bgeot_poly_composite.cc:32
bgeot::structured_mesh_for_convex
void structured_mesh_for_convex(pconvex_ref cvr, short_type k, pbasic_mesh &pm, pmesh_precomposite &pmp, bool force_simplexification)
This function returns a mesh in pm which contains a refinement of the convex cvr if force_simplexific...
Definition: bgeot_poly_composite.cc:499
bgeot::mesh_structure
Mesh structure definition.
Definition: bgeot_mesh_structure.h:73
bgeot
Basic Geometric Tools.
Definition: bgeot_convex_ref.cc:27
bgeot::pconvex_structure
std::shared_ptr< const convex_structure > pconvex_structure
Pointer on a convex structure description.
Definition: bgeot_convex_structure.h:54
dal::static_stored_object
base class for static stored objects
Definition: dal_static_stored_objects.h:206
bgeot::parallelepiped_structure
pconvex_structure parallelepiped_structure(dim_type nc, dim_type k)
Give a pointer on the structures of a parallelepiped of dimension d.
Definition: bgeot_convex_structure.cc:397
bgeot::pgeometric_trans
std::shared_ptr< const bgeot::geometric_trans > pgeometric_trans
pointer type for a geometric transformation
Definition: bgeot_geometric_trans.h:186
bgeot_comma_init.h
convenient initialization of vectors via overload of "operator,".
bgeot_poly_composite.h
Handle composite polynomials.
bgeot::mesh_structure::nb_points_of_convex
short_type nb_points_of_convex(size_type ic) const
Return the number of points of convex ic.
Definition: bgeot_mesh_structure.h:115
gmm::mult_add
void mult_add(const L1 &l1, const L2 &l2, L3 &l3)
*‍/
Definition: gmm_blas.h:1781
bgeot::mesh_structure::structure_of_convex
pconvex_structure structure_of_convex(size_type ic) const
Return the pconvex_structure of the convex ic.
Definition: bgeot_mesh_structure.h:112