Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
geometry_info.h
1 // ---------------------------------------------------------------------
2 // @f$Id: geometry_info.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 1998 - 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__geometry_info_h
18 #define __deal2__geometry_info_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/point.h>
24 
26 
27 
28 
42 template <int dim>
44 {
80  {
81  no_refinement= 0,
82 
83  isotropic_refinement = static_cast<unsigned char>(-1)
84  };
85 };
86 
87 
88 
99 template <>
101 {
137  {
138  no_refinement= 0,
139  cut_x = 1,
140 
141  isotropic_refinement = cut_x
142  };
143 };
144 
145 
146 
158 template <>
160 {
196  {
197  no_refinement= 0,
198  cut_x = 1,
199  cut_y = 2,
200  cut_xy = cut_x | cut_y,
201 
202  isotropic_refinement = cut_xy
203  };
204 };
205 
206 
207 
220 template <>
222 {
258  {
259  no_refinement= 0,
260  cut_x = 1,
261  cut_y = 2,
262  cut_xy = cut_x | cut_y,
263  cut_z = 4,
264  cut_xz = cut_x | cut_z,
265  cut_yz = cut_y | cut_z,
266  cut_xyz = cut_x | cut_y | cut_z,
267 
268  isotropic_refinement = cut_xyz
269  };
270 };
271 
272 
273 
284 template <int dim>
286 {
287 public:
292  RefinementCase ();
293 
301  RefinementCase (const typename RefinementPossibilities<dim>::Possibilities refinement_case);
302 
312  explicit RefinementCase (const unsigned char refinement_case);
313 
339  operator unsigned char () const;
340 
347  RefinementCase operator | (const RefinementCase &r) const;
348 
355  RefinementCase operator & (const RefinementCase &r) const;
356 
372  RefinementCase operator ~ () const;
373 
374 
384  static
385  RefinementCase cut_axis (const unsigned int i);
386 
392  static std::size_t memory_consumption ();
393 
398  template <class Archive>
399  void serialize(Archive &ar,
400  const unsigned int version);
401 
405  DeclException1 (ExcInvalidRefinementCase,
406  int,
407  << "The refinement flags given (" << arg1 << ") contain set bits that do not "
408  << "make sense for the space dimension of the object to which they are applied.");
409 
410 private:
417 unsigned char value :
418  (dim > 0 ? dim : 1);
419 };
420 
421 
422 namespace internal
423 {
424 
425 
447  template <int dim>
449  {
456  {
457  case_none = 0,
458 
459  case_isotropic = static_cast<unsigned char>(-1)
460  };
461  };
462 
463 
474  template <>
476  {
485  {
486  case_none = 0,
487 
488  case_isotropic = case_none
489  };
490  };
491 
492 
493 
505  template <>
507  {
518  {
519  case_none = 0,
520 
521  case_isotropic = case_none
522  };
523  };
524 
525 
526 
539  template <>
541  {
555  {
556  case_none = 0,
557  case_x = 1,
558 
559  case_isotropic = case_x
560  };
561  };
562 
563 
564 
660  template <>
662  {
674  {
675  case_none = 0,
676  case_x = 1,
677  case_x1y = 2,
678  case_x2y = 3,
679  case_x1y2y = 4,
680  case_y = 5,
681  case_y1x = 6,
682  case_y2x = 7,
683  case_y1x2x = 8,
684  case_xy = 9,
685 
686  case_isotropic = case_xy
687  };
688  };
689 
690 
691 
692 
701  template <int dim>
702  class SubfaceCase : public SubfacePossibilities<dim>
703  {
704  public:
713  SubfaceCase (const typename SubfacePossibilities<dim>::Possibilities subface_possibility);
714 
740  operator unsigned char () const;
741 
747  static std::size_t memory_consumption ();
748 
752  DeclException1 (ExcInvalidSubfaceCase,
753  int,
754  << "The subface case given (" << arg1 << ") does not make sense "
755  << "for the space dimension of the object to which they are applied.");
756 
757  private:
764  unsigned char value :
765  (dim == 3 ? 4 : 1);
766  };
767 
768 } // namespace internal
769 
770 
771 
772 template <int dim> struct GeometryInfo;
773 
774 
775 
776 
797 template <>
798 struct GeometryInfo<0>
799 {
800 
812  static const unsigned int max_children_per_cell = 1;
813 
817  static const unsigned int faces_per_cell = 0;
818 
830  static const unsigned int max_children_per_face = 0;
831 
840  static unsigned int n_children(const RefinementCase<0> &refinement_case);
841 
845  static const unsigned int vertices_per_cell = 1;
846 
856  static const unsigned int vertices_per_face = 0;
857 
861  static const unsigned int lines_per_face = 0;
862 
866  static const unsigned int quads_per_face = 0;
867 
871  static const unsigned int lines_per_cell = 0;
872 
877  static const unsigned int quads_per_cell = 0;
878 
883  static const unsigned int hexes_per_cell = 0;
884 };
885 
886 
887 
888 
889 
1422 template <int dim>
1423 struct GeometryInfo
1424 {
1425 
1437  static const unsigned int max_children_per_cell = 1 << dim;
1438 
1442  static const unsigned int faces_per_cell = 2 * dim;
1443 
1455  static const unsigned int max_children_per_face = GeometryInfo<dim-1>::max_children_per_cell;
1456 
1460  static const unsigned int vertices_per_cell = 1 << dim;
1461 
1466  static const unsigned int vertices_per_face = GeometryInfo<dim-1>::vertices_per_cell;
1467 
1471  static const unsigned int lines_per_face
1472  = GeometryInfo<dim-1>::lines_per_cell;
1473 
1477  static const unsigned int quads_per_face
1478  = GeometryInfo<dim-1>::quads_per_cell;
1479 
1492  static const unsigned int lines_per_cell
1493  = (2*GeometryInfo<dim-1>::lines_per_cell +
1494  GeometryInfo<dim-1>::vertices_per_cell);
1495 
1506  static const unsigned int quads_per_cell
1507  = (2*GeometryInfo<dim-1>::quads_per_cell +
1508  GeometryInfo<dim-1>::lines_per_cell);
1509 
1514  static const unsigned int hexes_per_cell
1515  = (2*GeometryInfo<dim-1>::hexes_per_cell +
1516  GeometryInfo<dim-1>::quads_per_cell);
1517 
1544  static const unsigned int ucd_to_deal[vertices_per_cell];
1545 
1565  static const unsigned int dx_to_deal[vertices_per_cell];
1566 
1584  static const unsigned int vertex_to_face[vertices_per_cell][dim];
1585 
1591  static
1592  unsigned int
1593  n_children(const RefinementCase<dim> &refinement_case);
1594 
1601  static
1602  unsigned int
1603  n_subfaces(const internal::SubfaceCase<dim> &subface_case);
1604 
1620  static
1621  double
1622  subface_ratio(const internal::SubfaceCase<dim> &subface_case,
1623  const unsigned int subface_no);
1624 
1633  static
1634  RefinementCase<dim-1>
1635  face_refinement_case (const RefinementCase<dim> &cell_refinement_case,
1636  const unsigned int face_no,
1637  const bool face_orientation = true,
1638  const bool face_flip = false,
1639  const bool face_rotation = false);
1640 
1649  static
1651  min_cell_refinement_case_for_face_refinement
1652  (const RefinementCase<dim-1> &face_refinement_case,
1653  const unsigned int face_no,
1654  const bool face_orientation = true,
1655  const bool face_flip = false,
1656  const bool face_rotation = false);
1657 
1665  static
1667  line_refinement_case(const RefinementCase<dim> &cell_refinement_case,
1668  const unsigned int line_no);
1669 
1676  static
1678  min_cell_refinement_case_for_line_refinement(const unsigned int line_no);
1679 
1745  static
1746  unsigned int
1747  child_cell_on_face (const RefinementCase<dim> &ref_case,
1748  const unsigned int face,
1749  const unsigned int subface,
1750  const bool face_orientation = true,
1751  const bool face_flip = false,
1752  const bool face_rotation = false,
1753  const RefinementCase<dim-1> &face_refinement_case
1755 
1777  static
1778  unsigned int
1779  line_to_cell_vertices (const unsigned int line,
1780  const unsigned int vertex);
1781 
1816  static
1817  unsigned int
1818  face_to_cell_vertices (const unsigned int face,
1819  const unsigned int vertex,
1820  const bool face_orientation = true,
1821  const bool face_flip = false,
1822  const bool face_rotation = false);
1823 
1843  static
1844  unsigned int
1845  face_to_cell_lines (const unsigned int face,
1846  const unsigned int line,
1847  const bool face_orientation = true,
1848  const bool face_flip = false,
1849  const bool face_rotation = false);
1850 
1864  static
1865  unsigned int
1866  standard_to_real_face_vertex (const unsigned int vertex,
1867  const bool face_orientation = true,
1868  const bool face_flip = false,
1869  const bool face_rotation = false);
1870 
1884  static
1885  unsigned int
1886  real_to_standard_face_vertex (const unsigned int vertex,
1887  const bool face_orientation = true,
1888  const bool face_flip = false,
1889  const bool face_rotation = false);
1890 
1904  static
1905  unsigned int
1906  standard_to_real_face_line (const unsigned int line,
1907  const bool face_orientation = true,
1908  const bool face_flip = false,
1909  const bool face_rotation = false);
1910 
1924  static
1925  unsigned int
1926  real_to_standard_face_line (const unsigned int line,
1927  const bool face_orientation = true,
1928  const bool face_flip = false,
1929  const bool face_rotation = false);
1930 
1938  static
1939  Point<dim>
1940  unit_cell_vertex (const unsigned int vertex);
1941 
1957  static
1958  unsigned int
1959  child_cell_from_point (const Point<dim> &p);
1960 
1975  static
1976  Point<dim>
1977  cell_to_child_coordinates (const Point<dim> &p,
1978  const unsigned int child_index,
1979  const RefinementCase<dim> refine_case
1981 
1990  static
1991  Point<dim>
1992  child_to_cell_coordinates (const Point<dim> &p,
1993  const unsigned int child_index,
1994  const RefinementCase<dim> refine_case
1996 
2002  static
2003  bool
2004  is_inside_unit_cell (const Point<dim> &p);
2005 
2028  static
2029  bool
2030  is_inside_unit_cell (const Point<dim> &p,
2031  const double eps);
2032 
2039  static
2040  Point<dim>
2041  project_to_unit_cell (const Point<dim> &p);
2042 
2051  static
2052  double
2053  distance_to_unit_cell (const Point<dim> &p);
2054 
2060  static
2061  double
2062  d_linear_shape_function (const Point<dim> &xi,
2063  const unsigned int i);
2064 
2070  static
2072  d_linear_shape_function_gradient (const Point<dim> &xi,
2073  const unsigned int i);
2074 
2171  template <int spacedim>
2172  static void
2173  alternating_form_at_vertices
2174 #ifndef DEAL_II_CONSTEXPR_BUG
2175  (const Point<spacedim> (&vertices)[vertices_per_cell],
2176  Tensor<spacedim-dim,spacedim> (&forms)[vertices_per_cell]);
2177 #else
2178  (const Point<spacedim> *vertices,
2179  Tensor<spacedim-dim,spacedim> *forms);
2180 #endif
2181 
2199  static const unsigned int unit_normal_direction[faces_per_cell];
2200 
2232  static const int unit_normal_orientation[faces_per_cell];
2233 
2241  static const unsigned int opposite_face[faces_per_cell];
2242 
2243 
2247  DeclException1 (ExcInvalidCoordinate,
2248  double,
2249  << "The coordinates must satisfy 0 <= x_i <= 1, "
2250  << "but here we have x_i=" << arg1);
2251 
2255  DeclException3 (ExcInvalidSubface,
2256  int, int, int,
2257  << "RefinementCase<dim> " << arg1 << ": face " << arg2
2258  << " has no subface " << arg3);
2259 };
2260 
2261 
2262 
2263 
2264 #ifndef DOXYGEN
2265 
2266 
2267 /* -------------- declaration of explicit specializations ------------- */
2268 
2269 #ifndef DEAL_II_MEMBER_ARRAY_SPECIALIZATION_BUG
2270 template <>
2271 const unsigned int GeometryInfo<1>::unit_normal_direction[faces_per_cell];
2272 template <>
2273 const unsigned int GeometryInfo<2>::unit_normal_direction[faces_per_cell];
2274 template <>
2275 const unsigned int GeometryInfo<3>::unit_normal_direction[faces_per_cell];
2276 template <>
2277 const unsigned int GeometryInfo<4>::unit_normal_direction[faces_per_cell];
2278 
2279 template <>
2280 const int GeometryInfo<1>::unit_normal_orientation[faces_per_cell];
2281 template <>
2282 const int GeometryInfo<2>::unit_normal_orientation[faces_per_cell];
2283 template <>
2284 const int GeometryInfo<3>::unit_normal_orientation[faces_per_cell];
2285 template <>
2286 const int GeometryInfo<4>::unit_normal_orientation[faces_per_cell];
2287 
2288 template <>
2289 const unsigned int GeometryInfo<1>::opposite_face[faces_per_cell];
2290 template <>
2291 const unsigned int GeometryInfo<2>::opposite_face[faces_per_cell];
2292 template <>
2293 const unsigned int GeometryInfo<3>::opposite_face[faces_per_cell];
2294 template <>
2295 const unsigned int GeometryInfo<4>::opposite_face[faces_per_cell];
2296 #endif
2297 
2298 
2299 template <>
2303  const unsigned int i);
2304 template <>
2308  const unsigned int i);
2309 template <>
2313  const unsigned int i);
2314 
2315 
2316 
2317 
2318 /* -------------- inline functions ------------- */
2319 
2320 namespace internal
2321 {
2322 
2323  template <int dim>
2324  inline
2325  SubfaceCase<dim>::SubfaceCase (const typename SubfacePossibilities<dim>::Possibilities subface_possibility)
2326  :
2327  value (subface_possibility)
2328  {}
2329 
2330 
2331  template <int dim>
2332  inline
2333  SubfaceCase<dim>::operator unsigned char () const
2334  {
2335  return value;
2336  }
2337 
2338 
2339 } // namespace internal
2340 
2341 
2342 template <int dim>
2343 inline
2345 RefinementCase<dim>::cut_axis (const unsigned int)
2346 {
2347  Assert (false, ExcInternalError());
2348  return static_cast<unsigned char>(-1);
2349 }
2350 
2351 
2352 template <>
2353 inline
2355 RefinementCase<1>::cut_axis (const unsigned int i)
2356 {
2357  const unsigned int dim = 1;
2358  Assert (i < dim, ExcIndexRange(i, 0, dim));
2359 
2360  static const RefinementCase options[dim] = { RefinementPossibilities<1>::cut_x };
2361  return options[i];
2362 }
2363 
2364 
2365 
2366 template <>
2367 inline
2369 RefinementCase<2>::cut_axis (const unsigned int i)
2370 {
2371  const unsigned int dim = 2;
2372  Assert (i < dim, ExcIndexRange(i, 0, dim));
2373 
2374  static const RefinementCase options[dim] = { RefinementPossibilities<2>::cut_x,
2376  };
2377  return options[i];
2378 }
2379 
2380 
2381 
2382 template <>
2383 inline
2385 RefinementCase<3>::cut_axis (const unsigned int i)
2386 {
2387  const unsigned int dim = 3;
2388  Assert (i < dim, ExcIndexRange(i, 0, dim));
2389 
2390  static const RefinementCase options[dim] = { RefinementPossibilities<3>::cut_x,
2393  };
2394  return options[i];
2395 }
2396 
2397 
2398 
2399 template <int dim>
2400 inline
2402  :
2403  value (RefinementPossibilities<dim>::no_refinement)
2404 {}
2405 
2406 
2407 
2408 template <int dim>
2409 inline
2411 RefinementCase (const typename RefinementPossibilities<dim>::Possibilities refinement_case)
2412  :
2413  value (refinement_case)
2414 {
2415  // check that only those bits of
2416  // the given argument are set that
2417  // make sense for a given space
2418  // dimension
2420  refinement_case,
2421  ExcInvalidRefinementCase (refinement_case));
2422 }
2423 
2424 
2425 
2426 template <int dim>
2427 inline
2428 RefinementCase<dim>::RefinementCase (const unsigned char refinement_case)
2429  :
2430  value (refinement_case)
2431 {
2432  // check that only those bits of
2433  // the given argument are set that
2434  // make sense for a given space
2435  // dimension
2437  refinement_case,
2438  ExcInvalidRefinementCase (refinement_case));
2439 }
2440 
2441 
2442 
2443 template <int dim>
2444 inline
2445 RefinementCase<dim>::operator unsigned char () const
2446 {
2447  return value;
2448 }
2449 
2450 
2451 
2452 template <int dim>
2453 inline
2456 {
2457  return RefinementCase<dim>(static_cast<unsigned char> (value | r.value));
2458 }
2459 
2460 
2461 
2462 template <int dim>
2463 inline
2466 {
2467  return RefinementCase<dim>(static_cast<unsigned char> (value & r.value));
2468 }
2469 
2470 
2471 
2472 template <int dim>
2473 inline
2476 {
2477  return RefinementCase<dim>(static_cast<unsigned char> (
2479 }
2480 
2481 
2482 
2483 
2484 template <int dim>
2485 inline
2486 std::size_t
2488 {
2489  return sizeof(RefinementCase<dim>);
2490 }
2491 
2492 
2493 
2494 template <int dim>
2495 template <class Archive>
2496 void RefinementCase<dim>::serialize (Archive &ar,
2497  const unsigned int)
2498 {
2499  // serialization can't deal with bitfields, so copy from/to a full sized
2500  // unsigned char
2501  unsigned char uchar_value = value;
2502  ar &uchar_value;
2503  value = uchar_value;
2504 }
2505 
2506 
2507 
2508 
2509 template <>
2510 inline
2511 Point<1>
2512 GeometryInfo<1>::unit_cell_vertex (const unsigned int vertex)
2513 {
2514  Assert (vertex < vertices_per_cell,
2515  ExcIndexRange (vertex, 0, vertices_per_cell));
2516 
2517  return Point<1>(static_cast<double>(vertex));
2518 }
2519 
2520 
2521 
2522 template <>
2523 inline
2524 Point<2>
2525 GeometryInfo<2>::unit_cell_vertex (const unsigned int vertex)
2526 {
2527  Assert (vertex < vertices_per_cell,
2528  ExcIndexRange (vertex, 0, vertices_per_cell));
2529 
2530  return Point<2>(vertex%2, vertex/2);
2531 }
2532 
2533 
2534 
2535 template <>
2536 inline
2537 Point<3>
2538 GeometryInfo<3>::unit_cell_vertex (const unsigned int vertex)
2539 {
2540  Assert (vertex < vertices_per_cell,
2541  ExcIndexRange (vertex, 0, vertices_per_cell));
2542 
2543  return Point<3>(vertex%2, vertex/2%2, vertex/4);
2544 }
2545 
2546 
2547 
2548 template <int dim>
2549 inline
2550 Point<dim>
2551 GeometryInfo<dim>::unit_cell_vertex (const unsigned int)
2552 {
2553  Assert(false, ExcNotImplemented());
2554 
2555  return Point<dim> ();
2556 }
2557 
2558 
2559 
2560 template <>
2561 inline
2562 unsigned int
2564 {
2565  Assert ((p[0] >= 0) && (p[0] <= 1), ExcInvalidCoordinate(p[0]));
2566 
2567  return (p[0] <= 0.5 ? 0 : 1);
2568 }
2569 
2570 
2571 
2572 template <>
2573 inline
2574 unsigned int
2576 {
2577  Assert ((p[0] >= 0) && (p[0] <= 1), ExcInvalidCoordinate(p[0]));
2578  Assert ((p[1] >= 0) && (p[1] <= 1), ExcInvalidCoordinate(p[1]));
2579 
2580  return (p[0] <= 0.5 ?
2581  (p[1] <= 0.5 ? 0 : 2) :
2582  (p[1] <= 0.5 ? 1 : 3));
2583 }
2584 
2585 
2586 
2587 template <>
2588 inline
2589 unsigned int
2591 {
2592  Assert ((p[0] >= 0) && (p[0] <= 1), ExcInvalidCoordinate(p[0]));
2593  Assert ((p[1] >= 0) && (p[1] <= 1), ExcInvalidCoordinate(p[1]));
2594  Assert ((p[2] >= 0) && (p[2] <= 1), ExcInvalidCoordinate(p[2]));
2595 
2596  return (p[0] <= 0.5 ?
2597  (p[1] <= 0.5 ?
2598  (p[2] <= 0.5 ? 0 : 4) :
2599  (p[2] <= 0.5 ? 2 : 6)) :
2600  (p[1] <= 0.5 ?
2601  (p[2] <= 0.5 ? 1 : 5) :
2602  (p[2] <= 0.5 ? 3 : 7)));
2603 }
2604 
2605 
2606 template <int dim>
2607 inline
2608 unsigned int
2610 {
2611  Assert(false, ExcNotImplemented());
2612 
2613  return 0;
2614 }
2615 
2616 
2617 
2618 template <>
2619 inline
2620 Point<1>
2622  const unsigned int child_index,
2623  const RefinementCase<1> refine_case)
2624 
2625 {
2626  Assert (child_index < 2,
2627  ExcIndexRange (child_index, 0, 2));
2628  Assert (refine_case==RefinementCase<1>::cut_x,
2629  ExcInternalError());
2630  (void)refine_case; // removes -Wunused-parameter warning in optimized mode
2631 
2632  return p*2.0-unit_cell_vertex(child_index);
2633 }
2634 
2635 
2636 
2637 template <>
2638 inline
2639 Point<2>
2641  const unsigned int child_index,
2642  const RefinementCase<2> refine_case)
2643 
2644 {
2645  Assert (child_index < GeometryInfo<2>::n_children(refine_case),
2646  ExcIndexRange (child_index, 0, GeometryInfo<2>::n_children(refine_case)));
2647 
2648  Point<2> point=p;
2649  switch (refine_case)
2650  {
2652  point[0]*=2.0;
2653  if (child_index==1)
2654  point[0]-=1.0;
2655  break;
2657  point[1]*=2.0;
2658  if (child_index==1)
2659  point[1]-=1.0;
2660  break;
2662  point*=2.0;
2663  point-=unit_cell_vertex(child_index);
2664  break;
2665  default:
2666  Assert(false, ExcInternalError());
2667  }
2668 
2669  return point;
2670 }
2671 
2672 
2673 
2674 template <>
2675 inline
2676 Point<3>
2678  const unsigned int child_index,
2679  const RefinementCase<3> refine_case)
2680 
2681 {
2682  Assert (child_index < GeometryInfo<3>::n_children(refine_case),
2683  ExcIndexRange (child_index, 0, GeometryInfo<3>::n_children(refine_case)));
2684 
2685  Point<3> point=p;
2686  // there might be a cleverer way to do
2687  // this, but since this function is called
2688  // in very few places for initialization
2689  // purposes only, I don't care at the
2690  // moment
2691  switch (refine_case)
2692  {
2694  point[0]*=2.0;
2695  if (child_index==1)
2696  point[0]-=1.0;
2697  break;
2699  point[1]*=2.0;
2700  if (child_index==1)
2701  point[1]-=1.0;
2702  break;
2704  point[2]*=2.0;
2705  if (child_index==1)
2706  point[2]-=1.0;
2707  break;
2709  point[0]*=2.0;
2710  point[1]*=2.0;
2711  if (child_index%2==1)
2712  point[0]-=1.0;
2713  if (child_index/2==1)
2714  point[1]-=1.0;
2715  break;
2717  // careful, this is slightly
2718  // different from xy and yz due to
2719  // different internal numbering of
2720  // children!
2721  point[0]*=2.0;
2722  point[2]*=2.0;
2723  if (child_index/2==1)
2724  point[0]-=1.0;
2725  if (child_index%2==1)
2726  point[2]-=1.0;
2727  break;
2729  point[1]*=2.0;
2730  point[2]*=2.0;
2731  if (child_index%2==1)
2732  point[1]-=1.0;
2733  if (child_index/2==1)
2734  point[2]-=1.0;
2735  break;
2737  point*=2.0;
2738  point-=unit_cell_vertex(child_index);
2739  break;
2740  default:
2741  Assert(false, ExcInternalError());
2742  }
2743 
2744  return point;
2745 }
2746 
2747 
2748 
2749 template <int dim>
2750 inline
2751 Point<dim>
2753  const unsigned int /*child_index*/,
2754  const RefinementCase<dim> /*refine_case*/)
2755 
2756 {
2757  AssertThrow (false, ExcNotImplemented());
2758  return Point<dim>();
2759 }
2760 
2761 
2762 
2763 template <>
2764 inline
2765 Point<1>
2767  const unsigned int child_index,
2768  const RefinementCase<1> refine_case)
2769 
2770 {
2771  Assert (child_index < 2,
2772  ExcIndexRange (child_index, 0, 2));
2773  Assert (refine_case==RefinementCase<1>::cut_x,
2774  ExcInternalError());
2775  (void)refine_case; // removes -Wunused-parameter warning in optimized mode
2776 
2777  return (p+unit_cell_vertex(child_index))*0.5;
2778 }
2779 
2780 
2781 
2782 template <>
2783 inline
2784 Point<3>
2786  const unsigned int child_index,
2787  const RefinementCase<3> refine_case)
2788 
2789 {
2790  Assert (child_index < GeometryInfo<3>::n_children(refine_case),
2791  ExcIndexRange (child_index, 0, GeometryInfo<3>::n_children(refine_case)));
2792 
2793  Point<3> point=p;
2794  // there might be a cleverer way to do
2795  // this, but since this function is called
2796  // in very few places for initialization
2797  // purposes only, I don't care at the
2798  // moment
2799  switch (refine_case)
2800  {
2802  if (child_index==1)
2803  point[0]+=1.0;
2804  point[0]*=0.5;
2805  break;
2807  if (child_index==1)
2808  point[1]+=1.0;
2809  point[1]*=0.5;
2810  break;
2812  if (child_index==1)
2813  point[2]+=1.0;
2814  point[2]*=0.5;
2815  break;
2817  if (child_index%2==1)
2818  point[0]+=1.0;
2819  if (child_index/2==1)
2820  point[1]+=1.0;
2821  point[0]*=0.5;
2822  point[1]*=0.5;
2823  break;
2825  // careful, this is slightly
2826  // different from xy and yz due to
2827  // different internal numbering of
2828  // children!
2829  if (child_index/2==1)
2830  point[0]+=1.0;
2831  if (child_index%2==1)
2832  point[2]+=1.0;
2833  point[0]*=0.5;
2834  point[2]*=0.5;
2835  break;
2837  if (child_index%2==1)
2838  point[1]+=1.0;
2839  if (child_index/2==1)
2840  point[2]+=1.0;
2841  point[1]*=0.5;
2842  point[2]*=0.5;
2843  break;
2845  point+=unit_cell_vertex(child_index);
2846  point*=0.5;
2847  break;
2848  default:
2849  Assert(false, ExcInternalError());
2850  }
2851 
2852  return point;
2853 }
2854 
2855 
2856 
2857 template <>
2858 inline
2859 Point<2>
2861  const unsigned int child_index,
2862  const RefinementCase<2> refine_case)
2863 {
2864  Assert (child_index < GeometryInfo<2>::n_children(refine_case),
2865  ExcIndexRange (child_index, 0, GeometryInfo<2>::n_children(refine_case)));
2866 
2867  Point<2> point=p;
2868  switch (refine_case)
2869  {
2871  if (child_index==1)
2872  point[0]+=1.0;
2873  point[0]*=0.5;
2874  break;
2876  if (child_index==1)
2877  point[1]+=1.0;
2878  point[1]*=0.5;
2879  break;
2881  point+=unit_cell_vertex(child_index);
2882  point*=0.5;
2883  break;
2884  default:
2885  Assert(false, ExcInternalError());
2886  }
2887 
2888  return point;
2889 }
2890 
2891 
2892 
2893 template <int dim>
2894 inline
2895 Point<dim>
2897  const unsigned int /*child_index*/,
2898  const RefinementCase<dim> /*refine_case*/)
2899 {
2900  AssertThrow (false, ExcNotImplemented());
2901  return Point<dim>();
2902 }
2903 
2904 
2905 
2906 template <>
2907 inline
2908 bool
2910 {
2911  return (p[0] >= 0.) && (p[0] <= 1.);
2912 }
2913 
2914 
2915 
2916 template <>
2917 inline
2918 bool
2920 {
2921  return (p[0] >= 0.) && (p[0] <= 1.) &&
2922  (p[1] >= 0.) && (p[1] <= 1.);
2923 }
2924 
2925 
2926 
2927 template <>
2928 inline
2929 bool
2931 {
2932  return (p[0] >= 0.) && (p[0] <= 1.) &&
2933  (p[1] >= 0.) && (p[1] <= 1.) &&
2934  (p[2] >= 0.) && (p[2] <= 1.);
2935 }
2936 
2937 template <>
2938 inline
2939 bool
2941  const double eps)
2942 {
2943  return (p[0] >= -eps) && (p[0] <= 1.+eps);
2944 }
2945 
2946 
2947 
2948 template <>
2949 inline
2950 bool
2952  const double eps)
2953 {
2954  const double l = -eps, u = 1+eps;
2955  return (p[0] >= l) && (p[0] <= u) &&
2956  (p[1] >= l) && (p[1] <= u);
2957 }
2958 
2959 
2960 
2961 template <>
2962 inline
2963 bool
2965  const double eps)
2966 {
2967  const double l = -eps, u = 1.0+eps;
2968  return (p[0] >= l) && (p[0] <= u) &&
2969  (p[1] >= l) && (p[1] <= u) &&
2970  (p[2] >= l) && (p[2] <= u);
2971 }
2972 
2973 
2974 #endif // DOXYGEN
2975 DEAL_II_NAMESPACE_CLOSE
2976 
2977 #endif
static Point< dim > child_to_cell_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)
UpdateFlags operator&(UpdateFlags f1, UpdateFlags f2)
static unsigned int child_cell_from_point(const Point< dim > &p)
RefinementCase operator|(const RefinementCase &r) const
static bool is_inside_unit_cell(const Point< dim > &p)
void serialize(Archive &ar, const unsigned int version)
static Point< dim > cell_to_child_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)
static Point< dim > unit_cell_vertex(const unsigned int vertex)
#define AssertThrow(cond, exc)
Definition: exceptions.h:362
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)
RefinementCase operator~() const
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
#define Assert(cond, exc)
Definition: exceptions.h:299
unsigned char value
UpdateFlags operator|(UpdateFlags f1, UpdateFlags f2)
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static std::size_t memory_consumption()
RefinementCase operator&(const RefinementCase &r) const
::ExceptionBase & ExcNotImplemented()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:535
::ExceptionBase & ExcInternalError()
static RefinementCase cut_axis(const unsigned int i)