Hermes2D  3.0
neighbor_search.cpp
1 // This file is part of Hermes2D.
2 //
3 // Hermes2D is free software: you can redistribute it and/or modify
4 // it under the terms of the GNU General Public License as published by
5 // the Free Software Foundation, either version 2 of the License, or
6 // (at your option) any later version.
7 //
8 // Hermes2D is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.See the
11 // GNU General Public License for more details.
12 //
13 // You should have received a copy of the GNU General Public License
14 // along with Hermes2D. If not, see <http://www.gnu.org/licenses/>.
15 
16 #include "neighbor_search.h"
17 #include <algorithm>
18 
19 namespace Hermes
20 {
21  namespace Hermes2D
22  {
23  template<typename Scalar>
25  supported_shapes(nullptr),
26  mesh(mesh),
27  central_el(el),
28  neighb_el(nullptr),
29  quad(&g_quad_2d_std),
30  central_transformations_size(0),
31  neighbor_transformations_size(0),
32  central_transformations_alloc_size(H2D_INITIAL_NEIGHBOR_NUMBER_GUESS),
33  neighbor_transformations_alloc_size(H2D_INITIAL_NEIGHBOR_NUMBER_GUESS)
34  {
35  if (central_el == nullptr || central_el->active != 1)
36  throw Exceptions::Exception("You must pass an active element to the NeighborSearch constructor.");
37  neighbors.reserve(2);
38  neighbor_edges.reserve(2);
39 
40  central_transformations = calloc_with_check<NeighborSearch<Scalar>, Transformations*>(H2D_INITIAL_NEIGHBOR_NUMBER_GUESS, this);
41  neighbor_transformations = calloc_with_check<NeighborSearch<Scalar>, Transformations*>(H2D_INITIAL_NEIGHBOR_NUMBER_GUESS, this);
42 
43  ignore_errors = false;
44  n_neighbors = 0;
45  neighborhood_type = H2D_DG_NOT_INITIALIZED;
47  }
48 
49  template<typename Scalar>
51  supported_shapes(nullptr),
52  mesh(ns.mesh),
53  central_el(ns.central_el),
54  neighb_el(nullptr),
55  neighbor_edge(ns.neighbor_edge),
56  active_segment(ns.active_segment)
57  {
58  neighbors.reserve(2);
59  neighbor_edges.reserve(2);
60 
61  central_transformations = calloc_with_check<NeighborSearch<Scalar>, Transformations*>(ns.central_transformations_alloc_size, this);
62  neighbor_transformations = calloc_with_check<NeighborSearch<Scalar>, Transformations*>(ns.neighbor_transformations_alloc_size, this);
63 
64  for (unsigned int j = 0; j < ns.central_transformations_alloc_size; j++)
65  {
66  if (ns.central_transformations[j])
67  this->central_transformations[j] = new Transformations(ns.central_transformations[j]);
68  }
69 
70  this->central_transformations_size = ns.central_transformations_size;
71  this->central_transformations_alloc_size = ns.central_transformations_alloc_size;
72 
73  for (unsigned int j = 0; j < ns.neighbor_transformations_alloc_size; j++)
74  {
75  if (ns.neighbor_transformations[j])
76  this->neighbor_transformations[j] = new Transformations(ns.neighbor_transformations[j]);
77  }
78 
79  this->neighbor_transformations_size = ns.neighbor_transformations_size;
80  this->neighbor_transformations_alloc_size = ns.neighbor_transformations_alloc_size;
81 
82  if (central_el == nullptr || central_el->active != 1)
83  throw Exceptions::Exception("You must pass an active element to the NeighborSearch constructor.");
84 
85  for (unsigned int i = 0; i < ns.neighbors.size(); i++)
86  this->neighbors.push_back(ns.neighbors[i]);
87  for (unsigned int i = 0; i < ns.neighbor_edges.size(); i++)
88  this->neighbor_edges.push_back(ns.neighbor_edges[i]);
89 
91  n_neighbors = ns.n_neighbors;
92  neighborhood_type = ns.neighborhood_type;
94  quad = (&g_quad_2d_std);
95  active_edge = ns.active_edge;
96  }
97 
98  template<typename Scalar>
100  {
101  this->free();
102  }
103 
104  template<typename Scalar>
106  {
107  neighbor_edges.clear();
108  neighbors.clear();
109  clear_supported_shapes();
110 
111  for (unsigned int i = 0; i < central_transformations_alloc_size; i++)
112  if (this->central_transformations[i])
113  delete this->central_transformations[i];
114  free_with_check(this->central_transformations);
115  for (unsigned int i = 0; i < neighbor_transformations_alloc_size; i++)
116  if (this->neighbor_transformations[i])
117  delete this->neighbor_transformations[i];
118  free_with_check(this->neighbor_transformations);
119  }
120 
121  template<typename Scalar>
123  {
124  return n_neighbors;
125  }
126 
127  template<typename Scalar>
128  const std::vector<Element*>* NeighborSearch<Scalar>::get_neighbors() const
129  {
130  return &neighbors;
131  }
132 
133  template<typename Scalar>
135  {
136  if (supported_shapes != nullptr)
137  delete supported_shapes;
138  supported_shapes = nullptr;
139  }
140 
141  template<typename Scalar>
143  {
144  this->ignore_errors = value;
145  }
146 
147  template<typename Scalar>
149  {
150  // Reset transformations.
151  for (unsigned int i = 0; i < this->central_transformations_alloc_size; i++)
152  {
153  if (this->central_transformations[i])
154  this->central_transformations[i]->reset();
155  }
156  for (unsigned int i = 0; i < this->neighbor_transformations_alloc_size; i++)
157  {
158  if (this->neighbor_transformations[i])
159  this->neighbor_transformations[i]->reset();
160  }
161 
162  // Reset information about the neighborhood's active state.
163  active_segment = 0;
164  active_edge = 0;
165  neighb_el = nullptr;
166  neighbor_edge.local_num_of_edge = 0;
167 
168  // Clear vectors with neighbor elements and their edge info for the active edge.
169  neighbor_edges.clear();
170  neighbors.clear();
171  n_neighbors = 0;
172 
173  neighborhood_type = H2D_DG_NOT_INITIALIZED;
174  }
175 
176  template<typename Scalar>
178  {
179  reset_neighb_info();
180  active_edge = edge;
181 
182  if (!this->ignore_errors && central_el->en[active_edge]->bnd == 1)
183  throw Hermes::Exceptions::Exception("The given edge isn't inner");
184  else
185  {
186  neighb_el = central_el->get_neighbor(active_edge);
187 
188  // First case : The neighboring element is of the same size as the central one.
189  if (neighb_el != nullptr)
190  {
191  //std::cout << "\t active neighbor el: " << neighb_el->id << std::endl;
192 
193  // Get local number of the edge used by the neighbor.
194  for (unsigned int j = 0; j < neighb_el->get_nvert(); j++)
195  {
196  if (central_el->en[active_edge] == neighb_el->en[j])
197  {
198  neighbor_edge.local_num_of_edge = j;
199  break;
200  }
201  }
202 
203  NeighborEdgeInfo local_edge_info;
204  local_edge_info.local_num_of_edge = neighbor_edge.local_num_of_edge;
205 
206  // Query the orientation of the neighbor edge relative to the central el.
207  int p1 = central_el->vn[active_edge]->id;
208  int p2 = central_el->vn[(active_edge + 1) % central_el->get_nvert()]->id;
209  local_edge_info.orientation = neighbor_edge_orientation(p1, p2, false);
210 
211  neighbor_edges.push_back(local_edge_info);
212 
213  // There is only one neighbor in this case.
214  n_neighbors = 1;
215  neighbors.push_back(neighb_el);
216 
217  // No need for transformation, since the neighboring element is of the same size.
218  neighborhood_type = H2D_DG_NO_TRANSF;
219  }
220  else
221  {
222  // Peek the vertex in the middle of the active edge (if there is none, vertex will be nullptr).
223  Node* vertex = mesh->peek_vertex_node(central_el->en[active_edge]->p1, central_el->en[active_edge]->p2);
224 
225  // Endpoints of the active edge.
226  int orig_vertex_id[2];
227  orig_vertex_id[0] = central_el->vn[active_edge]->id;
228  orig_vertex_id[1] = central_el->vn[(active_edge + 1) % central_el->get_nvert()]->id;
229 
230  if (vertex == nullptr)
231  {
232  neighborhood_type = H2D_DG_GO_UP;
233 
234  Element* parent = central_el->parent;
235 
236  // Array of middle-point vertices of the intermediate parent edges that we climb up to the correct parent element.
237  Node** par_mid_vertices = malloc_with_check<NeighborSearch<Scalar>, Node*>(Transformations::max_level, this);
238  // Number of visited intermediate parents.
239  int n_parents = 0;
240 
241  for (unsigned int j = 0; j < (unsigned)Transformations::max_level; j++)
242  par_mid_vertices[j] = nullptr;
243 
244  find_act_elem_up(parent, orig_vertex_id, par_mid_vertices, n_parents);
245 
246  free_with_check(par_mid_vertices);
247  }
248  else
249  {
250  neighborhood_type = H2D_DG_GO_DOWN;
251 
252  // array of virtual sons of the central el. visited on the way down to the neighbor
253  int sons[Transformations::max_level];
254  // number of used transformations
255  int n_sons = 0;
256 
257  // Start the search by going down to the first son.
258  find_act_elem_down(vertex, orig_vertex_id, sons, n_sons + 1);
259 
260  //debug_log("number of neighbors on the way down: %d ", n_neighbors);
261  }
262  }
263  }
264  }
265 
266  template<typename Scalar>
268  {
269  std::vector<unsigned int> transformations = get_transforms(original_central_el_transform);
270  // Inter-element edge.
271  if (is_inter_edge(edge, transformations))
272  {
273  set_active_edge(edge);
274  update_according_to_sub_idx(transformations);
275  return true;
276  }
277  // Intra-element edge.
278  else
279  {
280  neighb_el = central_el;
281 
282  assert(!neighbor_transformations[0]);
283  this->add_neighbor_transformations(new Transformations(transformations), 0);
284 
285  neighbor_edge.local_num_of_edge = active_edge = edge;
286  NeighborEdgeInfo local_edge_info;
287  local_edge_info.local_num_of_edge = neighbor_edge.local_num_of_edge;
288  // The "opposite" view of the same edge has the same orientation.
289  local_edge_info.orientation = 0;
290  neighbor_edges.push_back(local_edge_info);
291 
292  n_neighbors = 1;
293  neighbors.push_back(neighb_el);
294  neighborhood_type = H2D_DG_NO_TRANSF;
295  return false;
296  }
297  }
298 
299  template<typename Scalar>
300  std::vector<unsigned int> NeighborSearch<Scalar>::get_transforms(uint64_t sub_idx) const
301  {
302  std::vector<unsigned int> transformations_backwards;
303  while (sub_idx > 0)
304  {
305  transformations_backwards.push_back((sub_idx - 1) & 7);
306  sub_idx = (sub_idx - 1) >> 3;
307  }
308  std::vector<unsigned int> transformations;
309  for (unsigned int i = 0; i < transformations_backwards.size(); i++)
310  transformations.push_back(transformations_backwards[transformations_backwards.size() - 1 - i]);
311 
312  return transformations;
313  }
314 
315  template<typename Scalar>
316  bool NeighborSearch<Scalar>::is_inter_edge(const int& edge, const std::vector<unsigned int>& transformations) const
317  {
318  // No subelements => of course this edge is an inter-element one.
319  if (transformations.size() == 0)
320  return true;
321 
322  // Triangles.
323  for (unsigned int i = 0; i < transformations.size(); i++)
324  if (central_el->get_mode() == HERMES_MODE_TRIANGLE)
325  {
326  if ((edge == 0 && (transformations[i] == 2 || transformations[i] == 3)) ||
327  (edge == 1 && (transformations[i] == 0 || transformations[i] == 3)) ||
328  (edge == 2 && (transformations[i] == 1 || transformations[i] == 3)))
329  return false;
330  }
331  // Quads.
332  else
333  {
334  if ((edge == 0 && (transformations[i] == 2 || transformations[i] == 3 || transformations[i] == 5)) ||
335  (edge == 1 && (transformations[i] == 0 || transformations[i] == 3 || transformations[i] == 6)) ||
336  (edge == 2 && (transformations[i] == 0 || transformations[i] == 1 || transformations[i] == 4)) ||
337  (edge == 3 && (transformations[i] == 1 || transformations[i] == 2 || transformations[i] == 7)))
338  return false;
339  }
340  return true;
341  }
342 
343  template<typename Scalar>
344  void NeighborSearch<Scalar>::add_central_transformations(Transformations* to_add, int position, bool replace)
345  {
346  if (position >= this->central_transformations_alloc_size)
347  {
348  int new_alloc_size = std::min(position + 1, (int)(position * 1.5));
349  this->central_transformations = (Transformations**)realloc(this->central_transformations, new_alloc_size * sizeof(Transformations*));
350  memset(this->central_transformations + this->central_transformations_alloc_size, 0, sizeof(Transformations*)* (new_alloc_size - this->central_transformations_alloc_size));
351  this->central_transformations_alloc_size = new_alloc_size;
352  }
353 
354  if (this->central_transformations[position])
355  {
356  if (replace)
357  {
358  delete this->central_transformations[position];
359  this->central_transformations[position] = to_add;
360  }
361  }
362  else
363  {
364  this->central_transformations[position] = to_add;
365  this->central_transformations_size++;
366  }
367  }
368 
369  template<typename Scalar>
370  void NeighborSearch<Scalar>::add_neighbor_transformations(Transformations* to_add, int position, bool replace)
371  {
372  if (position >= this->neighbor_transformations_alloc_size)
373  {
374  int new_alloc_size = std::min(position + 1, (int)(position * 1.5));
375  this->neighbor_transformations = (Transformations**)realloc(this->neighbor_transformations, new_alloc_size * sizeof(Transformations*));
376  memset(this->neighbor_transformations + this->neighbor_transformations_alloc_size, 0, sizeof(Transformations*)* (new_alloc_size - this->neighbor_transformations_alloc_size));
377  this->neighbor_transformations_alloc_size = new_alloc_size;
378  }
379 
380  if (this->neighbor_transformations[position])
381  {
382  if (replace)
383  {
384  delete this->neighbor_transformations[position];
385  this->neighbor_transformations[position] = to_add;
386  }
387  }
388  else
389  {
390  this->neighbor_transformations[position] = to_add;
391  this->neighbor_transformations_size++;
392  }
393  }
394 
395  template<typename Scalar>
396  void NeighborSearch<Scalar>::update_according_to_sub_idx(const std::vector<unsigned int>& transformations)
397  {
398  if (neighborhood_type == H2D_DG_NO_TRANSF && transformations.size() == 0)
399  return;
400 
401  if (neighborhood_type == H2D_DG_NO_TRANSF || neighborhood_type == H2D_DG_GO_UP)
402  {
403  if (!neighbor_transformations[0])
404  {
405  this->add_neighbor_transformations(new Transformations(transformations), 0);
406  }
407 
408  Transformations *tr = neighbor_transformations[0];
409 
410  for (unsigned int i = 0; i < transformations.size(); i++)
411  // Triangles.
412  if (central_el->get_mode() == HERMES_MODE_TRIANGLE)
413  if ((active_edge == 0 && transformations[i] == 0) ||
414  (active_edge == 1 && transformations[i] == 1) ||
415  (active_edge == 2 && transformations[i] == 2))
416  tr->transf[tr->num_levels++] = (!neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 3);
417  else
418  tr->transf[tr->num_levels++] = (neighbor_edges[0].orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 3);
419  // Quads.
420  else
421  if ((active_edge == 0 && (transformations[i] == 0 || transformations[i] == 6)) ||
422  (active_edge == 1 && (transformations[i] == 1 || transformations[i] == 4)) ||
423  (active_edge == 2 && (transformations[i] == 2 || transformations[i] == 7)) ||
424  (active_edge == 3 && (transformations[i] == 3 || transformations[i] == 5)))
425  tr->transf[tr->num_levels++] = (!neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 4);
426  else if ((active_edge == 0 && (transformations[i] == 1 || transformations[i] == 7)) ||
427  (active_edge == 1 && (transformations[i] == 2 || transformations[i] == 5)) ||
428  (active_edge == 2 && (transformations[i] == 3 || transformations[i] == 6)) ||
429  (active_edge == 3 && (transformations[i] == 0 || transformations[i] == 4)))
430  tr->transf[tr->num_levels++] = (neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 4);
431  }
432  else handle_sub_idx_way_down(transformations);
433  }
434 
435  template<typename Scalar>
436  void NeighborSearch<Scalar>::handle_sub_idx_way_down(const std::vector<unsigned int>& transformations)
437  {
438  std::vector<unsigned int> neighbors_to_be_deleted;
439  std::vector<unsigned int> neighbors_not_to_be_deleted;
440 
441  std::vector<unsigned int> updated_transformations;
442  for (unsigned short i = 0; i < transformations.size(); i++)
443  {
444  if (!((active_edge == 0 && transformations[i] == 4) || (active_edge == 1 && transformations[i] == 7) || (active_edge == 2 && transformations[i] == 5) || (active_edge == 3 && transformations[i] == 6)))
445  {
446  if (active_edge == 0 && transformations[i] == 6)
447  updated_transformations.push_back(0);
448  else if (active_edge == 0 && transformations[i] == 7)
449  updated_transformations.push_back(1);
450  else if (active_edge == 1 && transformations[i] == 4)
451  updated_transformations.push_back(1);
452  else if (active_edge == 1 && transformations[i] == 5)
453  updated_transformations.push_back(2);
454  else if (active_edge == 2 && transformations[i] == 6)
455  updated_transformations.push_back(3);
456  else if (active_edge == 2 && transformations[i] == 7)
457  updated_transformations.push_back(2);
458  else if (active_edge == 3 && transformations[i] == 4)
459  updated_transformations.push_back(0);
460  else if (active_edge == 3 && transformations[i] == 5)
461  updated_transformations.push_back(3);
462  else
463  updated_transformations.push_back(transformations[i]);
464  }
465  }
466 
467  // We basically identify the neighbors that are not compliant with the current sub-element mapping on the central element.
468  for (unsigned int neighbor_i = 0; neighbor_i < n_neighbors; neighbor_i++)
469  {
470  bool deleted = false;
471 
472  Transformations* current_transforms = central_transformations[neighbor_i];
473 
474  for (unsigned int level = 0; level < std::min((unsigned int)updated_transformations.size(), current_transforms->num_levels); level++)
475  {
476  // If the found neighbor is not a neighbor of this subelement.
477  if (!compatible_transformations(current_transforms->transf[level], updated_transformations[level], active_edge))
478  {
479  deleted = true;
480  break;
481  }
482  }
483 
484  /*
485  // If there were more sub-element transformation from the assembling than from the neighbor search.
486  if(!deleted)
487  {
488  if((unsigned int)updated_transformations.size() > current_transforms->num_levels)
489  {
490  for(unsigned int level = current_transforms->num_levels; level < (unsigned int)updated_transformations.size(); level++)
491  {
492  // If the found neighbor is not a neighbor of this subelement.
493  if(!compatible_transformations(current_transforms->transf[current_transforms->num_levels - 1], updated_transformations[level], active_edge))
494  {
495  deleted = true;
496  break;
497  }
498  }
499  }
500  }
501  */
502  if (deleted)
503  neighbors_to_be_deleted.push_back(neighbor_i);
504  else
505  neighbors_not_to_be_deleted.push_back(neighbor_i);
506  }
507 
508  // Now we truly delete (in the reverse order) the neighbors.
509  if (neighbors_to_be_deleted.size() > 0)
510  for (unsigned int neighbors_to_be_deleted_i = neighbors_to_be_deleted.size(); neighbors_to_be_deleted_i >= 1; neighbors_to_be_deleted_i--)
511  delete_neighbor(neighbors_to_be_deleted[neighbors_to_be_deleted_i - 1]);
512  }
513 
514  template<typename Scalar>
515  bool NeighborSearch<Scalar>::compatible_transformations(unsigned int a, unsigned int b, int edge) const
516  {
517  if (a == b)
518  return true;
519  if (edge == 0)
520  {
521  if ((a == 0 && (b == 6 || b == 4)) ||
522  (a == 1 && (b == 7 || b == 4)))
523  return true;
524  else
525  return false;
526  }
527  if (edge == 1)
528  {
529  if ((a == 1 && (b == 4 || b == 7)) ||
530  (a == 2 && (b == 5 || b == 7)))
531  return true;
532  else
533  return false;
534  }
535  if (edge == 2)
536  {
537  if ((a == 2 && (b == 7 || b == 5)) ||
538  (a == 3 && (b == 6 || b == 5)))
539  return true;
540  else
541  return false;
542  }
543  if (edge == 3)
544  {
545  if ((a == 3 && (b == 5 || b == 6)) ||
546  (a == 0 && (b == 4 || b == 6)))
547  return true;
548  else
549  return false;
550  }
551  return false;
552  }
553 
554  template<typename Scalar>
556  {
557  if (neighborhood_type != H2D_DG_GO_DOWN)
558  return;
559  // Obtain the transformations sequence.
560  std::vector<unsigned int> transformations = get_transforms(original_central_el_transform);
561 
562  std::vector<unsigned int> updated_transformations;
563  for (unsigned short i = 0; i < transformations.size(); i++)
564  {
565  if (!((active_edge == 0 && transformations[i] == 4) || (active_edge == 1 && transformations[i] == 7) || (active_edge == 2 && transformations[i] == 5) || (active_edge == 3 && transformations[i] == 6)))
566  {
567  if (active_edge == 0 && transformations[i] == 6)
568  updated_transformations.push_back(0);
569  else if (active_edge == 0 && transformations[i] == 7)
570  updated_transformations.push_back(1);
571  else if (active_edge == 1 && transformations[i] == 4)
572  updated_transformations.push_back(1);
573  else if (active_edge == 1 && transformations[i] == 5)
574  updated_transformations.push_back(2);
575  else if (active_edge == 2 && transformations[i] == 6)
576  updated_transformations.push_back(3);
577  else if (active_edge == 2 && transformations[i] == 7)
578  updated_transformations.push_back(2);
579  else if (active_edge == 3 && transformations[i] == 4)
580  updated_transformations.push_back(0);
581  else if (active_edge == 3 && transformations[i] == 5)
582  updated_transformations.push_back(3);
583  else
584  updated_transformations.push_back(transformations[i]);
585  }
586  }
587 
588  // Test for active element.
589  if (updated_transformations.size() == 0)
590  return;
591 
592  for (unsigned int i = 0; i < n_neighbors; i++)
593  {
594  // Find the index where the additional subelement mapping (on top of the initial one from assembling) starts.
595  unsigned int j = 0;
596  // Note that we do not have to test if central_transformations is empty or how long it is, because it has to be
597  // longer than transformations (and that is tested).
598  // Also the function compatible_transformations() does not have to be used, as now the array central_transformations
599  // has been adjusted so that it contains the array transformations.
600  while (central_transformations[i]->transf[j] == updated_transformations[j])
601  if (++j > updated_transformations.size() - 1)
602  break;
603  if (j > central_transformations[i]->num_levels)
604  j = central_transformations[i]->num_levels;
605 
606  for (unsigned int level = central_transformations[i]->num_levels; level < updated_transformations.size(); level++)
607  {
608  if ((i >= neighbor_transformations_alloc_size) || !neighbor_transformations[i])
609  this->add_neighbor_transformations(new Transformations, i);
610 
611  Transformations* neighbor_transforms = neighbor_transformations[i];
612 
613  // Triangles.
614  if (central_el->get_mode() == HERMES_MODE_TRIANGLE)
615  if ((active_edge == 0 && updated_transformations[level] == 0) ||
616  (active_edge == 1 && updated_transformations[level] == 1) ||
617  (active_edge == 2 && updated_transformations[level] == 2))
618  neighbor_transforms->transf[neighbor_transforms->num_levels++] = (!neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 3);
619  else
620  neighbor_transforms->transf[neighbor_transforms->num_levels++] = (neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 3);
621  // Quads.
622  else
623  if ((active_edge == 0 && (updated_transformations[level] == 0 || updated_transformations[level] == 6)) ||
624  (active_edge == 1 && (updated_transformations[level] == 1 || updated_transformations[level] == 4)) ||
625  (active_edge == 2 && (updated_transformations[level] == 2 || updated_transformations[level] == 7)) ||
626  (active_edge == 3 && (updated_transformations[level] == 3 || updated_transformations[level] == 5)))
627  neighbor_transforms->transf[neighbor_transforms->num_levels++] = (!neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 4);
628  else if ((active_edge == 0 && (updated_transformations[level] == 1 || updated_transformations[level] == 7)) ||
629  (active_edge == 1 && (updated_transformations[level] == 2 || updated_transformations[level] == 5)) ||
630  (active_edge == 2 && (updated_transformations[level] == 3 || updated_transformations[level] == 6)) ||
631  (active_edge == 3 && (updated_transformations[level] == 0 || updated_transformations[level] == 4)))
632  neighbor_transforms->transf[neighbor_transforms->num_levels++] = (neighbor_edge.orientation ? neighbor_edge.local_num_of_edge : (neighbor_edge.local_num_of_edge + 1) % 4);
633  }
634 
635  central_transformations[i]->strip_initial_transformations(j);
636  }
637  }
638 
639  template<typename Scalar>
641  {
642  // Create a new_ clear array for the remaining transformations.
643  unsigned int shifted_trfs[max_level];
644  memset(shifted_trfs, 0, max_level * sizeof(unsigned int));
645  // Move the old one to the new_ one.
646  for (unsigned int k = number_of_stripped; k < num_levels; k++)
647  shifted_trfs[k - number_of_stripped] = transf[k];
648  // Point to the new_ one.
649  memcpy(transf, shifted_trfs, max_level*sizeof(unsigned int));
650  // We also have to store the information about length of the transformation array for this neighbor.
651  num_levels -= number_of_stripped;
652  }
653 
654  template<typename Scalar>
655  void NeighborSearch<Scalar>::delete_neighbor(unsigned int position)
656  {
657  for (unsigned int i = position; i < n_neighbors - 1; i++)
658  if (i + 1 < this->central_transformations_alloc_size)
659  central_transformations[i]->copy_from(central_transformations[i + 1]);
660 
661  if ((n_neighbors - 1 < this->central_transformations_alloc_size) && central_transformations[n_neighbors - 1])
662  central_transformations[n_neighbors - 1]->reset();
663 
664  for (unsigned int i = position; i < n_neighbors - 1; i++)
665  {
666  if ((i + 1 < this->neighbor_transformations_alloc_size) && neighbor_transformations[i + 1])
667  {
668  if (!neighbor_transformations[i])
669  this->add_neighbor_transformations(new Transformations, i);
670 
671  neighbor_transformations[i]->copy_from(neighbor_transformations[i + 1]);
672  }
673  }
674  if ((n_neighbors - 1 < this->neighbor_transformations_alloc_size) && neighbor_transformations[n_neighbors - 1])
675  neighbor_transformations[n_neighbors - 1]->reset();
676 
677  neighbor_edges.erase(neighbor_edges.begin() + position);
678  neighbors.erase(neighbors.begin() + position);
679  n_neighbors--;
680  }
681 
682  template<typename Scalar>
683  void NeighborSearch<Scalar>::find_act_elem_up(Element* elem, int* orig_vertex_id, Node** par_mid_vertices, int n_parents)
684  {
685  Node* edge = nullptr;
686  Node* vertex = nullptr;
687 
688  assert(n_parents <= (int)Transformations::max_level);
689 
690  // IDs of vertices bounding the current intermediate parent edge.
691  int p1 = elem->vn[active_edge]->id;
692  int p2 = elem->vn[(active_edge + 1) % elem->get_nvert()]->id;
693 
694  int id_of_par_orient_1 = p1;
695  int id_of_par_orient_2 = p2;
696 
697  // Find if p1 and p2 bound a used edge (used by the neighbor element).
698  edge = mesh->peek_edge_node(p1, p2);
699 
700  // Add the vertex in the middle of the parent edge to the array of intermediate parent vertices. This is for
701  // consequent transformation of functions on neighbor element.
702  vertex = mesh->peek_vertex_node(p1, p2);
703  if (vertex != nullptr)
704  {
705  if (n_parents == 0)
706  par_mid_vertices[n_parents++] = vertex;
707  else
708  if (n_parents == Transformations::max_level - 1)
709  throw Hermes::Exceptions::Exception("Maximum number of intermediate parents exceeded in NeighborSearch<Scalar>::finding_act_elem_up");
710  else
711  if (par_mid_vertices[n_parents - 1]->id != vertex->id)
712  par_mid_vertices[n_parents++] = vertex;
713  }
714 
715  if ((edge == nullptr) || (central_el->en[active_edge]->id == edge->id))
716  {
717  // We have not yet found the parent of the central element completely adjacent to the neighbor.
718  find_act_elem_up(elem->parent, orig_vertex_id, par_mid_vertices, n_parents);
719  }
720  else
721  {
722  for (int i = 0; i < 2; i++)
723  {
724  // Get a pointer to the active neighbor element.
725  if ((edge->elem[i] != nullptr) && (edge->elem[i]->active == 1))
726  {
727  //debug_log("way up neighbor: %d", neighb_el->id);
728  neighb_el = edge->elem[i];
729 
730  // Get local number of the edge used by the neighbor.
731  neighbor_edge.local_num_of_edge = -1;
732  for (unsigned int j = 0; j < neighb_el->get_nvert(); j++)
733  if (neighb_el->en[j] == edge)
734  {
735  neighbor_edge.local_num_of_edge = j;
736  break;
737  }
738  if (neighbor_edge.local_num_of_edge == -1) throw Hermes::Exceptions::Exception("Neighbor edge wasn't found");
739 
740  // Add to the array of neighbor_transformations one that transforms central el. to its parent completely
741  // adjacent to the single big neighbor.
742  assert(n_neighbors == 0);
743 
744  if (n_neighbors >= this->neighbor_transformations_alloc_size || !neighbor_transformations[n_neighbors])
745  this->add_neighbor_transformations(new Transformations, n_neighbors);
746 
747  Transformations *neighbor_transforms = neighbor_transformations[n_neighbors];
748 
749  neighbor_transforms->num_levels = n_parents;
750 
751  // Go back through the intermediate inactive parents down to the central element and stack corresponding
752  // neighbor_transformations into the array 'neighbor_transformations'.
753  for (int j = n_parents - 1; j > 0; j--)
754  {
755  Node* n = mesh->peek_vertex_node(par_mid_vertices[j]->id, p1);
756  if (n == nullptr)
757  {
758  neighbor_transforms->transf[n_parents - j - 1] = neighbor_edge.local_num_of_edge;
759  p1 = par_mid_vertices[j]->id;
760  }
761  else
762  {
763  if (n->id == par_mid_vertices[j - 1]->id)
764  {
765  neighbor_transforms->transf[n_parents - j - 1] = (neighbor_edge.local_num_of_edge + 1) % neighb_el->get_nvert();
766  p2 = par_mid_vertices[j]->id;
767  }
768  else
769  {
770  neighbor_transforms->transf[n_parents - j - 1] = neighbor_edge.local_num_of_edge;
771  p1 = par_mid_vertices[j]->id;
772  }
773  }
774  }
775 
776  // Final transformation to the central element itself.
777  if (orig_vertex_id[0] == par_mid_vertices[0]->id)
778  neighbor_transforms->transf[n_parents - 1] = neighbor_edge.local_num_of_edge;
779  else
780  neighbor_transforms->transf[n_parents - 1] = (neighbor_edge.local_num_of_edge + 1) % neighb_el->get_nvert();
781 
782  NeighborEdgeInfo local_edge_info;
783  local_edge_info.local_num_of_edge = neighbor_edge.local_num_of_edge;
784 
785  // Query the orientation of the neighbor edge relative to the central el.
786  local_edge_info.orientation = neighbor_edge_orientation(id_of_par_orient_1, id_of_par_orient_2, false);
787 
788  neighbor_edges.push_back(local_edge_info);
789 
790  // There is only one neighbor, ...
791  n_neighbors = 1;
792 
793  // ...add it to the vector of neighbors.
794  neighbors.push_back(neighb_el);
795  }
796  }
797  }
798  }
799 
800  template<typename Scalar>
801  void NeighborSearch<Scalar>::find_act_elem_down(Node* vertex, int* bounding_verts_id, int* sons, unsigned int n_sons)
802  {
803  // ID of vertex in between vertices from par_vertex_id.
804  int mid_vert = vertex->id;
805  int bnd_verts[2];
806  bnd_verts[0] = bounding_verts_id[0];
807  bnd_verts[1] = bounding_verts_id[1];
808 
809  assert(n_sons < (unsigned)Transformations::max_level);
810 
811  for (int i = 0; i < 2; i++)
812  {
813  sons[n_sons - 1] = (active_edge + i) % central_el->get_nvert();
814 
815  // Try to get a pointer to the edge between the middle vertex and one of the vertices bounding the previously
816  // tested segment.
817  Node* edge = mesh->peek_edge_node(mid_vert, bnd_verts[i]);
818 
819  if (edge == nullptr) // The edge is not used, i.e. there is no active element on either side.
820  {
821  // Get the middle vertex of this edge and try again on the segments into which this vertex splits the edge.
822  Node * n = mesh->peek_vertex_node(mid_vert, bnd_verts[i]);
823  if (n == nullptr)
824  throw Hermes::Exceptions::Exception("wasn't able to find middle vertex");
825  else
826  {
827  // Make sure the next visited segment has the same orientation as the original central element's active edge.
828  if (i == 0)
829  bounding_verts_id[1] = mid_vert;
830  else
831  bounding_verts_id[0] = mid_vert;
832 
833  find_act_elem_down(n, bounding_verts_id, sons, n_sons + 1);
834 
835  bounding_verts_id[0] = bnd_verts[0];
836  bounding_verts_id[1] = bnd_verts[1];
837  }
838  }
839  else // We have found a used edge, the active neighbor we are looking for is on one of its sides.
840  {
841  for (int j = 0; j < 2; j++)
842  {
843  if ((edge->elem[j] != nullptr) && (edge->elem[j]->active == 1))
844  {
845  //debug_log("way down neighbor: %d", edge->elem[j]->id);
846  neighb_el = mesh->get_element(edge->elem[j]->id);
847 
848  // Get local number of the edge used by the neighbor.
849  neighbor_edge.local_num_of_edge = -1;
850  for (unsigned int k = 0; k < neighb_el->get_nvert(); k++)
851  if (neighb_el->en[k] == edge)
852  {
853  neighbor_edge.local_num_of_edge = k;
854  break;
855  }
856 
857  if (neighbor_edge.local_num_of_edge == -1) throw Hermes::Exceptions::Exception("Neighbor edge wasn't found");
858 
859  if ((n_neighbors >= this->central_transformations_alloc_size) || !central_transformations[n_neighbors])
860  this->add_central_transformations(new Transformations, n_neighbors);
861 
862  Transformations *tr = central_transformations[n_neighbors];
863 
864  // Construct the transformation path to the current neighbor.
865  for (unsigned int k = 0; k < n_sons; k++)
866  tr->transf[k] = sons[k];
867  tr->num_levels = n_sons;
868 
869  NeighborEdgeInfo local_edge_info;
870  local_edge_info.local_num_of_edge = neighbor_edge.local_num_of_edge;
871  // Query the orientation of the neighbor edge relative to the central el.
872  local_edge_info.orientation = neighbor_edge_orientation(bnd_verts[0], bnd_verts[1], i);
873 
874  neighbor_edges.push_back(local_edge_info);
875 
876  // Append the new_ neighbor.
877  n_neighbors++;
878  neighbors.push_back(neighb_el);
879  }
880  }
881  }
882  }
883  }
884 
885  template<typename Scalar>
886  bool NeighborSearch<Scalar>::neighbor_edge_orientation(int bounding_vert1, int bounding_vert2, int segment) const
887  {
888  if (segment == 0)
889  {
890  // neighbor edge goes from parent1 to middle vertex
891  if (neighb_el->vn[neighbor_edge.local_num_of_edge]->id != bounding_vert1)
892  // orientation reversed
893  return true;
894  }
895  else
896  {
897  // neighbor edge goes from middle vertex to parent2
898  if (neighb_el->vn[neighbor_edge.local_num_of_edge]->id == bounding_vert2)
899  // orientation reversed
900  return true;
901  }
902  return false;
903  }
904 
905  template<typename Scalar>
906  NeighborSearch<Scalar>::ExtendedShapeset::ExtendedShapeset(const ExtendedShapeset & other)
907  {
908  this->central_al = new AsmList<Scalar>(*other.central_al);
909  this->cnt = other.cnt;
910  this->dof = other.dof;
911  this->neighbor_al = new AsmList<Scalar>(*other.neighbor_al);
912  this->combine_assembly_lists();
913  }
914 
915  template<typename Scalar>
917  {
918  ExtendedShapeset* new_supp_shapes = new ExtendedShapeset(this, al, space);
919 
920  return new_supp_shapes;
921  }
922 
923  template<typename Scalar>
925  {
926  if (supported_shapes != nullptr)
927  delete supported_shapes;
928 
929  supported_shapes = new ExtendedShapeset(this, al, space);
930 
931  return new ExtendedShapeset(*supported_shapes);
932  }
933 
934  template<typename Scalar>
936  {
937  neighb_quad_order = quad->get_edge_points(neighbor_edge.local_num_of_edge, order, neighbors[active_segment]->get_mode());
938  central_quad_order = quad->get_edge_points(active_edge, order, central_el->get_mode());
939  }
940 
941  template<typename Scalar>
942  int NeighborSearch<Scalar>::get_quad_eo(bool on_neighbor) const
943  {
944  if (on_neighbor)
945  return neighb_quad_order;
946  else
947  return central_quad_order;
948  }
949 
950  template<typename Scalar>
952  {
953  return this->active_segment;
954  }
955 
956  template<typename Scalar>
958  {
959  if (index >= n_neighbors)
960  throw Hermes::Exceptions::Exception("NeighborSearch<Scalar>::set_active_segment() called with an incorrect index.");
961 
962  this->active_segment = index;
963  this->neighb_el = this->neighbors[index];
964  this->neighbor_edge = this->neighbor_edges[index];
965  }
966 
967  template<typename Scalar>
969  {
970  return this->neighb_el;
971  }
972 
973  template<typename Scalar>
975  {
976  return this->neighbor_edge;
977  }
978 
979  template<typename Scalar>
980  unsigned int NeighborSearch<Scalar>::get_central_n_trans(unsigned int index) const
981  {
982  if (central_transformations[index] != nullptr)
983  return this->central_transformations[index]->num_levels;
984  else
985  return 0;
986  }
987 
988  template<typename Scalar>
989  unsigned int NeighborSearch<Scalar>::get_central_transformations(unsigned int index_1, unsigned int index_2) const
990  {
991  if (central_transformations[index_1] == nullptr)
992  throw Hermes::Exceptions::Exception("Out of bounds of central_transformations.");
993  if (index_2 >= (unsigned)Transformations::max_level)
994  throw Hermes::Exceptions::Exception("Trying to access transformation deeper than allowed.");
995  return this->central_transformations[index_1]->transf[index_2];
996  }
997 
998  template<typename Scalar>
999  unsigned int NeighborSearch<Scalar>::get_neighbor_n_trans(unsigned int index) const
1000  {
1001  if (neighbor_transformations[index] != nullptr)
1002  return this->neighbor_transformations[index]->num_levels;
1003  else
1004  return 0;
1005  }
1006 
1007  template<typename Scalar>
1008  unsigned int NeighborSearch<Scalar>::get_neighbor_transformations(unsigned int index_1, unsigned int index_2) const
1009  {
1010  if (neighbor_transformations[index_1] == nullptr)
1011  throw Hermes::Exceptions::Exception("Out of bounds of neighbor_transformations.");
1012  if (index_2 >= (unsigned)Transformations::max_level)
1013  throw Hermes::Exceptions::Exception("Trying to access transformation deeper than allowed.");
1014  return this->neighbor_transformations[index_1]->transf[index_2];
1015  }
1016 
1017  template<typename Scalar>
1019  {
1020  Func<Scalar>* fn_central = init_fn(fu, get_quad_eo(false));
1021 
1022  uint64_t original_transform = fu->get_transform();
1023 
1024  // Change the active element of the function. Note that this also resets the transformations on the function.
1025  fu->set_active_element(neighbors[active_segment]);
1026 
1027  if (active_segment < this->neighbor_transformations_alloc_size && neighbor_transformations[active_segment])
1028  neighbor_transformations[active_segment]->apply_on(fu);
1029 
1030  Func<Scalar>* fn_neighbor = init_fn(fu, get_quad_eo(true));
1031 
1032  // Restore the original function.
1034  fu->set_transform(original_transform);
1035 
1036  return new DiscontinuousFunc<Scalar>(fn_central, fn_neighbor, (neighbor_edge.orientation == 1));
1037 
1038  //NOTE: This function is not very efficient, since it sets the active elements and possibly pushes transformations
1039  // for each mesh function in each cycle of the innermost assembly loop. This is neccessary because only in
1040  // this innermost cycle (in function DiscreteProblem::eval_form), we know the quadrature order (dependent on
1041  // the actual basis and test function), which is needed for creating the Func<Scalar> objects via init_fn.
1042  // The reason for storing the central and neighbor values of any given function in these objects is that otherwise
1043  // we would have to have one independent copy of the function for each of the neighboring elements. However, it
1044  // could unify the way PrecalcShapesets and MeshFunctions are treated in NeighborSearch and maybe these additional
1045  // deep memory copying, performed only after setting the active edge part (before the nested loops over basis and
1046  // test functions), would be actually more efficient than this. This would require implementing copy for Filters.
1047  }
1048 
1049  template<typename Scalar>
1051  central_al(central_al)
1052  {
1053  neighbor_al = new AsmList<Scalar>();
1054  space->get_boundary_assembly_list(neighborhood->neighb_el, neighborhood->neighbor_edge.local_num_of_edge, neighbor_al);
1056  }
1057 
1058  template<typename Scalar>
1060  {
1061  assert(central_al != nullptr && neighbor_al != nullptr);
1062  cnt = central_al->cnt + neighbor_al->cnt;
1063  dof = malloc_with_check<int>(cnt);
1064  memcpy(dof, central_al->dof, sizeof(int)*central_al->cnt);
1065  memcpy(dof + central_al->cnt, neighbor_al->dof, sizeof(int)*neighbor_al->cnt);
1066  }
1067 
1068  template<typename Scalar>
1070  {
1071  free_with_check(dof);
1072  delete neighbor_al;
1073  }
1074 
1075  template<typename Scalar>
1077  {
1078  delete central_al;
1079  }
1080 
1081  template<typename Scalar>
1083  {
1084  free_with_check(this->dof);
1085  space->get_boundary_assembly_list(neighborhood->neighb_el, neighborhood->neighbor_edge.local_num_of_edge, neighbor_al);
1086  combine_assembly_lists();
1087  }
1088 
1089  template<typename Scalar>
1091  {
1092  return (index >= central_al->cnt);
1093  }
1094 
1095  template<typename Scalar>
1097  {
1098  }
1099 
1100  template<typename Scalar>
1101  NeighborSearch<Scalar>::Transformations::Transformations(const Transformations* t)
1102  {
1103  copy_from(t);
1104  }
1105 
1106  template<typename Scalar>
1107  void NeighborSearch<Scalar>::Transformations::operator=(const Transformations* t)
1108  {
1109  copy_from(t);
1110  }
1111 
1112  template<typename Scalar>
1113  NeighborSearch<Scalar>::Transformations::Transformations(const std::vector<unsigned int>& t)
1114  {
1115  copy_from(t);
1116  }
1117 
1118  template<typename Scalar>
1119  void NeighborSearch<Scalar>::Transformations::copy_from(const std::vector<unsigned int>& t)
1120  {
1121  num_levels = std::min<unsigned int>(t.size(), max_level);
1122  std::copy(t.begin(), t.begin() + num_levels, transf);
1123  }
1124 
1125  template<typename Scalar>
1126  void NeighborSearch<Scalar>::Transformations::copy_from(const Transformations* t)
1127  {
1128  num_levels = t->num_levels;
1129  memcpy(transf, t->transf, max_level * sizeof(unsigned int));
1130  }
1131 
1132  template<typename Scalar>
1133  void NeighborSearch<Scalar>::Transformations::copy_to(std::vector<unsigned int>* t)
1134  {
1135  t->assign(transf, transf + num_levels);
1136  }
1137 
1138  template<typename Scalar>
1139  void NeighborSearch<Scalar>::Transformations::reset()
1140  {
1141  memset(transf, 0, num_levels * sizeof(unsigned int));
1142  num_levels = 0;
1143  }
1144 
1145  template<typename Scalar>
1146  void NeighborSearch<Scalar>::Transformations::apply_on(Transformable* tr) const
1147  {
1148  for (unsigned int i = 0; i < num_levels; i++)
1149  tr->push_transform(transf[i]);
1150  }
1151 
1152  template<typename Scalar>
1153  void NeighborSearch<Scalar>::Transformations::apply_on(const std::vector<Transformable*>& tr) const
1154  {
1155  for (std::vector<Transformable*>::const_iterator it = tr.begin(); it != tr.end(); ++it)
1156  for (unsigned int i = 0; i < num_levels; i++)
1157  (*it)->push_transform(transf[i]);
1158  }
1159 
1160  template class HERMES_API NeighborSearch < double > ;
1161  template class HERMES_API NeighborSearch < std::complex<double> > ;
1162  }
1163 }
DiscontinuousFunc< Scalar > * init_ext_fn(MeshFunction< Scalar > *fu)
std::vector< unsigned int > get_transforms(uint64_t sub_idx) const
Extract transformations in the correct direction from the provided sub_idx.
void update_according_to_sub_idx(const std::vector< unsigned int > &transformations)
Update according to the subelement mapping of the central element.
ExtendedShapeset * supported_shapes
Object allowing to set/get a particular shape function from the extended.
int id
node id number
Definition: element.h:48
Definition: adapt.h:24
Stores one element of a mesh.
Definition: element.h:107
virtual void set_active_element(Element *e)
Element * elem[2]
elements sharing the edge node
Definition: element.h:70
This class represents a function with jump discontinuity on an interface of two elements.
Definition: forms.h:335
Element * neighb_el
Currently selected neighbor element (on the other side of active segment).
unsigned int get_central_n_trans(unsigned int index) const
Returns the number(depth) of the current central transformation according to the current active segme...
void set_active_segment(unsigned int index)
Sets the active segment, neighbor element, and neighbor edge accordingly.
static const int max_level
Number of allowed transformations (or equiv. number of neighbors.
ExtendedShapeset(NeighborSearch< Scalar > *neighborhood, AsmList< Scalar > *central_al, SpaceSharedPtr< Scalar >space)
Represents a function defined on a mesh.
Definition: mesh_function.h:56
Element * central_el
Central (currently assembled) element.
Used to pass the instances of Space around.
Definition: space.h:34
bool compatible_transformations(unsigned int a, unsigned int b, int edge) const
void set_ignore_errors(bool value)
Function that sets the variable ignore_errors. See the variable description.
NeighborEdgeInfo get_neighbor_edge() const
Returns the current active neighbor edge according to the current active segment. ...
Stores one node of a mesh.
Definition: element.h:45
HERMES_API Func< double > * init_fn(PrecalcShapeset *fu, RefMap *rm, const int order)
Init the shape function for the evaluation of the volumetric/surface integral (transformation of valu...
Definition: forms.cpp:469
const std::vector< Element * > * get_neighbors() const
Transformations ** central_transformations
Array of transformations of the central element to each neighbor.
uint64_t original_central_el_transform
Sub-element transformation of any function that comes from the.
unsigned int num_levels
Number of transformation levels actually used in transf.
void delete_neighbor(unsigned int position)
In case we determine a neighbor is not correct due to a subelement mapping, we delete it...
Element * get_neighb_el() const
Returns the current neighbor element according to the current active segment.
int get_active_segment() const
Returns the current active segment.
virtual void set_transform(uint64_t idx)
Definition: function.cpp:155
Node * en[H2D_MAX_NUMBER_EDGES]
edge node pointers
Definition: element.h:138
Transformations of an element to one of its neighbors.
NeighborSearch(Element *el, MeshSharedPtr mesh)
int local_num_of_edge
Local number of the edge on neighbor element.
unsigned int transf[max_level]
Array holding the transformations at subsequent levels.
void update(NeighborSearch *neighborhood, SpaceSharedPtr< Scalar > space)
int get_quad_eo(bool on_neighbor=false) const
bool is_inter_edge(const int &edge, const std::vector< unsigned int > &transformations) const
Gives an info if edge is an intra- or inter- element edge.
Transformations ** neighbor_transformations
Array of transformations of the neighbor to the central element (go-up).
unsigned int get_neighbor_n_trans(unsigned int index) const
Returns the number(depth) of the current neighbor transformation according to the current active segm...
ExtendedShapeset * create_extended_asmlist(SpaceSharedPtr< Scalar > space, AsmList< Scalar > *al)
void handle_sub_idx_way_down(const std::vector< unsigned int > &transformations)
Special function for handling subelement transformations in the case of more than one neighboring act...
bool set_active_edge_multimesh(const int &edge)
Enhancement of set_active_edge for multimesh assembling.
unsigned int get_neighbor_transformations(unsigned int index_1, unsigned int index_2) const
Returns the current neighbor transformations according to the current active segment.
This class characterizes a neighborhood of a given edge in terms of adjacent elements and provides me...
bool active
0 = active, no sons; 1 = inactive (refined), has sons
Definition: element.h:114
AsmList< Scalar > * neighbor_al
Assembly list for the currently assembled edge on the neighbor elem.
void clear_supported_shapes()
Frees the memory occupied by the extended shapeset.
void combine_assembly_lists()
Create assembly list for the extended shapeset by joining central and neighbor element's assembly lis...
bool orientation
Relative orientation of the neighbor edge with respect to the active edge.
Element * parent
pointer to the parent element for the current son
Definition: element.h:118
Node * vn[H2D_MAX_NUMBER_VERTICES]
vertex node pointers
Definition: element.h:134
Structure containing all the needed information about the active edge from the neighbor's side...
unsigned int get_central_transformations(unsigned int index_1, unsigned int index_2) const
Returns the current central transformations according to the current active segment.