Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
tria_accessor.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: tria_accessor.templates.h 30040 2013-07-18 17:06:48Z maier @f$
3 //
4 // Copyright (C) 1999 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__tria_accessor_templates_h
18 #define __deal2__tria_accessor_templates_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/geometry_info.h>
23 #include <deal.II/base/template_constraints.h>
24 #include <deal.II/grid/tria.h>
25 #include <deal.II/grid/tria_levels.h>
26 #include <deal.II/grid/tria_faces.h>
27 #include <deal.II/grid/tria_iterator.h>
28 #include <deal.II/grid/tria_accessor.h>
29 #include <deal.II/grid/tria_iterator.templates.h>
30 #include <deal.II/distributed/tria.h>
31 
32 #include <cmath>
33 
35 
36 namespace parallel
37 {
38  namespace distributed
39  {
40  template <int, int> class Triangulation;
41  }
42 }
43 
44 
45 /*------------------------ Functions: TriaAccessorBase ---------------------------*/
46 
47 template <int structdim, int dim, int spacedim>
48 inline
50  const Triangulation<dim,spacedim> *tria,
51  const int level,
52  const int index,
53  const AccessorData *)
54  :
55  present_level((structdim==dim) ? level : 0),
56  present_index (index),
57  tria (tria)
58 {
59 
60  // non-cells have no level, so a 0
61  // should have been passed, or a -1
62  // for an end-iterator, or -2 for
63  // an invalid (default constructed)
64  // iterator
65  if (structdim != dim)
66  {
67  Assert ((level == 0) || (level == -1) || (level == -2),
69  }
70 }
71 
72 
73 template <int structdim, int dim, int spacedim>
74 inline
76  :
77  present_level(a.present_level),
78  present_index(a.present_index),
79  tria(a.tria)
80 {}
81 
82 
83 template <int structdim, int dim, int spacedim>
84 inline
85 void
87 {
88  present_level = a.present_level;
89  present_index = a.present_index;
90  tria = a.tria;
91 
92  if (structdim != dim)
93  {
94  Assert ((present_level == 0) || (present_level == -1) || (present_level == -2),
96  }
97 }
98 
99 
100 
101 template <int structdim, int dim, int spacedim>
102 inline
105 {
106  present_level = a.present_level;
107  present_index = a.present_index;
108  tria = a.tria;
109 
110  if (structdim != dim)
111  {
112  Assert ((present_level == 0) || (present_level == -1) || (present_level == -2),
113  ExcInternalError());
114  }
115  return *this;
116 }
117 
118 
119 
120 template <int structdim, int dim, int spacedim>
121 inline
122 bool
124 {
125  Assert (tria == a.tria, TriaAccessorExceptions::ExcCantCompareIterators());
126  return ((present_level == a.present_level) &&
127  (present_index == a.present_index));
128 }
129 
130 
131 
132 template <int structdim, int dim, int spacedim>
133 inline
134 bool
136 {
137  Assert (tria == a.tria, TriaAccessorExceptions::ExcCantCompareIterators());
138  return ((present_level != a.present_level) ||
139  (present_index != a.present_index));
140 }
141 
142 
143 
144 template <int structdim, int dim, int spacedim>
145 inline
146 bool
148 {
149  Assert (tria == other.tria, TriaAccessorExceptions::ExcCantCompareIterators());
150 
151  if (present_level != other.present_level)
152  return (present_level < other.present_level);
153 
154  return (present_index < other.present_index);
155 
156 }
157 
158 
159 
160 template <int structdim, int dim, int spacedim>
161 inline
162 int
164 {
165  // This is always zero or invalid
166  // if the object is not a cell
167  return present_level;
168 }
169 
170 
171 
172 template <int structdim, int dim, int spacedim>
173 inline
174 int
176 {
177  return present_index;
178 }
179 
180 
181 
182 template <int structdim, int dim, int spacedim>
183 inline
186 {
187  if ((present_level>=0) && (present_index>=0))
188  return IteratorState::valid;
189  else if (present_index==-1)
191  else
192  return IteratorState::invalid;
193 }
194 
195 
196 
197 template <int structdim, int dim, int spacedim>
198 inline
201 {
202  return *tria;
203 }
204 
205 
206 
207 template <int structdim, int dim, int spacedim>
208 inline
209 void
211 {
212  // this iterator is used for
213  // objects without level
214  ++this->present_index;
215 
216  if (structdim != dim)
217  {
218  // is index still in the range of
219  // the vector? (note that we don't
220  // have to set the level, since
221  // dim!=1 and the object therefore
222  // has no level)
223  if (this->present_index
224  >=
225  static_cast<int>(objects().cells.size()))
226  this->present_index = -1;
227  }
228  else
229  {
230  while (this->present_index
231  >=
232  static_cast<int>(this->tria->levels[this->present_level]->cells.cells.size()))
233  {
234  // no -> go one level up until we find
235  // one with more than zero cells
236  ++this->present_level;
237  this->present_index = 0;
238  // highest level reached?
239  if (this->present_level >= static_cast<int>(this->tria->levels.size()))
240  {
241  // return with past the end pointer
242  this->present_level = this->present_index = -1;
243  return;
244  }
245  }
246  }
247 }
248 
249 
250 template <int structdim, int dim, int spacedim>
251 inline
252 void
254 {
255  // same as operator++
256  --this->present_index;
257 
258  if (structdim != dim)
259  {
260  if (this->present_index < 0)
261  this->present_index = -1;
262  }
263  else
264  {
265  while (this->present_index < 0)
266  {
267  // no -> go one level down
268  --this->present_level;
269  // lowest level reached?
270  if (this->present_level == -1)
271  {
272  // return with past the end pointer
273  this->present_level = this->present_index = -1;
274  return;
275  }
276  // else
277  this->present_index = this->tria->levels[this->present_level]->cells.cells.size()-1;
278  }
279  }
280 }
281 
282 
283 namespace internal
284 {
285  namespace TriaAccessorBase
286  {
292  template <int dim>
293  inline
294  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<1> > *
295  get_objects (::internal::Triangulation::TriaFaces<dim> *faces,
296  const ::internal::int2type<1>)
297  {
298  return &faces->lines;
299  }
300 
301 
302  template <int dim>
303  inline
304  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<2> > *
305  get_objects (::internal::Triangulation::TriaFaces<dim> *faces,
306  const ::internal::int2type<2>)
307  {
308  return &faces->quads;
309  }
310 
311  inline
312  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<1> > *
313  get_objects (::internal::Triangulation::TriaFaces<1> *,
314  const ::internal::int2type<1>)
315  {
316  Assert (false, ExcInternalError());
317  return 0;
318  }
319 
320  inline
321  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<2> > *
322  get_objects (::internal::Triangulation::TriaFaces<2> *,
323  const ::internal::int2type<2>)
324  {
325  Assert (false, ExcInternalError());
326  return 0;
327  }
328 
329  inline
330  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<3> > *
331  get_objects (::internal::Triangulation::TriaFaces<3> *,
332  const ::internal::int2type<3>)
333  {
334  Assert (false, ExcInternalError());
335  return 0;
336  }
337 
343  template <int dim>
344  inline
345  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<3> > *
346  get_objects (::internal::Triangulation::TriaFaces<dim> *faces,
347  const ::internal::int2type<3>)
348  {
349  Assert (false, ExcInternalError());
350  return 0;
351  }
352 
357  template <int structdim, int dim>
358  inline
359  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<structdim> > *
361  const ::internal::int2type<structdim>)
362  {
363  Assert (false, ExcInternalError());
364  return 0;
365  }
366 
367  template <int dim>
368  inline
369  ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<dim> > *
371  const ::internal::int2type<dim>)
372  {
373  return cells;
374  }
375  }
376 }
377 
378 
379 
380 template <int structdim, int dim, int spacedim>
381 inline
382 ::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject<structdim> > &
384 {
385  if (structdim != dim)
386  // get sub-objects. note that the
387  // current class is only used for
388  // objects that are *not* cells
389  return *::internal::TriaAccessorBase::get_objects (this->tria->faces,
391  else
392  return *::internal::TriaAccessorBase::get_objects (&this->tria->levels[this->present_level]->cells,
394 }
395 
396 
397 
398 /*------------------------ Functions: InvalidAccessor ---------------------------*/
399 
400 template <int structdim, int dim, int spacedim>
403  const int ,
404  const int ,
405  const AccessorData *)
406 {
407  Assert (false,
408  ExcMessage ("You are attempting an illegal conversion between "
409  "iterator/accessor types. The constructor you call "
410  "only exists to make certain template constructs "
411  "easier to write as dimension independent code but "
412  "the conversion is not valid in the current context."));
413 }
414 
415 
416 
417 template <int structdim, int dim, int spacedim>
420  :
421  TriaAccessorBase<structdim,dim,spacedim> (static_cast<const TriaAccessorBase<structdim,dim,spacedim>&>(i))
422 {
423  Assert (false,
424  ExcMessage ("You are attempting an illegal conversion between "
425  "iterator/accessor types. The constructor you call "
426  "only exists to make certain template constructs "
427  "easier to write as dimension independent code but "
428  "the conversion is not valid in the current context."));
429 }
430 
431 
432 
433 template <int structdim, int dim, int spacedim>
434 void
437 {
438  // nothing to do here. we could
439  // throw an exception but we can't
440  // get here without first creating
441  // an object which would have
442  // already thrown
443 }
444 
445 
446 
447 template <int structdim, int dim, int spacedim>
448 bool
451 {
452  // nothing to do here. we could
453  // throw an exception but we can't
454  // get here without first creating
455  // an object which would have
456  // already thrown
457  return false;
458 }
459 
460 
461 
462 template <int structdim, int dim, int spacedim>
463 bool
465 operator != (const InvalidAccessor &) const
466 {
467  // nothing to do here. we could
468  // throw an exception but we can't
469  // get here without first creating
470  // an object which would have
471  // already thrown
472  return true;
473 }
474 
475 
476 
477 template <int structdim, int dim, int spacedim>
478 bool
480 {
481  // nothing to do here. we could
482  // throw an exception but we can't
483  // get here without first creating
484  // an object which would have
485  // already thrown
486  return false;
487 }
488 
489 
490 
491 template <int structdim, int dim, int spacedim>
492 bool
494 {
495  // nothing to do here. we could
496  // throw an exception but we can't
497  // get here without first creating
498  // an object which would have
499  // already thrown
500  return false;
501 }
502 
503 
504 
505 template <int structdim, int dim, int spacedim>
506 void
508 {}
509 
510 
511 
512 template <int structdim, int dim, int spacedim>
513 void
515 {}
516 
517 
518 
519 /*------------------------ Functions: TriaAccessor ---------------------------*/
520 
521 
522 namespace internal
523 {
524  namespace TriaAccessor
525  {
526  // make sure that if in the following we
527  // write TriaAccessor
528  // we mean the *class*
529  // ::TriaAccessor, not the
530  // enclosing namespace
531  // ::internal::TriaAccessor
532  using ::TriaAccessor;
533 
539  {
545  template <int dim, int spacedim>
546  static
547  unsigned int
549  const unsigned int)
550  {
551  Assert (false,
552  ExcMessage ("You can't ask for the index of a line bounding "
553  "a one-dimensional cell because it is not "
554  "bounded by lines."));
556  }
557 
558 
559  template <int dim, int spacedim>
560  static
561  unsigned int
563  const unsigned int i)
564  {
565  return accessor.objects().cells[accessor.present_index].face(i);
566  }
567 
568 
569  template <int dim, int spacedim>
570  static
571  unsigned int
573  const unsigned int i)
574  {
575  // get the line index by asking the
576  // quads. first assume standard orientation
577  //
578  // set up a table that for each
579  // line describes a) from which
580  // quad to take it, b) which line
581  // therein it is if the face is
582  // oriented correctly
583  static const unsigned int lookup_table[12][2] =
584  {
585  { 4, 0 }, // take first four lines from bottom face
586  { 4, 1 },
587  { 4, 2 },
588  { 4, 3 },
589 
590  { 5, 0 }, // second four lines from top face
591  { 5, 1 },
592  { 5, 2 },
593  { 5, 3 },
594 
595  { 0, 0 }, // the rest randomly
596  { 1, 0 },
597  { 0, 1 },
598  { 1, 1 }
599  };
600 
601  // respect non-standard faces by calling the
602  // reordering function from GeometryInfo
603 
604  const unsigned int quad_index=lookup_table[i][0];
605  const unsigned int std_line_index=lookup_table[i][1];
606 
608  std_line_index,
609  accessor.face_orientation(quad_index),
610  accessor.face_flip(quad_index),
611  accessor.face_rotation(quad_index));
612 
613  return (accessor.quad(quad_index)->line_index(line_index));
614  }
615 
616 
617 
623  template <int structdim, int dim, int spacedim>
624  static
625  unsigned int
627  const unsigned int)
628  {
629  Assert (false,
630  ExcMessage ("You can't ask for the index of a quad bounding "
631  "a one- or two-dimensional cell because it is not "
632  "bounded by quads."));
634  }
635 
636 
637  template <int dim, int spacedim>
638  static
639  unsigned int
641  const unsigned int i)
642  {
645  return accessor.tria->levels[accessor.present_level]
646  ->cells.cells[accessor.present_index].face(i);
647  }
648 
649 
650 
654  template <int structdim, int dim, int spacedim>
655  static
656  bool
658  const unsigned int)
659  {
660  /*
661  * Default implementation used in 1d and 2d
662  *
663  * In 1d and 2d, face_orientation is always true
664  */
665 
666  return true;
667  }
668 
669 
670  template <int dim, int spacedim>
671  static
672  bool
674  const unsigned int face)
675  {
676  return (accessor.tria->levels[accessor.present_level]
677  ->cells.face_orientation(accessor.present_index, face));
678  }
679 
680 
681 
687  template <int structdim, int dim, int spacedim>
688  static
689  bool
691  const unsigned int)
692  {
693  /*
694  * Default implementation used in 1d and 2d
695  *
696  * In 1d, face_flip is always false as there is no such concept as
697  * "flipped" faces in 1d.
698  *
699  * In 2d, we currently only support meshes where all faces are in
700  * standard orientation, so the result is also false. This also
701  * matches the fact that one can *always* orient faces in 2d in such a
702  * way that the don't need to be flipped
703  */
704  return false;
705 
706  }
707 
708 
709  template <int dim, int spacedim>
710  static
711  bool
712  face_flip (const TriaAccessor<3, dim, spacedim> &accessor,
713  const unsigned int face)
714  {
718  < accessor.tria->levels[accessor.present_level]
719  ->cells.face_flips.size(),
720  ExcInternalError());
721 
722  return (accessor.tria->levels[accessor.present_level]
723  ->cells.face_flips[accessor.present_index *
725  + face]);
726  }
727 
728 
729 
735  template <int structdim, int dim, int spacedim>
736  static
737  bool
739  const unsigned int)
740  {
741  /*
742  * Default implementation used in 1d and 2d
743  *
744  * In 1d and 2d, face_rotation is always false as there is no such
745  * concept as "rotated" faces in 1d and 2d.
746  */
747  return false;
748  }
749 
750 
751  template <int dim, int spacedim>
752  static
753  bool
755  const unsigned int face)
756  {
760  < accessor.tria->levels[accessor.present_level]
761  ->cells.face_rotations.size(),
762  ExcInternalError());
763 
764  return (accessor.tria->levels[accessor.present_level]
765  ->cells.face_rotations[accessor.present_index *
767  + face]);
768  }
769 
775  template <int dim, int spacedim>
776  static
777  bool
779  const unsigned int)
780  {
781  return true;
782  }
783 
784 
785  template <int spacedim>
786  static
787  bool
789  const unsigned int)
790  {
791  // quads in 2d have no non-standard orientation
792  return true;
793  }
794 
795 
796  template <int spacedim>
797  static
798  bool
800  const unsigned int line)
801  {
802  // quads as part of 3d hexes can have non-standard orientation
803 
804  //TODO: why is this face_orientation, not line_orientation as in the setter function?
805  return accessor.tria->faces->quads.face_orientation(accessor.present_index, line);
806  }
807 
808 
809  template <int dim, int spacedim>
810  static
811  bool
813  const unsigned int line)
814  {
815  Assert (accessor.used(), TriaAccessorExceptions::ExcCellNotUsed());
818 
819  // get the line index by asking the
820  // quads. first assume standard orientation
821  //
822  // set up a table that for each
823  // line describes a) from which
824  // quad to take it, b) which line
825  // therein it is if the face is
826  // oriented correctly
827  static const unsigned int lookup_table[12][2] =
828  {
829  { 4, 0 }, // take first four lines from bottom face
830  { 4, 1 },
831  { 4, 2 },
832  { 4, 3 },
833 
834  { 5, 0 }, // second four lines from top face
835  { 5, 1 },
836  { 5, 2 },
837  { 5, 3 },
838 
839  { 0, 0 }, // the rest randomly
840  { 1, 0 },
841  { 0, 1 },
842  { 1, 1 }
843  };
844 
845  const unsigned int quad_index=lookup_table[line][0];
846  const unsigned int std_line_index=lookup_table[line][1];
847 
848  const unsigned int line_index=GeometryInfo<dim>::standard_to_real_face_line(
849  std_line_index,
850  accessor.face_orientation(quad_index),
851  accessor.face_flip(quad_index),
852  accessor.face_rotation(quad_index));
853 
854  // now we got to the correct line and ask
855  // the quad for its line_orientation. however, if
856  // the face is rotated, it might be possible,
857  // that a standard orientation of the line
858  // with respect to the face corrsponds to a
859  // non-standard orientation for the line with
860  // respect to the cell.
861  //
862  // set up a table indicating if the two
863  // standard orientations coincide
864  //
865  // first index: two pairs of lines 0(lines
866  // 0/1) and 1(lines 2/3)
867  //
868  // second index: face_orientation; 0:
869  // opposite normal, 1: standard
870  //
871  // third index: face_flip; 0: standard, 1:
872  // face rotated by 180 degrees
873  //
874  // forth index: face_rotation: 0: standard,
875  // 1: face rotated by 90 degrees
876 
877  static const bool bool_table[2][2][2][2] =
878  {
879  { { { true, false }, // lines 0/1, face_orientation=false, face_flip=false, face_rotation=false and true
880  { false, true }
881  }, // lines 0/1, face_orientation=false, face_flip=true, face_rotation=false and true
882  { { true, true }, // lines 0/1, face_orientation=true, face_flip=false, face_rotation=false and true
883  { false, false }
884  }
885  },// linea 0/1, face_orientation=true, face_flip=true, face_rotation=false and true
886 
887  { { { true, true }, // lines 2/3 ...
888  { false, false }
889  },
890  { { true, false },
891  { false, true }
892  }
893  }
894  };
895 
896 
897  return (accessor.quad(quad_index)
898  ->line_orientation(line_index)
899  == bool_table[std_line_index/2]
900  [accessor.face_orientation(quad_index)]
901  [accessor.face_flip(quad_index)]
902  [accessor.face_rotation(quad_index)]);
903  }
904 
905 
906 
912  template <int structdim, int dim, int spacedim>
913  static
914  void
916  const unsigned int,
917  const bool)
918  {
919  Assert (false, ExcInternalError());
920  }
921 
922 
923  template <int dim, int spacedim>
924  static
925  void
927  const unsigned int face,
928  const bool value)
929  {
930  Assert (accessor.used(), TriaAccessorExceptions::ExcCellNotUsed());
934  < accessor.tria->levels[accessor.present_level]
935  ->cells.face_orientations.size(),
936  ExcInternalError());
937  accessor.tria->levels[accessor.present_level]
938  ->cells.face_orientations[accessor.present_index *
940  +
941  face] = value;
942  }
943 
944 
945 
951  template <int structdim, int dim, int spacedim>
952  static
953  void
955  const unsigned int,
956  const bool)
957  {
958  Assert (false, ExcInternalError());
959  }
960 
961 
962  template <int dim, int spacedim>
963  static
964  void
966  const unsigned int face,
967  const bool value)
968  {
972  < accessor.tria->levels[accessor.present_level]
973  ->cells.face_flips.size(),
974  ExcInternalError());
975 
976  accessor.tria->levels[accessor.present_level]
977  ->cells.face_flips[accessor.present_index *
979  + face] = value;
980  }
981 
982 
983 
989  template <int structdim, int dim, int spacedim>
990  static
991  void
993  const unsigned int,
994  const bool)
995  {
996  Assert (false, ExcInternalError());
997  }
998 
999 
1000  template <int dim, int spacedim>
1001  static
1002  void
1004  const unsigned int face,
1005  const bool value)
1006  {
1010  < accessor.tria->levels[accessor.present_level]
1011  ->cells.face_rotations.size(),
1012  ExcInternalError());
1013 
1014  accessor.tria->levels[accessor.present_level]
1015  ->cells.face_rotations[accessor.present_index *
1017  + face] = value;
1018  }
1019 
1025  template <int dim, int spacedim>
1026  static
1027  void
1029  const unsigned int,
1030  const bool)
1031  {
1032  Assert (false, ExcInternalError());
1033  }
1034 
1035 
1036  template <int spacedim>
1037  static
1038  void
1040  const unsigned int,
1041  const bool)
1042  {
1043  // quads in 2d have no
1044  // non-standard orientation
1045  Assert (false, ExcInternalError());
1046  }
1047 
1048 
1049  template <int spacedim>
1050  static
1051  void
1053  const unsigned int line,
1054  const bool value)
1055  {
1056  Assert (accessor.used(), TriaAccessorExceptions::ExcCellNotUsed());
1060  < accessor.tria->faces->quads.line_orientations.size(),
1061  ExcInternalError());
1062  // quads as part of 3d hexes
1063  // can have non-standard
1064  // orientation
1065  accessor.tria->faces->quads.line_orientations[accessor.present_index *
1067  + line]
1068  = value;
1069  }
1070 
1071 
1072  template <int dim, int spacedim>
1073  static
1074  void
1076  const unsigned int,
1077  const bool)
1078  {
1079  // it seems like we don't need this
1080  // one
1081  Assert (false, ExcNotImplemented());
1082  }
1083 
1084 
1089  template <int dim, int spacedim>
1090  static
1091  unsigned int
1093  const unsigned int corner)
1094  {
1095  return accessor.objects().cells[accessor.present_index].face (corner);
1096  }
1097 
1098 
1099  template <int dim, int spacedim>
1100  static
1101  unsigned int
1102  vertex_index (const TriaAccessor<2,dim,spacedim> &accessor,
1103  const unsigned int corner)
1104  {
1105  // table used to switch the vertices, if the
1106  // line orientation is wrong,
1107  //
1108  // first index: line orientation 0: false or
1109  // 1: true=standard
1110  //
1111  // second index: vertex index to be switched
1112  // (or not)
1113 
1114  static const unsigned int switch_table[2][2]= {{1,0},{0,1}};
1115 
1116  return accessor.line(corner%2)
1117  ->vertex_index(switch_table[accessor.line_orientation(corner%2)][corner/2]);
1118  }
1119 
1120 
1121 
1122  template <int dim, int spacedim>
1123  static
1124  unsigned int
1125  vertex_index (const TriaAccessor<3,dim,spacedim> &accessor,
1126  const unsigned int corner)
1127  {
1128  // get the corner indices by asking either
1129  // the bottom or the top face for its
1130  // vertices. handle non-standard faces by
1131  // calling the vertex reordering function
1132  // from GeometryInfo
1133 
1134  // bottom face (4) for first four vertices,
1135  // top face (5) for the rest
1136  const unsigned int face_index=4+corner/4;
1137 
1138  return accessor.quad(face_index)
1139  ->vertex_index(GeometryInfo<dim>
1140  ::standard_to_real_face_vertex(corner%4,
1141  accessor.face_orientation(face_index),
1142  accessor.face_flip(face_index),
1143  accessor.face_rotation(face_index)));
1144  }
1145  };
1146  }
1147 }
1148 
1149 
1150 
1151 
1152 template <int structdim, int dim, int spacedim>
1153 inline
1156  const int level,
1157  const int index,
1158  const AccessorData *local_data)
1159  :
1160  TriaAccessorBase<structdim,dim,spacedim> (parent, level, index, local_data)
1161 {}
1162 
1163 
1164 
1165 template <int structdim, int dim, int spacedim>
1166 inline
1167 bool
1169 {
1170  Assert (this->state() == IteratorState::valid,
1171  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1172  return this->objects().used[this->present_index];
1173 }
1174 
1175 
1176 
1177 template <int structdim, int dim, int spacedim>
1178 unsigned int
1180 vertex_index (const unsigned int corner) const
1181 {
1184 
1185  return ::internal::TriaAccessor::Implementation::vertex_index (*this, corner);
1186 }
1187 
1188 
1189 
1190 template <int structdim, int dim, int spacedim>
1193 {
1194  return const_cast<Point<spacedim> &> (this->tria->vertices[vertex_index(i)]);
1195 }
1196 
1197 
1198 
1199 template <int structdim, int dim, int spacedim>
1200 inline
1201 typename ::internal::Triangulation::Iterators<dim,spacedim>::line_iterator
1203 {
1204  // checks happen in line_index
1205  return typename ::internal::Triangulation::Iterators<dim,spacedim>::line_iterator
1206  (this->tria, 0, line_index (i));
1207 }
1208 
1209 
1210 
1211 template <int structdim, int dim, int spacedim>
1212 inline
1213 unsigned int
1215 {
1218 
1219  return ::internal::TriaAccessor::Implementation::line_index (*this, i);
1220 }
1221 
1222 
1223 
1224 
1225 template <int structdim, int dim, int spacedim>
1226 inline
1227 typename ::internal::Triangulation::Iterators<dim,spacedim>::quad_iterator
1229 {
1230  // checks happen in quad_index
1231  return typename ::internal::Triangulation::Iterators<dim,spacedim>::quad_iterator
1232  (this->tria, 0, quad_index (i));
1233 }
1234 
1235 
1236 
1237 template <int structdim, int dim, int spacedim>
1238 inline
1239 unsigned int
1241 {
1242  return ::internal::TriaAccessor::Implementation::quad_index (*this, i);
1243 }
1244 
1245 
1246 
1247 template <int structdim, int dim, int spacedim>
1248 inline
1249 bool
1251 {
1252  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1253 
1254  return ::internal::TriaAccessor::Implementation::face_orientation (*this, face);
1255 }
1256 
1257 
1258 
1259 template <int structdim, int dim, int spacedim>
1260 inline
1261 bool
1263 {
1264  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1265 
1266  return ::internal::TriaAccessor::Implementation::face_flip (*this, face);
1267 }
1268 
1269 
1270 template <int structdim, int dim, int spacedim>
1271 inline
1272 bool
1274 {
1275  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1276 
1277  return ::internal::TriaAccessor::Implementation::face_rotation (*this, face);
1278 }
1279 
1280 
1281 
1282 template <int structdim, int dim, int spacedim>
1283 inline
1284 bool
1286 {
1287  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1290 
1291  return ::internal::TriaAccessor::Implementation::line_orientation (*this, line);
1292 }
1293 
1294 
1295 
1296 template <int structdim, int dim, int spacedim>
1297 inline
1298 void
1300  const bool value) const
1301 {
1302  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1303 
1305 }
1306 
1307 
1308 
1309 template <int structdim, int dim, int spacedim>
1310 inline
1311 void
1313  const bool value) const
1314 {
1315  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1316 
1318 }
1319 
1320 
1321 template <int structdim, int dim, int spacedim>
1322 inline
1323 void
1325  const bool value) const
1326 {
1327  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1328 
1330 }
1331 
1332 
1333 
1334 template <int structdim, int dim, int spacedim>
1335 inline
1336 void
1338  const bool value) const
1339 {
1340  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
1343 
1345 }
1346 
1347 
1348 
1349 template <int structdim, int dim, int spacedim>
1351 {
1352  Assert (this->state() == IteratorState::valid,
1353  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1354  this->objects().used[this->present_index] = true;
1355 }
1356 
1357 
1358 
1359 template <int structdim, int dim, int spacedim>
1360 void
1362 {
1363  Assert (this->state() == IteratorState::valid,
1364  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1365  this->objects().used[this->present_index] = false;
1366 }
1367 
1368 
1369 template <int structdim, int dim, int spacedim>
1370 int
1373 {
1374  Assert (this->present_level > 0, TriaAccessorExceptions::ExcCellHasNoParent ());
1375 
1376  // the parent of two consecutive cells
1377  // is stored only once, since it is
1378  // the same
1379  return this->tria->levels[this->present_level]->parents[this->present_index / 2];
1380 }
1381 
1382 
1383 template <int structdim, int dim, int spacedim>
1384 int
1386 child_index (const unsigned int i) const
1387 {
1388  Assert (has_children(), TriaAccessorExceptions::ExcCellHasNoChildren());
1389  Assert (i<n_children(),
1390  ExcIndexRange (i, 0, n_children()));
1391 
1392  // each set of two children are stored
1393  // consecutively, so we only have to find
1394  // the location of the set of children
1395  const unsigned int n_sets_of_two = GeometryInfo<structdim>::max_children_per_cell/2;
1396  return this->objects().children[n_sets_of_two*this->present_index+i/2]+i%2;
1397 }
1398 
1399 
1400 
1401 template <int structdim, int dim, int spacedim>
1402 int
1404 isotropic_child_index (const unsigned int i) const
1405 {
1408 
1409  switch (structdim)
1410  {
1411  case 1:
1412  return child_index (i);
1413  case 2:
1414  {
1415  const RefinementCase<2>
1416  this_refinement_case (static_cast<unsigned char>(refinement_case()));
1417 
1418  Assert (this_refinement_case != RefinementCase<2>::no_refinement,
1419  TriaAccessorExceptions::ExcCellHasNoChildren());
1420 
1421  if (this_refinement_case == RefinementCase<2>::cut_xy)
1422  return child_index(i);
1423  else if ((this_refinement_case == RefinementCase<2>::cut_x)
1424  &&
1425  (child(i%2)->refinement_case()==RefinementCase<2>::cut_y))
1426  return child(i%2)->child_index(i/2);
1427  else if ((this_refinement_case == RefinementCase<2>::cut_y)
1428  &&
1429  (child(i/2)->refinement_case()==RefinementCase<2>::cut_x))
1430  return child(i/2)->child_index(i%2);
1431  else
1432  Assert(false,
1433  ExcMessage("This cell has no grandchildren equivalent to isotropic refinement"));
1434  }
1435 
1436  case 3:
1437  Assert (false, ExcNotImplemented());
1438  }
1439  return -1;
1440 }
1441 
1442 
1443 
1444 template <int structdim, int dim, int spacedim>
1447 {
1448  Assert (this->state() == IteratorState::valid,
1449  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1450 
1451  switch (structdim)
1452  {
1453  case 1:
1455  (this->objects().children[this->present_index] != -1
1456  ?
1457  // cast the branches
1458  // here first to uchar
1459  // and then (above) to
1460  // RefinementCase<structdim>
1461  // so that the
1462  // conversion is valid
1463  // even for the case
1464  // structdim>1 (for
1465  // which this part of
1466  // the code is dead
1467  // anyway)
1468  static_cast<unsigned char>(RefinementCase<1>::cut_x) :
1469  static_cast<unsigned char>(RefinementCase<1>::no_refinement)));
1470 
1471  default:
1472  Assert (static_cast<unsigned int> (this->present_index) <
1473  this->objects().refinement_cases.size(),
1474  ExcIndexRange(this->present_index, 0,
1475  this->objects().refinement_cases.size()));
1476 
1477  return (static_cast<RefinementCase<structdim> >
1478  (this->objects().refinement_cases[this->present_index]));
1479  }
1480 }
1481 
1482 
1483 
1484 
1485 template <int structdim, int dim, int spacedim>
1486 inline
1489 
1490 {
1491  // checking of 'i' happens in child_index
1493  q (this->tria,
1494  (dim == structdim ? this->level() + 1 : 0),
1495  child_index (i));
1496 
1497  Assert ((q.state() == IteratorState::past_the_end) || q->used(),
1498  TriaAccessorExceptions::ExcUnusedCellAsChild());
1499 
1500  return q;
1501 }
1502 
1503 
1504 
1505 template <int structdim, int dim, int spacedim>
1506 inline
1509 isotropic_child (const unsigned int i) const
1510 {
1511  // checking of 'i' happens in child() or
1512  // child_index() called below
1513  switch (structdim)
1514  {
1515  case 1:
1516  // no anisotropic refinement in 1D
1517  return child(i);
1518 
1519  case 2:
1520  {
1521  const RefinementCase<2>
1522  this_refinement_case (static_cast<unsigned char>(refinement_case()));
1523 
1524  Assert (this_refinement_case != RefinementCase<2>::no_refinement,
1525  TriaAccessorExceptions::ExcCellHasNoChildren());
1526 
1527  if (this_refinement_case == RefinementCase<2>::cut_xy)
1528  return child(i);
1529  else if ((this_refinement_case == RefinementCase<2>::cut_x)
1530  &&
1531  (child(i%2)->refinement_case()==RefinementCase<2>::cut_y))
1532  return child(i%2)->child(i/2);
1533  else if ((this_refinement_case == RefinementCase<2>::cut_y)
1534  &&
1535  (child(i/2)->refinement_case()==RefinementCase<2>::cut_x))
1536  return child(i/2)->child(i%2);
1537  else
1538  Assert(false,
1539  ExcMessage("This cell has no grandchildren equivalent to isotropic refinement"));
1540  }
1541 
1542  default:
1543  Assert (false, ExcNotImplemented());
1544  }
1545  // we don't get here but have to return
1546  // something...
1547  return child(0);
1548 }
1549 
1550 
1551 
1552 template <int structdim, int dim, int spacedim>
1553 inline
1554 bool
1556 {
1557  Assert (this->state() == IteratorState::valid,
1558  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1559 
1560  // each set of two children are stored
1561  // consecutively, so we only have to find
1562  // the location of the set of children
1563  const unsigned int n_sets_of_two = GeometryInfo<structdim>::max_children_per_cell/2;
1564  return (this->objects().children[n_sets_of_two * this->present_index] != -1);
1565 }
1566 
1567 
1568 
1569 
1570 template <int structdim, int dim, int spacedim>
1571 inline
1572 unsigned int
1574 {
1575  return GeometryInfo<structdim>::n_children(refinement_case());
1576 }
1577 
1578 
1579 
1580 template <int structdim, int dim, int spacedim>
1581 inline
1582 void
1584 set_refinement_case (const RefinementCase<structdim> &refinement_case) const
1585 {
1586  Assert (this->state() == IteratorState::valid,
1587  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1588  Assert (static_cast<unsigned int> (this->present_index) <
1589  this->objects().refinement_cases.size(),
1590  ExcIndexRange(this->present_index, 0,
1591  this->objects().refinement_cases.size()));
1592 
1593  this->objects().refinement_cases[this->present_index] = refinement_case;
1594 }
1595 
1596 
1597 template <int structdim, int dim, int spacedim>
1598 inline
1599 void
1601 {
1602  Assert (this->state() == IteratorState::valid,
1603  TriaAccessorExceptions::ExcDereferenceInvalidObject());
1604  Assert (static_cast<unsigned int> (this->present_index) <
1605  this->objects().refinement_cases.size(),
1606  ExcIndexRange(this->present_index, 0,
1607  this->objects().refinement_cases.size()));
1608 
1609  this->objects().refinement_cases[this->present_index]
1611 }
1612 
1613 
1614 
1615 
1616 
1617 template <int structdim, int dim, int spacedim>
1618 void
1620  const int index) const
1621 {
1622  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1623  Assert (i%2==0, TriaAccessorExceptions::ExcSetOnlyEvenChildren(i));
1624 
1625  // each set of two children are stored
1626  // consecutively, so we only have to find
1627  // the location of the set of children
1628  const unsigned int n_sets_of_two = GeometryInfo<structdim>::max_children_per_cell/2;
1629 
1630  Assert ((index==-1) ||
1631  (i==0 && !this->has_children() && (index>=0)) ||
1632  (i>0 && this->has_children() && (index>=0) &&
1633  this->objects().children[n_sets_of_two*this->present_index+i/2] == -1),
1634  TriaAccessorExceptions::ExcCantSetChildren(index));
1635 
1636  this->objects().children[n_sets_of_two*this->present_index+i/2] = index;
1637 }
1638 
1639 
1640 
1641 template <int structdim, int dim, int spacedim>
1642 void
1644 {
1645  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1646  Assert (this->present_level > 0, TriaAccessorExceptions::ExcCellHasNoParent ());
1647  this->tria->levels[this->present_level]->parents[this->present_index / 2]
1648  = parent_index;
1649 }
1650 
1651 
1652 
1653 template <int structdim, int dim, int spacedim>
1654 void
1656 {
1657  // each set of two children are stored
1658  // consecutively, so we only have to find
1659  // the location of the set of children
1660  const unsigned int n_sets_of_two = GeometryInfo<structdim>::max_children_per_cell/2;
1661 
1662  for (unsigned int i=0; i<n_sets_of_two; ++i)
1663  set_children (2*i,-1);
1664 }
1665 
1666 
1667 
1668 template <int structdim, int dim, int spacedim>
1669 inline
1670 bool
1672 {
1673  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1674  return this->objects().user_flags[this->present_index];
1675 }
1676 
1677 
1678 
1679 template <int structdim, int dim, int spacedim>
1680 inline
1681 void
1683 {
1684  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1685  this->objects().user_flags[this->present_index] = true;
1686 }
1687 
1688 
1689 
1690 template <int structdim, int dim, int spacedim>
1691 inline
1692 void
1694 {
1695  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1696  this->objects().user_flags[this->present_index] = false;
1697 }
1698 
1699 
1700 
1701 template <int structdim, int dim, int spacedim>
1703 {
1704  set_user_flag ();
1705 
1706  if (this->has_children())
1707  for (unsigned int c=0; c<this->n_children(); ++c)
1708  this->child(c)->recursively_set_user_flag ();
1709 }
1710 
1711 
1712 
1713 template <int structdim, int dim, int spacedim>
1715 {
1716  clear_user_flag ();
1717 
1718  if (this->has_children())
1719  for (unsigned int c=0; c<this->n_children(); ++c)
1720  this->child(c)->recursively_clear_user_flag ();
1721 }
1722 
1723 
1724 
1725 template <int structdim, int dim, int spacedim>
1727 {
1728  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1729  this->objects().clear_user_data(this->present_index);
1730 }
1731 
1732 
1733 
1734 template <int structdim, int dim, int spacedim>
1736 {
1737  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1738  this->objects().user_pointer(this->present_index) = p;
1739 }
1740 
1741 
1742 
1743 template <int structdim, int dim, int spacedim>
1745 {
1746  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1747  this->objects().user_pointer(this->present_index) = 0;
1748 }
1749 
1750 
1751 
1752 template <int structdim, int dim, int spacedim>
1754 {
1755  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1756  return this->objects().user_pointer(this->present_index);
1757 }
1758 
1759 
1760 
1761 template <int structdim, int dim, int spacedim>
1762 void
1764 {
1765  set_user_pointer (p);
1766 
1767  if (this->has_children())
1768  for (unsigned int c=0; c<this->n_children(); ++c)
1769  this->child(c)->recursively_set_user_pointer (p);
1770 }
1771 
1772 
1773 
1774 template <int structdim, int dim, int spacedim>
1775 void
1777 {
1778  clear_user_pointer ();
1779 
1780  if (this->has_children())
1781  for (unsigned int c=0; c<this->n_children(); ++c)
1782  this->child(c)->recursively_clear_user_pointer ();
1783 }
1784 
1785 
1786 
1787 template <int structdim, int dim, int spacedim>
1789 {
1790  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1791  this->objects().user_index(this->present_index) = p;
1792 }
1793 
1794 
1795 
1796 template <int structdim, int dim, int spacedim>
1798 {
1799  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1800  this->objects().user_index(this->present_index) = 0;
1801 }
1802 
1803 
1804 
1805 template <int structdim, int dim, int spacedim>
1807 {
1808  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1809  return this->objects().user_index(this->present_index);
1810 }
1811 
1812 
1813 
1814 template <int structdim, int dim, int spacedim>
1815 void
1817 {
1818  set_user_index (p);
1819 
1820  if (this->has_children())
1821  for (unsigned int c=0; c<this->n_children(); ++c)
1822  this->child(c)->recursively_set_user_index (p);
1823 }
1824 
1825 
1826 
1827 template <int structdim, int dim, int spacedim>
1828 void
1830 {
1831  clear_user_index ();
1832 
1833  if (this->has_children())
1834  for (unsigned int c=0; c<this->n_children(); ++c)
1835  this->child(c)->recursively_clear_user_index ();
1836 }
1837 
1838 
1839 
1840 template <int structdim, int dim, int spacedim>
1841 inline
1842 unsigned int
1844 {
1845  if (!this->has_children())
1846  return 0;
1847 
1848  unsigned int max_depth = 1;
1849  for (unsigned int c=0; c<n_children(); ++c)
1850  max_depth = std::max (max_depth,
1851  child(c)->max_refinement_depth() + 1);
1852  return max_depth;
1853 }
1854 
1855 
1856 
1857 template <int structdim, int dim, int spacedim>
1858 unsigned int
1860 {
1861  if (!this->has_children())
1862  return 1;
1863  else
1864  {
1865  unsigned int sum = 0;
1866  for (unsigned int c=0; c<n_children(); ++c)
1867  sum += this->child(c)->number_of_children();
1868  return sum;
1869  }
1870 }
1871 
1872 
1873 
1874 template <int structdim, int dim, int spacedim>
1877 {
1878  Assert (structdim<dim, ExcImpossibleInDim(dim));
1879  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1880 
1881  return this->objects().boundary_or_material_id[this->present_index].boundary_id;
1882 }
1883 
1884 
1885 
1886 template <int structdim, int dim, int spacedim>
1887 void
1889 set_boundary_indicator (const types::boundary_id boundary_ind) const
1890 {
1891  Assert (structdim<dim, ExcImpossibleInDim(dim));
1892  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1893 
1894  this->objects().boundary_or_material_id[this->present_index].boundary_id = boundary_ind;
1895 }
1896 
1897 
1898 
1899 template <int structdim, int dim, int spacedim>
1900 void
1903 {
1904  set_boundary_indicator (boundary_ind);
1905 
1906  switch (structdim)
1907  {
1908  case 1:
1909  // 1d objects have no sub-objects
1910  // where we have to do anything
1911  break;
1912 
1913  case 2:
1914  // for boundary quads also set
1915  // boundary_id of bounding lines
1916  for (unsigned int i=0; i<4; ++i)
1917  this->line(i)->set_boundary_indicator (boundary_ind);
1918  break;
1919 
1920  default:
1921  Assert (false, ExcNotImplemented());
1922  }
1923 }
1924 
1925 
1926 
1927 template <int structdim, int dim, int spacedim>
1928 bool
1930 {
1931  // error checking is done
1932  // in boundary_indicator()
1933  return (boundary_indicator() != numbers::internal_face_boundary_id);
1934 }
1935 
1936 
1937 
1938 template <int structdim, int dim, int spacedim>
1939 const Boundary<dim,spacedim> &
1941 {
1942  Assert (structdim<dim, ExcImpossibleInDim(dim));
1943  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
1944 
1945  return this->tria->get_boundary(this->objects()
1946  .boundary_or_material_id[this->present_index].boundary_id);
1947 }
1948 
1949 
1950 
1951 template <int structdim, int dim, int spacedim>
1952 double
1954 {
1955  switch (structdim)
1956  {
1957  case 1:
1958  return std::sqrt((this->vertex(1)-this->vertex(0)).square());
1959  case 2:
1960  return std::sqrt(std::max((this->vertex(3)-this->vertex(0)).square(),
1961  (this->vertex(2)-this->vertex(1)).square()));
1962  case 3:
1963  return std::sqrt(std::max( std::max((this->vertex(7)-this->vertex(0)).square(),
1964  (this->vertex(6)-this->vertex(1)).square()),
1965  std::max((this->vertex(2)-this->vertex(5)).square(),
1966  (this->vertex(3)-this->vertex(4)).square()) ));
1967  default:
1968  Assert (false, ExcNotImplemented());
1969  return -1e10;
1970  }
1971 }
1972 
1973 
1974 
1975 template <int structdim, int dim, int spacedim>
1976 double
1978 {
1979  switch (structdim)
1980  {
1981  case 1:
1982  return std::sqrt((this->vertex(1)-this->vertex(0)).square());
1983  case 2:
1984  case 3:
1985  {
1986  double min = std::numeric_limits<double>::max();
1987  for (unsigned int i=0; i<GeometryInfo<structdim>::vertices_per_cell; ++i)
1988  for (unsigned int j=i+1; j<GeometryInfo<structdim>::vertices_per_cell; ++j)
1989  min = std::min(min, (this->vertex(i)-this->vertex(j)).square());
1990  return std::sqrt(min);
1991  }
1992  default:
1993  Assert (false, ExcNotImplemented());
1994  return -1e10;
1995  }
1996 }
1997 
1998 
1999 
2000 template <int structdim, int dim, int spacedim>
2003 {
2004  Point<spacedim> p;
2005  for (unsigned int v=0; v<GeometryInfo<structdim>::vertices_per_cell; ++v)
2006  p += vertex(v);
2008 }
2009 
2010 
2011 template <int structdim, int dim, int spacedim>
2012 bool
2015 {
2016  // go through the vertices and check... The
2017  // cell is a translation of the previous
2018  // one in case the distance between the
2019  // individual vertices in the two cell is
2020  // the same for all the vertices. So do the
2021  // check by first getting the distance on
2022  // the first vertex, and then checking
2023  // whether all others have the same down to
2024  // rounding errors (we have to be careful
2025  // here because the calculation of the
2026  // displacement between one cell and the
2027  // next can already result in the loss of
2028  // one or two digits), so we choose 1e-12
2029  // times the distance between the zeroth
2030  // vertices here.
2031  bool is_translation = true;
2032  const Point<spacedim> dist = o->vertex(0) - this->vertex(0);
2033  const double tol_square = 1e-24 * dist.norm_square();
2034  for (unsigned int i=1; i<GeometryInfo<structdim>::vertices_per_cell; ++i)
2035  {
2036  const Point<spacedim> dist_new = (o->vertex(i) - this->vertex(i)) - dist;
2037  if (dist_new.norm_square() > tol_square)
2038  {
2039  is_translation = false;
2040  break;
2041  }
2042  }
2043  return is_translation;
2044 }
2045 
2046 
2047 /*------------------------ Functions: TriaAccessor<0,1,spacedim> -----------------------*/
2048 
2049 template <int spacedim>
2050 inline
2053  const VertexKind vertex_kind,
2054  const unsigned int vertex_index)
2055  :
2056  tria (tria),
2057  vertex_kind (vertex_kind),
2058  global_vertex_index (vertex_index)
2059 {}
2060 
2061 
2062 
2063 template <int spacedim>
2064 inline
2067  const int level,
2068  const int index,
2069  const AccessorData *)
2070  :
2071  tria (tria),
2072  vertex_kind (interior_vertex),
2073  global_vertex_index (numbers::invalid_unsigned_int)
2074 {
2075  // in general, calling this constructor should yield an error -- users should
2076  // instead call the one immediately above. however, if you create something
2077  // like Triangulation<1>::face_iterator() then this calls the default constructor
2078  // of the iterator which calls the accessor with argument list (0,-2,-2,0), so
2079  // in this particular case accept this call and create an object that corresponds
2080  // to the default constructed (invalid) vertex accessor
2081  Assert ((level == -2) && (index == -2), ExcInternalError());
2082 }
2083 
2084 
2085 
2086 template <int spacedim>
2087 template <int structdim2, int dim2, int spacedim2>
2088 inline
2091  :
2092  tria (0),
2093  vertex_kind (interior_vertex),
2094  global_vertex_index (numbers::invalid_unsigned_int)
2095 {
2096  Assert(false, ExcImpossibleInDim(0));
2097 }
2098 
2099 
2100 
2101 template <int spacedim>
2102 template <int structdim2, int dim2, int spacedim2>
2103 inline
2106  :
2107  tria (0),
2108  vertex_kind (interior_vertex),
2109  global_vertex_index (numbers::invalid_unsigned_int)
2110 {
2111  Assert(false, ExcImpossibleInDim(0));
2112 }
2113 
2114 
2115 
2116 template <int spacedim>
2117 inline
2118 void
2120 {
2121  tria = t.tria;
2122  vertex_kind = t.vertex_kind;
2123  global_vertex_index = t.global_vertex_index;
2124 }
2125 
2126 
2127 
2128 template <int spacedim>
2129 inline
2132 {
2133  return IteratorState::valid;
2134 }
2135 
2136 
2137 template <int spacedim>
2138 inline
2139 int
2141 {
2142  return 0;
2143 }
2144 
2145 
2146 
2147 template <int spacedim>
2148 inline
2149 int
2151 {
2152  return global_vertex_index;
2153 }
2154 
2155 
2156 
2157 template <int spacedim>
2158 inline
2159 void
2161 {
2162  Assert (false, ExcNotImplemented());
2163 }
2164 
2165 
2166 template <int spacedim>
2167 inline
2168 void
2170 {
2171  Assert (false, ExcNotImplemented());
2172 }
2173 
2174 
2175 
2176 template <int spacedim>
2177 inline
2178 bool
2180 {
2181  const bool result = ((tria == t.tria)
2182  &&
2183  (global_vertex_index == t.global_vertex_index));
2184  // if we point to the same vertex,
2185  // make sure we know the same about
2186  // it
2187  if (result == true)
2188  Assert (vertex_kind == t.vertex_kind, ExcInternalError());
2189 
2190  return result;
2191 }
2192 
2193 
2194 
2195 template <int spacedim>
2196 inline
2197 bool
2199 {
2200  return !(*this==t);
2201 }
2202 
2203 
2204 
2205 
2206 template <int spacedim>
2207 inline
2208 int
2210 {
2211  return -1;
2212 }
2213 
2214 
2215 template <int spacedim>
2216 inline
2217 unsigned int
2219 {
2220  Assert(i==0, ExcIndexRange(i, 0, 1));
2221  (void)i;
2222  return global_vertex_index;
2223 }
2224 
2225 
2226 
2227 template <int spacedim>
2228 inline
2230 TriaAccessor<0, 1, spacedim>::vertex (const unsigned int i) const
2231 {
2232  Assert(i==0, ExcIndexRange(i, 0, 1));
2233  (void)i;
2234  return const_cast<Point<spacedim> &> (this->tria->vertices[global_vertex_index]);
2235 }
2236 
2237 
2238 
2239 template <int spacedim>
2240 inline
2243 {
2244  return this->tria->vertices[global_vertex_index];
2245 }
2246 
2247 
2248 
2249 template <int spacedim>
2250 inline
2251 typename ::internal::Triangulation::Iterators<1,spacedim>::line_iterator
2253 {
2254  return typename ::internal::Triangulation::Iterators<1,spacedim>::line_iterator();
2255 }
2256 
2257 
2258 template <int spacedim>
2259 inline
2260 unsigned int
2262 {
2263  Assert(false, ExcImpossibleInDim(0));
2265 }
2266 
2267 
2268 template <int spacedim>
2269 inline
2270 typename ::internal::Triangulation::Iterators<1,spacedim>::quad_iterator
2272 {
2273  return typename ::internal::Triangulation::Iterators<1,spacedim>::quad_iterator();
2274 }
2275 
2276 
2277 
2278 template <int spacedim>
2279 inline
2280 unsigned int
2282 {
2283  Assert(false, ExcImpossibleInDim(0));
2285 }
2286 
2287 
2288 template <int spacedim>
2289 inline
2290 bool
2292 {
2293  return vertex_kind != interior_vertex;
2294 }
2295 
2296 
2297 template <int spacedim>
2298 inline
2301 {
2302  switch (vertex_kind)
2303  {
2304  case left_vertex:
2305  case right_vertex:
2306  {
2307  Assert (tria->vertex_to_boundary_id_map_1d->find (this->vertex_index())
2308  != tria->vertex_to_boundary_id_map_1d->end(),
2309  ExcInternalError());
2310 
2311  return (*tria->vertex_to_boundary_id_map_1d)[this->vertex_index()];
2312  }
2313 
2314  default:
2316  }
2317 
2318 }
2319 
2320 
2321 template <int spacedim>
2322 inline
2324 {
2325  return false;
2326 }
2327 
2328 
2329 
2330 template <int spacedim>
2331 inline
2332 bool TriaAccessor<0, 1, spacedim>::face_flip (const unsigned int face)
2333 {
2334  return false;
2335 }
2336 
2337 
2338 
2339 template <int spacedim>
2340 inline
2341 bool TriaAccessor<0, 1, spacedim>::face_rotation (const unsigned int face)
2342 {
2343  return false;
2344 }
2345 
2346 
2347 
2348 template <int spacedim>
2349 inline
2351 {
2352  return false;
2353 }
2354 
2355 
2356 
2357 template <int spacedim>
2358 inline
2360 {
2361  return false;
2362 }
2363 
2364 
2365 
2366 template <int spacedim>
2367 inline
2369 {
2370  return 0;
2371 }
2372 
2373 
2374 
2375 template <int spacedim>
2376 inline
2378 {
2379  return 0;
2380 }
2381 
2382 
2383 
2384 template <int spacedim>
2385 inline
2387 {
2388  return 0;
2389 }
2390 
2391 
2392 template <int spacedim>
2393 inline
2396 {
2398 }
2399 
2400 
2401 template <int spacedim>
2402 inline
2405 {
2407 }
2408 
2409 
2410 template <int spacedim>
2411 inline
2413 {
2415 }
2416 
2417 template <int spacedim>
2418 inline
2420 {
2421  return -1;
2422 }
2423 
2424 
2425 template <int spacedim>
2426 inline
2428 {
2429  return -1;
2430 }
2431 
2432 
2433 template <int spacedim>
2434 inline
2435 void
2437 {
2438  Assert (tria->vertex_to_boundary_id_map_1d->find (this->vertex_index())
2439  != tria->vertex_to_boundary_id_map_1d->end(),
2440  ExcInternalError());
2441 
2442  (*tria->vertex_to_boundary_id_map_1d)[this->vertex_index()] = b;
2443 }
2444 
2445 
2446 
2447 template <int spacedim>
2448 inline
2450 {
2452 }
2453 
2454 
2455 
2456 template <int spacedim>
2457 inline
2459 {
2460  return tria->vertex_used(global_vertex_index);
2461 }
2462 
2463 /*------------------------ Functions: CellAccessor<dim,spacedim> -----------------------*/
2464 
2465 
2466 template <int dim, int spacedim>
2467 inline
2470  const int level,
2471  const int index,
2472  const AccessorData *local_data)
2473  :
2474  TriaAccessor<dim, dim, spacedim> (parent, level, index, local_data)
2475 {}
2476 
2477 
2478 
2479 template <int dim, int spacedim>
2480 inline
2482  :
2483  TriaAccessor<dim, dim, spacedim> (static_cast<const TriaAccessor<dim, dim, spacedim>&>(cell_accessor))
2484 {}
2485 
2486 
2487 
2488 namespace internal
2489 {
2490  namespace CellAccessor
2491  {
2492  template <int spacedim>
2493  inline
2494  ::TriaIterator<::TriaAccessor<0, 1, spacedim> >
2495  get_face (const ::CellAccessor<1,spacedim> &cell,
2496  const unsigned int i)
2497  {
2499  a (&cell.get_triangulation(),
2500  ((i == 0) && cell.at_boundary(0)
2501  ?
2503  :
2504  ((i == 1) && cell.at_boundary(1)
2505  ?
2507  :
2509  cell.vertex_index(i));
2510  return ::TriaIterator<::TriaAccessor<0, 1, spacedim> >(a);
2511  }
2512 
2513 
2514  template <int spacedim>
2515  inline
2516  ::TriaIterator<::TriaAccessor<1, 2, spacedim> >
2517  get_face (const ::CellAccessor<2,spacedim> &cell,
2518  const unsigned int i)
2519  {
2520  return cell.line(i);
2521  }
2522 
2523 
2524  template <int spacedim>
2525  inline
2526  ::TriaIterator<::TriaAccessor<2, 3, spacedim> >
2527  get_face (const ::CellAccessor<3,spacedim> &cell,
2528  const unsigned int i)
2529  {
2530  return cell.quad(i);
2531  }
2532  }
2533 }
2534 
2535 
2536 
2537 template <int dim, int spacedim>
2538 inline
2539 TriaIterator<TriaAccessor<dim-1, dim, spacedim> >
2540 CellAccessor<dim,spacedim>::face (const unsigned int i) const
2541 {
2542  return ::internal::CellAccessor::get_face (*this, i);
2543 }
2544 
2545 
2546 
2547 template <int dim, int spacedim>
2548 inline
2549 unsigned int
2550 CellAccessor<dim,spacedim>::face_index (const unsigned int i) const
2551 {
2552  switch (dim)
2553  {
2554  case 1:
2555  {
2556  return this->vertex_index(i);
2557  }
2558 
2559  case 2:
2560  return this->line_index(i);
2561 
2562  case 3:
2563  return this->quad_index(i);
2564 
2565  default:
2567  }
2568 }
2569 
2570 
2571 
2572 template <int dim, int spacedim>
2573 inline
2574 int
2576 {
2578  return this->tria->levels[this->present_level]->
2579  neighbors[this->present_index*GeometryInfo<dim>::faces_per_cell+i].second;
2580 }
2581 
2582 
2583 
2584 template <int dim, int spacedim>
2585 inline
2586 int
2588 {
2590  return this->tria->levels[this->present_level]->
2591  neighbors[this->present_index*GeometryInfo<dim>::faces_per_cell+i].first;
2592 }
2593 
2594 
2595 
2596 template <int dim, int spacedim>
2597 inline
2600 {
2601  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
2602  // cells flagged for refinement must be active
2603  // (the @p set_refine_flag function checks this,
2604  // but activity may change when refinement is
2605  // executed and for some reason the refine
2606  // flag is not cleared).
2607  Assert (this->active() || !this->tria->levels[this->present_level]->refine_flags[this->present_index],
2608  ExcRefineCellNotActive());
2609  return RefinementCase<dim>(this->tria->levels[this->present_level]->refine_flags[this->present_index]);
2610 }
2611 
2612 
2613 
2614 template <int dim, int spacedim>
2615 inline
2616 void
2618 {
2619  Assert (this->used() && this->active(), ExcRefineCellNotActive());
2620  Assert (!coarsen_flag_set(),
2621  ExcCellFlaggedForCoarsening());
2622 
2623  this->tria->levels[this->present_level]->refine_flags[this->present_index] = refinement_case;
2624 }
2625 
2626 
2627 
2628 template <int dim, int spacedim>
2629 inline
2630 void
2632 {
2633  Assert (this->used() && this->active(), ExcRefineCellNotActive());
2634  this->tria->levels[this->present_level]->refine_flags[this->present_index] =
2636 }
2637 
2638 
2639 
2640 template <int dim, int spacedim>
2641 inline
2642 bool
2644  const RefinementCase<dim-1> &face_refinement_case) const
2645 {
2646  Assert (dim>1, ExcImpossibleInDim(dim));
2649  Assert (face_refinement_case < RefinementCase<dim>::isotropic_refinement+1,
2650  ExcIndexRange(face_refinement_case,0,RefinementCase<dim>::isotropic_refinement+1));
2651 
2652  // the new refinement case is a combination
2653  // of the minimum required one for the given
2654  // face refinement and the already existing
2655  // flagged refinement case
2656  RefinementCase<dim> old_ref_case = refine_flag_set();
2658  new_ref_case = (old_ref_case
2660  face_no,
2661  this->face_orientation(face_no),
2662  this->face_flip(face_no),
2663  this->face_rotation(face_no)));
2664  set_refine_flag(new_ref_case);
2665  // return, whether we had to change the
2666  // refinement flag
2667  return new_ref_case != old_ref_case;
2668 }
2669 
2670 
2671 
2672 template <int dim, int spacedim>
2673 inline
2674 bool
2676 {
2677  Assert (dim>1, ExcImpossibleInDim(dim));
2680 
2681  // the new refinement case is a combination
2682  // of the minimum required one for the given
2683  // line refinement and the already existing
2684  // flagged refinement case
2685  RefinementCase<dim> old_ref_case=refine_flag_set(),
2686  new_ref_case=old_ref_case
2688  set_refine_flag(new_ref_case);
2689  // return, whether we had to change the
2690  // refinement flag
2691  return new_ref_case != old_ref_case;
2692 }
2693 
2694 
2695 
2696 template <>
2697 inline
2698 ::internal::SubfaceCase<1>
2699 CellAccessor<1>::subface_case(const unsigned int) const
2700 {
2701  return ::internal::SubfaceCase<1>::case_none;
2702 }
2703 
2704 template <>
2705 inline
2706 ::internal::SubfaceCase<1>
2707 CellAccessor<1,2>::subface_case(const unsigned int) const
2708 {
2709  return ::internal::SubfaceCase<1>::case_none;
2710 }
2711 
2712 
2713 template <>
2714 inline
2715 ::internal::SubfaceCase<1>
2716 CellAccessor<1,3>::subface_case(const unsigned int) const
2717 {
2718  return ::internal::SubfaceCase<1>::case_none;
2719 }
2720 
2721 
2722 template <>
2723 inline
2724 ::internal::SubfaceCase<2>
2725 CellAccessor<2>::subface_case(const unsigned int face_no) const
2726 {
2727  Assert(active(), TriaAccessorExceptions::ExcCellNotActive());
2730  return ((face(face_no)->has_children()) ?
2733 }
2734 
2735 template <>
2736 inline
2737 ::internal::SubfaceCase<2>
2738 CellAccessor<2,3>::subface_case(const unsigned int face_no) const
2739 {
2740  Assert(active(), TriaAccessorExceptions::ExcCellNotActive());
2743  return ((face(face_no)->has_children()) ?
2746 }
2747 
2748 
2749 template <>
2750 inline
2751 ::internal::SubfaceCase<3>
2752 CellAccessor<3>::subface_case(const unsigned int face_no) const
2753 {
2754  Assert(active(), TriaAccessorExceptions::ExcCellNotActive());
2757  switch (static_cast<unsigned char> (face(face_no)->refinement_case()))
2758  {
2760  return ::internal::SubfaceCase<3>::case_none;
2761  break;
2763  if (face(face_no)->child(0)->has_children())
2764  {
2765  Assert(face(face_no)->child(0)->refinement_case()==RefinementCase<2>::cut_y,
2766  ExcInternalError());
2767  if (face(face_no)->child(1)->has_children())
2768  {
2769  Assert(face(face_no)->child(1)->refinement_case()==RefinementCase<2>::cut_y,
2770  ExcInternalError());
2771  return ::internal::SubfaceCase<3>::case_x1y2y;
2772  }
2773  else
2774  return ::internal::SubfaceCase<3>::case_x1y;
2775  }
2776  else
2777  {
2778  if (face(face_no)->child(1)->has_children())
2779  {
2780  Assert(face(face_no)->child(1)->refinement_case()==RefinementCase<2>::cut_y,
2781  ExcInternalError());
2782  return ::internal::SubfaceCase<3>::case_x2y;
2783  }
2784  else
2785  return ::internal::SubfaceCase<3>::case_x;
2786  }
2787  break;
2789  if (face(face_no)->child(0)->has_children())
2790  {
2791  Assert(face(face_no)->child(0)->refinement_case()==RefinementCase<2>::cut_x,
2792  ExcInternalError());
2793  if (face(face_no)->child(1)->has_children())
2794  {
2795  Assert(face(face_no)->child(1)->refinement_case()==RefinementCase<2>::cut_x,
2796  ExcInternalError());
2797  return ::internal::SubfaceCase<3>::case_y1x2x;
2798  }
2799  else
2800  return ::internal::SubfaceCase<3>::case_y1x;
2801  }
2802  else
2803  {
2804  if (face(face_no)->child(1)->has_children())
2805  {
2806  Assert(face(face_no)->child(1)->refinement_case()==RefinementCase<2>::cut_x,
2807  ExcInternalError());
2808  return ::internal::SubfaceCase<3>::case_y2x;
2809  }
2810  else
2811  return ::internal::SubfaceCase<3>::case_y;
2812  }
2813  break;
2815  return ::internal::SubfaceCase<3>::case_xy;
2816  break;
2817  default:
2818  Assert(false, ExcInternalError());
2819  }
2820  // we should never get here
2821  return ::internal::SubfaceCase<3>::case_none;
2822 }
2823 
2824 
2825 
2826 template <int dim, int spacedim>
2827 inline
2828 bool
2830 {
2831  Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
2832  // cells flagged for coarsening must be active
2833  // (the @p set_refine_flag function checks this,
2834  // but activity may change when refinement is
2835  // executed and for some reason the refine
2836  // flag is not cleared).
2837  Assert (this->active() || !this->tria->levels[this->present_level]->coarsen_flags[this->present_index],
2838  ExcRefineCellNotActive());
2839  return this->tria->levels[this->present_level]->coarsen_flags[this->present_index];
2840 }
2841 
2842 
2843 
2844 template <int dim, int spacedim>
2845 inline
2846 void
2848 {
2849  Assert (this->used() && this->active(), ExcRefineCellNotActive());
2850  Assert (!refine_flag_set(), ExcCellFlaggedForRefinement());
2851 
2852  this->tria->levels[this->present_level]->coarsen_flags[this->present_index] = true;
2853 }
2854 
2855 
2856 
2857 template <int dim, int spacedim>
2858 inline
2859 void
2861 {
2862  Assert (this->used() && this->active(), ExcRefineCellNotActive());
2863  this->tria->levels[this->present_level]->coarsen_flags[this->present_index] = false;
2864 }
2865 
2866 
2867 
2868 template <int dim, int spacedim>
2869 inline
2871 CellAccessor<dim,spacedim>::neighbor (const unsigned int i) const
2872 {
2874  q (this->tria, neighbor_level (i), neighbor_index (i));
2875 
2876  Assert ((q.state() == IteratorState::past_the_end) || q->used(),
2877  TriaAccessorExceptions::ExcUnusedCellAsNeighbor());
2878 
2879  return q;
2880 }
2881 
2882 
2883 
2884 template <int dim, int spacedim>
2885 inline
2887 CellAccessor<dim,spacedim>::child (const unsigned int i) const
2888 {
2890  q (this->tria, this->present_level+1, this->child_index (i));
2891 
2892  Assert ((q.state() == IteratorState::past_the_end) || q->used(),
2893  TriaAccessorExceptions::ExcUnusedCellAsChild());
2894 
2895  return q;
2896 }
2897 
2898 
2899 
2900 template <int dim, int spacedim>
2901 inline
2902 bool
2904 {
2905  return !this->has_children();
2906 }
2907 
2908 
2909 
2910 template <int dim, int spacedim>
2911 inline
2912 bool
2914 {
2915  Assert (this->active(), ExcMessage("is_locally_owned() only works on active cells!"));
2916 #ifndef DEAL_II_WITH_P4EST
2917  return true;
2918 #else
2919  if (is_artificial())
2920  return false;
2921 
2923  = dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim> *>(this->tria);
2924 
2925  if (pdt == 0)
2926  return true;
2927  else
2928  return (this->subdomain_id() == pdt->locally_owned_subdomain());
2929 #endif
2930 }
2931 
2932 
2933 template <int dim, int spacedim>
2934 inline
2935 bool
2937 {
2938 #ifndef DEAL_II_WITH_P4EST
2939  return true;
2940 #else
2942  = dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim> *>(this->tria);
2943 
2944  if (pdt == 0)
2945  return true;
2946  else
2947  return (this->level_subdomain_id() == pdt->locally_owned_subdomain());
2948 #endif
2949 }
2950 
2951 
2952 template <int dim, int spacedim>
2953 inline
2954 bool
2956 {
2957  Assert (this->active(), ExcMessage("is_ghost() only works on active cells!"));
2958 #ifndef DEAL_II_WITH_P4EST
2959  return false;
2960 #else
2961  if (is_artificial() || this->has_children())
2962  return false;
2963 
2965  = dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim> *>(this->tria);
2966 
2967  if (pdt == 0)
2968  return false;
2969  else
2970  return (this->subdomain_id() != pdt->locally_owned_subdomain());
2971 #endif
2972 }
2973 
2974 
2975 
2976 template <int dim, int spacedim>
2977 inline
2978 bool
2980 {
2981  Assert (this->active(), ExcMessage("is_artificial() only works on active cells!"));
2982 #ifndef DEAL_II_WITH_P4EST
2983  return false;
2984 #else
2985  return this->subdomain_id() == numbers::artificial_subdomain_id;
2986 #endif
2987 }
2988 
2989 
2990 
2991 template <int dim, int spacedim>
2992 inline
2993 unsigned int
2994 CellAccessor<dim,spacedim>::neighbor_face_no (const unsigned int neighbor) const
2995 {
2996  const unsigned int n2=neighbor_of_neighbor_internal(neighbor);
2998  // return this value as the
2999  // neighbor is not coarser
3000  return n2;
3001  else
3002  // the neighbor is coarser
3003  return neighbor_of_coarser_neighbor(neighbor).first;
3004 }
3005 
3006 
3007 template <int dim, int spacedim>
3008 inline
3009 bool
3011 {
3012  Assert (this->tria == other.tria, TriaAccessorExceptions::ExcCantCompareIterators());
3013 
3014  if (level_subdomain_id() != other.level_subdomain_id())
3015  return (level_subdomain_id() < other.level_subdomain_id());
3016 
3017  if (active() && other.active() &&
3018  (subdomain_id() != other.subdomain_id()))
3019  return (subdomain_id() < other.subdomain_id());
3020 
3022 }
3023 
3024 
3025 
3026 DEAL_II_NAMESPACE_CLOSE
3027 
3028 #endif
bool flag_for_line_refinement(const unsigned int line_no) const
types::boundary_id boundary_indicator() const
void set_face_orientation(const unsigned int face, const bool orientation) const
void clear_user_flag() const
const Triangulation< dim, spacedim > & get_triangulation() const
static const unsigned int invalid_unsigned_int
Definition: types.h:191
unsigned int user_index() const
bool face_flip(const unsigned int face) const
static unsigned int standard_to_real_face_line(const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
Iterator points to a valid object.
void recursively_clear_user_index() const
bool operator<(const TriaAccessorBase &other) const
TriaIterator< CellAccessor< dim, spacedim > > child(const unsigned int i) const
void set_refinement_case(const RefinementCase< structdim > &ref_case) const
bool user_flag_set() const
void set_user_pointer(void *p) const
static unsigned int vertex_index(const TriaAccessor< 1, dim, spacedim > &accessor, const unsigned int corner)
TriaAccessorBase & operator=(const TriaAccessorBase &)
::ExceptionBase & ExcMessage(std::string arg1)
#define AssertIndexRange(index, range)
Definition: exceptions.h:888
unsigned int neighbor_face_no(const unsigned int neighbor) const
TriaIterator< TriaAccessor< dim-1, dim, spacedim > > face(const unsigned int i) const
void copy_from(const TriaAccessorBase &)
::internal::SubfaceCase< dim > subface_case(const unsigned int face_no) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor(const unsigned int i) const
void recursively_set_user_index(const unsigned int p) const
void recursively_clear_user_flag() const
bool face_rotation(const unsigned int face) const
int neighbor_level(const unsigned int i) const
bool face_orientation(const unsigned int face) const
types::subdomain_id level_subdomain_id() const
bool is_translation_of(const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &o) const
void clear_user_data() const
unsigned int vertex_index(const unsigned int i) const
void copy_from(const InvalidAccessor &)
bool at_boundary() const
void clear_user_index() const
InvalidAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
bool is_locally_owned() const
void set_used_flag() const
static void set_face_rotation(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int, const bool)
unsigned int line_index(const unsigned int i) const
Iterator reached end of container.
static unsigned int quad_index(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int)
static void set_face_orientation(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int, const bool)
void clear_refinement_case() const
typename::internal::Triangulation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
Iterator is invalid, probably due to an error.
unsigned int n_children() const
typename::internal::Triangulation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
int child_index(const unsigned int i) const
static bool face_rotation(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int)
static bool face_flip(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int)
VectorizedArray< Number > min(const ::VectorizedArray< Number > &x, const ::VectorizedArray< Number > &y)
void clear_children() const
void recursively_set_user_pointer(void *p) const
#define Assert(cond, exc)
Definition: exceptions.h:299
void set_boundary_indicator(const types::boundary_id) const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child(const unsigned int i) const
void clear_used_flag() const
bool line_orientation(const unsigned int line) const
void set_children(const unsigned int i, const int index) const
types::subdomain_id locally_owned_subdomain() const
real_type norm_square() const
const Triangulation< dim, spacedim > * tria
bool is_locally_owned_on_level() const
::internal::Triangulation::TriaObjects<::internal::Triangulation::TriaObject< structdim > > & objects() const
bool is_artificial() const
const Boundary< dim, spacedim > & get_boundary() const
bool coarsen_flag_set() const
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void set_all_boundary_indicators(const types::boundary_id) const
void set_face_flip(const unsigned int face, const bool flip) const
RefinementCase< structdim > refinement_case() const
typename::internal::TriaAccessor::PresentLevelType< structdim, dim >::type present_level
bool operator<(const CellAccessor< dim, spacedim > &other) const
unsigned int number_of_children() const
::ExceptionBase & ExcImpossibleInDim(int arg1)
void clear_refine_flag() const
void set_refine_flag(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
bool operator==(const TriaAccessorBase &) const
unsigned int face_index(const unsigned int i) const
static void set_line_orientation(const TriaAccessor< 1, dim, spacedim > &, const unsigned int, const bool)
bool operator==(const InvalidAccessor &) const
void clear_coarsen_flag() const
void set_face_rotation(const unsigned int face, const bool rotation) const
int neighbor_index(const unsigned int i) const
static RefinementCase< dim > min_cell_refinement_case_for_line_refinement(const unsigned int line_no)
RefinementCase< dim > refine_flag_set() const
void set_user_flag() const
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
const types::subdomain_id artificial_subdomain_id
Definition: types.h:286
int isotropic_child_index(const unsigned int i) const
IteratorState::IteratorStates state() const
static bool line_orientation(const TriaAccessor< 1, dim, spacedim > &, const unsigned int)
static void set_face_flip(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int, const bool)
bool flag_for_face_refinement(const unsigned int face_no, const RefinementCase< dim-1 > &face_refinement_case=RefinementCase< dim-1 >::isotropic_refinement) const
void clear_user_pointer() const
Point< spacedim > & vertex(const unsigned int i) const
unsigned int quad_index(const unsigned int i) const
IteratorState::IteratorStates state() const
void set_coarsen_flag() const
CellAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
void set_parent(const unsigned int parent_index)
unsigned int max_refinement_depth() const
void recursively_clear_user_pointer() const
void set_line_orientation(const unsigned int line, const bool orientation) const
::ExceptionBase & ExcNotImplemented()
TriaAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
double diameter() const
Point< spacedim > center() const
unsigned char boundary_id
Definition: types.h:128
types::subdomain_id subdomain_id() const
const types::boundary_id internal_face_boundary_id
Definition: types.h:253
void * user_pointer() const
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *=0)
double minimum_vertex_distance() const
void recursively_set_user_flag() const
::ExceptionBase & ExcInternalError()
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
static bool face_orientation(const TriaAccessor< structdim, dim, spacedim > &, const unsigned int)
void set_user_index(const unsigned int p) const
static RefinementCase< dim > min_cell_refinement_case_for_face_refinement(const RefinementCase< dim-1 > &face_refinement_case, const unsigned int face_no, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
bool has_children() const
static unsigned int line_index(const TriaAccessor< 1, dim, spacedim > &, const unsigned int)
bool operator!=(const TriaAccessorBase &) const