21 #include "libmesh/fe_interface.h" 23 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 24 #include "libmesh/fe_interface_macros.h" 25 #include "libmesh/inf_fe.h" 28 #include "libmesh/elem.h" 29 #include "libmesh/fe.h" 30 #include "libmesh/fe_compute_data.h" 31 #include "libmesh/dof_map.h" 32 #include "libmesh/enum_fe_family.h" 33 #include "libmesh/enum_order.h" 34 #include "libmesh/enum_elem_type.h" 35 #include "libmesh/enum_to_string.h" 44 libmesh_error_msg(
"ERROR: Do not define an object of this type.");
49 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS 80 #endif //ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS 82 #define fe_family_case_func(family, dim, func_and_args, prefix, suffix) \ 84 prefix FE<dim,family>::func_and_args suffix \ 86 #define fe_family_case(family) \ 89 #define fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \ 90 fe_family_case_func(CLOUGH, dim, func_and_args, prefix, suffix) \ 91 fe_family_case_func(HERMITE, dim, func_and_args, prefix, suffix) \ 92 fe_family_case_func(HIERARCHIC, dim, func_and_args, prefix, suffix) \ 93 fe_family_case_func(L2_HIERARCHIC, dim, func_and_args, prefix, suffix) \ 94 fe_family_case_func(SIDE_HIERARCHIC, dim, func_and_args, prefix, suffix) \ 95 fe_family_case_func(LAGRANGE, dim, func_and_args, prefix, suffix) \ 96 fe_family_case_func(L2_LAGRANGE, dim, func_and_args, prefix, suffix) \ 97 fe_family_case_func(MONOMIAL, dim, func_and_args, prefix, suffix) \ 98 fe_family_case_func(SCALAR, dim, func_and_args, prefix, suffix) \ 99 fe_family_case_func(XYZ, dim, func_and_args, prefix, suffix) \ 100 fe_family_case_func(SUBDIVISION, 2, func_and_args, \ 101 libmesh_assert_equal_to (dim, 2); prefix, suffix) 103 #define fe_family_scalar_case() \ 104 fe_family_case(CLOUGH) \ 105 fe_family_case(HERMITE) \ 106 fe_family_case(HIERARCHIC) \ 107 fe_family_case(L2_HIERARCHIC) \ 108 fe_family_case(SIDE_HIERARCHIC) \ 109 fe_family_case(LAGRANGE) \ 110 fe_family_case(L2_LAGRANGE) \ 111 fe_family_case(MONOMIAL) \ 112 fe_family_case(SCALAR) \ 113 fe_family_case(XYZ) \ 114 fe_family_case(SUBDIVISION) \ 116 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES 118 #define fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \ 119 fe_family_case_func(BERNSTEIN, dim, func_and_args, prefix, suffix) \ 120 fe_family_case_func(SZABAB, dim, func_and_args, prefix, suffix) \ 121 fe_family_case_func(RATIONAL_BERNSTEIN, dim, func_and_args, prefix, suffix) \ 123 #define fe_family_horder_case() \ 124 fe_family_case(BERNSTEIN) \ 125 fe_family_case(SZABAB) \ 126 fe_family_case(RATIONAL_BERNSTEIN) \ 130 #define fe_family_horder_case_func(dim, func_and_args, prefix, suffix) 131 #define fe_family_horder_case() 135 #define fe_family_vector_case_func(dim, func_and_args, prefix, suffix) \ 136 fe_family_case_func(HIERARCHIC_VEC, dim, func_and_args, prefix, suffix) \ 137 fe_family_case_func(L2_HIERARCHIC_VEC, dim, func_and_args, prefix, suffix) \ 138 fe_family_case_func(LAGRANGE_VEC, dim, func_and_args, prefix, suffix) \ 139 fe_family_case_func(L2_LAGRANGE_VEC, dim, func_and_args, prefix, suffix) \ 140 fe_family_case_func(MONOMIAL_VEC, dim, func_and_args, prefix, suffix) \ 141 fe_family_case_func(NEDELEC_ONE, dim, func_and_args, prefix, suffix) \ 142 fe_family_case_func(RAVIART_THOMAS, dim, func_and_args, prefix, suffix) \ 143 fe_family_case_func(L2_RAVIART_THOMAS, dim, func_and_args, prefix, suffix) \ 145 #define fe_family_vector_case() \ 146 fe_family_case(HIERARCHIC_VEC) \ 147 fe_family_case(L2_HIERARCHIC_VEC) \ 148 fe_family_case(LAGRANGE_VEC) \ 149 fe_family_case(L2_LAGRANGE_VEC) \ 150 fe_family_case(MONOMIAL_VEC) \ 151 fe_family_case(NEDELEC_ONE) \ 152 fe_family_case(RAVIART_THOMAS) \ 153 fe_family_case(L2_RAVIART_THOMAS) \ 155 #define fe_family_switch(dim, func_and_args, prefix, suffix) \ 157 switch (fe_t.family) \ 159 fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \ 160 fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \ 162 libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \ 166 #define fe_family_with_vec_switch(dim, func_and_args, prefix, suffix) \ 168 switch (fe_t.family) \ 170 fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \ 171 fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \ 172 fe_family_vector_case_func(dim, func_and_args, prefix, suffix) \ 174 libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \ 178 #define fe_scalar_vec_error_switch(dim, func_and_args, prefix, suffix) \ 180 switch (fe_t.family) \ 182 fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \ 183 fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \ 184 fe_family_vector_case() \ 185 libmesh_error_msg("Error: Can only request scalar valued elements for Real FEInterface::func_and_args"); \ 187 libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \ 192 #define fe_vector_scalar_error_switch(dim, func_and_args, prefix, suffix) \ 194 switch (fe_t.family) \ 196 fe_family_vector_case_func(dim, func_and_args, prefix, suffix) \ 197 fe_family_scalar_case() \ 198 fe_family_horder_case() \ 199 libmesh_error_msg("Error: Can only request vector valued elements for RealGradient FEInterface::func_and_args"); \ 201 libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \ 206 #define fe_switch(func_and_args) \ 212 fe_family_switch (0, func_and_args, return, ;); \ 215 fe_family_switch (1, func_and_args, return, ;); \ 218 fe_family_switch (2, func_and_args, return, ;); \ 221 fe_family_switch (3, func_and_args, return, ;); \ 223 libmesh_error_msg("Invalid dim = " << dim); \ 227 #define fe_with_vec_switch(func_and_args) \ 233 fe_family_with_vec_switch (0, func_and_args, return, ;); \ 236 fe_family_with_vec_switch (1, func_and_args, return, ;); \ 239 fe_family_with_vec_switch (2, func_and_args, return, ;); \ 242 fe_family_with_vec_switch (3, func_and_args, return, ;); \ 244 libmesh_error_msg("Invalid dim = " << dim); \ 248 #define void_fe_with_vec_switch(func_and_args) \ 254 fe_family_with_vec_switch (0, func_and_args, , ; return;); \ 257 fe_family_with_vec_switch (1, func_and_args, , ; return;); \ 260 fe_family_with_vec_switch (2, func_and_args, , ; return;); \ 263 fe_family_with_vec_switch (3, func_and_args, , ; return;); \ 265 libmesh_error_msg("Invalid dim = " << dim); \ 271 #ifdef LIBMESH_ENABLE_DEPRECATED 276 libmesh_deprecated();
278 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 294 #endif // LIBMESH_ENABLE_DEPRECATED 301 const bool add_p_level)
306 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 319 auto total_order = fe_t.
order + add_p_level*elem->
p_level();
321 fe_with_vec_switch(
n_dofs(elem, total_order));
328 const int extra_order,
334 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 347 auto total_order = fe_t.
order + extra_order;
349 fe_with_vec_switch(
n_dofs(elem, total_order));
354 #ifdef LIBMESH_ENABLE_DEPRECATED 359 libmesh_deprecated();
361 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 370 fe_with_vec_switch(
n_dofs(t, o));
381 libmesh_deprecated();
385 #endif // LIBMESH_ENABLE_DEPRECATED 392 const bool add_p_level)
397 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 419 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 428 auto total_order = fe_t.
order + extra_order;
430 fe_with_vec_switch(
n_dofs(elem, total_order));
435 #ifdef LIBMESH_ENABLE_DEPRECATED 439 const unsigned int n)
441 libmesh_deprecated();
443 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 461 libmesh_deprecated();
465 #endif // LIBMESH_ENABLE_DEPRECATED 484 const unsigned int n,
485 const bool add_p_level)
490 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 498 auto total_order = fe_t.
order + add_p_level*elem->
p_level();
507 const int extra_order,
509 const unsigned int n)
514 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 522 auto total_order = fe_t.
order + extra_order;
529 #ifdef LIBMESH_ENABLE_DEPRECATED 534 libmesh_deprecated();
536 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 547 #endif // LIBMESH_ENABLE_DEPRECATED 554 const bool add_p_level)
559 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 567 auto total_order = fe_t.
order + add_p_level*elem->
p_level();
576 const int extra_order,
582 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 590 auto total_order = fe_t.
order + extra_order;
598 const unsigned int dim,
601 std::vector<unsigned int> & di,
602 const bool add_p_level)
606 void_fe_with_vec_switch(
dofs_on_side(elem, o, s, di, add_p_level));
612 const unsigned int dim,
615 std::vector<unsigned int> & di,
616 const bool add_p_level)
620 void_fe_with_vec_switch(
dofs_on_edge(elem, o, e, di, add_p_level));
628 const std::vector<Number> & elem_soln,
629 std::vector<Number> & nodal_soln,
630 const bool add_p_level,
631 const unsigned int vdim)
633 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 652 const unsigned int side,
653 const std::vector<Number> & elem_soln,
654 std::vector<Number> & nodal_soln,
655 const bool add_p_level,
656 const unsigned int vdim)
658 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 662 libmesh_not_implemented();
669 const unsigned int dim = elem->
dim();
676 #ifdef LIBMESH_ENABLE_DEPRECATED 682 libmesh_deprecated();
683 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 687 fe_with_vec_switch(
map(elem, p));
696 const Real tolerance,
699 libmesh_deprecated();
700 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 707 fe_with_vec_switch(
inverse_map(elem, p, tolerance, secure));
715 const std::vector<Point> & physical_points,
716 std::vector<Point> & reference_points,
717 const Real tolerance,
720 libmesh_deprecated();
722 const std::size_t n_pts = physical_points.size();
725 reference_points.resize(n_pts);
729 libMesh::err <<
"WARNING: empty vector physical_points!" 735 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 746 void_fe_with_vec_switch(
inverse_map(elem, physical_points, reference_points, tolerance, secure));
763 const unsigned int i,
766 libmesh_deprecated();
768 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 777 fe_switch(
shape(t,o,i,p));
783 const unsigned int i,
786 libmesh_deprecated();
788 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 797 fe_switch(
shape(elem,o,i,p));
799 #endif // LIBMESH_ENABLE_DEPRECATED 806 const unsigned int i,
808 const bool add_p_level)
813 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 825 fe_switch(
shape(elem, fe_t.
order, i, p, add_p_level));
834 const unsigned int i,
840 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 854 auto total_order = fe_t.
order + extra_order;
856 fe_switch(
shape(elem, total_order, i, p,
false));
861 #ifdef LIBMESH_ENABLE_DEPRECATED 863 void FEInterface::shape<Real>(
const unsigned int dim,
866 const unsigned int i,
870 libmesh_deprecated();
872 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 874 if (is_InfFE_elem(t))
876 phi = ifem_shape(
dim, fe_t, t, i, p);
882 const Order o = fe_t.order;
887 fe_scalar_vec_error_switch(0, shape(t,o,i,p), phi = , ;
break;);
890 fe_scalar_vec_error_switch(1, shape(t,o,i,p), phi = , ;
break;);
893 fe_scalar_vec_error_switch(2, shape(t,o,i,p), phi = , ;
break;);
896 fe_scalar_vec_error_switch(3, shape(t,o,i,p), phi = , ;
break;);
899 libmesh_error_msg(
"Invalid dimension = " <<
dim);
908 void FEInterface::shape<Real>(
const unsigned int dim,
911 const unsigned int i,
915 libmesh_deprecated();
917 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 919 if (elem && is_InfFE_elem(elem->type()))
921 phi = ifem_shape(fe_t, elem, i, p);
926 const Order o = fe_t.order;
931 fe_scalar_vec_error_switch(0, shape(elem,o,i,p), phi = , ;
break;);
934 fe_scalar_vec_error_switch(1, shape(elem,o,i,p), phi = , ;
break;);
937 fe_scalar_vec_error_switch(2, shape(elem,o,i,p), phi = , ;
break;);
940 fe_scalar_vec_error_switch(3, shape(elem,o,i,p), phi = , ;
break;);
943 libmesh_error_msg(
"Invalid dimension = " <<
dim);
948 #endif // LIBMESH_ENABLE_DEPRECATED 953 void FEInterface::shape<Real>(
const FEType & fe_t,
955 const unsigned int i,
960 auto dim = elem->dim();
962 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 964 if (is_InfFE_elem(elem->type()))
966 phi = ifem_shape(fe_t, elem, i, p);
977 fe_scalar_vec_error_switch(0, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
980 fe_scalar_vec_error_switch(1, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
983 fe_scalar_vec_error_switch(2, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
986 fe_scalar_vec_error_switch(3, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
989 libmesh_error_msg(
"Invalid dimension = " <<
dim);
996 void FEInterface::shape<Real>(
const FEType & fe_t,
999 const unsigned int i,
1004 auto dim = elem->dim();
1006 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1008 if (is_InfFE_elem(elem->type()))
1010 phi = ifem_shape(fe_t, elem, i, p);
1017 auto total_order = fe_t.order + extra_order;
1028 fe_scalar_vec_error_switch(0, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1031 fe_scalar_vec_error_switch(1, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1034 fe_scalar_vec_error_switch(2, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1037 fe_scalar_vec_error_switch(3, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1040 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1047 void FEInterface::shapes<Real>(
const unsigned int dim,
1050 const unsigned int i,
1051 const std::vector<Point> & p,
1052 std::vector<Real> & phi,
1053 const bool add_p_level)
1055 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1057 if (elem && is_InfFE_elem(elem->type()))
1060 elevated.
order = fe_t.order + add_p_level * elem->p_level();
1062 phi[qpi] = ifem_shape(elevated, elem, i, p[qpi]);
1067 const Order o = fe_t.order;
1072 fe_scalar_vec_error_switch(0, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1075 fe_scalar_vec_error_switch(1, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1078 fe_scalar_vec_error_switch(2, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1081 fe_scalar_vec_error_switch(3, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1084 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1092 void FEInterface::all_shapes<Real>(
const unsigned int dim,
1095 const std::vector<Point> & p,
1096 std::vector<std::vector<Real>> & phi,
1097 const bool add_p_level)
1099 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1101 if (elem && is_InfFE_elem(elem->type()))
1104 FEInterface::shapes<Real>(
dim, fe_t, elem, i, p, phi[i], add_p_level);
1109 const Order o = fe_t.order;
1114 fe_scalar_vec_error_switch(0, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1117 fe_scalar_vec_error_switch(1, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1120 fe_scalar_vec_error_switch(2, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1123 fe_scalar_vec_error_switch(3, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1126 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1134 #ifdef LIBMESH_ENABLE_DEPRECATED 1136 void FEInterface::shape<RealGradient>(
const unsigned int dim,
1139 const unsigned int i,
1143 libmesh_deprecated();
1146 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1147 if (is_InfFE_elem(t))
1149 libmesh_not_implemented();
1152 const Order o = fe_t.order;
1157 fe_vector_scalar_error_switch(0, shape(t,o,i,p), phi = , ;
break;);
1160 fe_vector_scalar_error_switch(1, shape(t,o,i,p), phi = , ;
break;);
1163 fe_vector_scalar_error_switch(2, shape(t,o,i,p), phi = , ;
break;);
1166 fe_vector_scalar_error_switch(3, shape(t,o,i,p), phi = , ;
break;);
1169 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1172 #endif // LIBMESH_ENABLE_DEPRECATED 1177 void FEInterface::shape<RealGradient>(
const FEType & fe_t,
1179 const unsigned int i,
1184 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1185 if (is_InfFE_elem(elem->type()))
1187 libmesh_not_implemented();
1191 auto dim = elem->dim();
1203 fe_vector_scalar_error_switch(0, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
1206 fe_vector_scalar_error_switch(1, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
1209 fe_vector_scalar_error_switch(2, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
1212 fe_vector_scalar_error_switch(3, shape(elem, fe_t.order, i, p,
true), phi = , ;
break;);
1215 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1222 void FEInterface::shape<RealGradient>(
const FEType & fe_t,
1225 const unsigned int i,
1230 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1231 if (is_InfFE_elem(elem->type()))
1233 libmesh_not_implemented();
1237 auto dim = elem->dim();
1246 auto total_order = fe_t.order + extra_order;
1251 fe_vector_scalar_error_switch(0, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1254 fe_vector_scalar_error_switch(1, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1257 fe_vector_scalar_error_switch(2, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1260 fe_vector_scalar_error_switch(3, shape(elem, total_order, i, p,
false), phi = , ;
break;);
1263 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1270 void FEInterface::shapes<RealGradient>(
const unsigned int dim,
1273 const unsigned int i,
1274 const std::vector<Point> & p,
1275 std::vector<RealGradient> & phi,
1276 const bool add_p_level)
1279 if (elem->infinite())
1280 libmesh_not_implemented();
1282 const Order o = fe_t.order;
1287 fe_vector_scalar_error_switch(0, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1290 fe_vector_scalar_error_switch(1, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1293 fe_vector_scalar_error_switch(2, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1296 fe_vector_scalar_error_switch(3, shapes(elem,o,i,p,phi,add_p_level), , ;
return;);
1299 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1308 void FEInterface::all_shapes<RealGradient>(
const unsigned int dim,
1311 const std::vector<Point> & p,
1312 std::vector<std::vector<RealGradient>> & phi,
1313 const bool add_p_level)
1316 if (elem->infinite())
1317 libmesh_not_implemented();
1319 const Order o = fe_t.order;
1324 fe_vector_scalar_error_switch(0, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1327 fe_vector_scalar_error_switch(1, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1330 fe_vector_scalar_error_switch(2, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1333 fe_vector_scalar_error_switch(3, all_shapes(elem,o,p,phi,add_p_level), , ;
return;);
1336 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1343 #ifdef LIBMESH_ENABLE_DEPRECATED 1349 libmesh_deprecated();
1351 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1354 inf_fe_switch(
shape);
1359 #endif // LIBMESH_ENABLE_DEPRECATED 1370 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1372 inf_fe_switch(
shape);
1379 #ifdef LIBMESH_ENABLE_DEPRECATED 1383 const unsigned int i,
1384 const unsigned int j,
1387 libmesh_deprecated();
1389 libmesh_assert_greater (
dim,j);
1390 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1407 const unsigned int i,
1408 const unsigned int j,
1411 libmesh_deprecated();
1413 libmesh_assert_greater (
dim,j);
1414 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1427 fe_family_switch (0,
shape_deriv(elem, o, i, j, p),
return , ;);
1429 fe_family_switch (1,
shape_deriv(elem, o, i, j, p),
return , ;);
1431 fe_family_switch (2,
shape_deriv(elem, o, i, j, p),
return , ;);
1433 fe_family_switch (3,
shape_deriv(elem, o, i, j, p),
return , ;);
1435 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1439 #endif // LIBMESH_ENABLE_DEPRECATED 1445 const unsigned int i,
1446 const unsigned int j,
1451 libmesh_assert_greater (
dim, j);
1452 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1474 const unsigned int i,
1475 const unsigned int j,
1480 libmesh_assert_greater (
dim, j);
1481 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1491 auto total_order = fe_t.
order + extra_order;
1498 fe_family_switch (0,
shape_deriv(elem, total_order, i, j, p,
false),
return , ;);
1500 fe_family_switch (1,
shape_deriv(elem, total_order, i, j, p,
false),
return , ;);
1502 fe_family_switch (2,
shape_deriv(elem, total_order, i, j, p,
false),
return , ;);
1504 fe_family_switch (3,
shape_deriv(elem, total_order, i, j, p,
false),
return , ;);
1506 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1516 void FEInterface::shape_derivs<Real>(
const FEType & fe_t,
1518 const unsigned int i,
1519 const unsigned int j,
1520 const std::vector<Point> & p,
1521 std::vector<Real> & dphi,
1522 const bool add_p_level)
1524 const Order o = fe_t.order;
1529 fe_scalar_vec_error_switch(0, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1532 fe_scalar_vec_error_switch(1, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1535 fe_scalar_vec_error_switch(2, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1538 fe_scalar_vec_error_switch(3, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1541 libmesh_error_msg(
"Invalid dimension = " << elem->dim());
1552 void FEInterface::all_shape_derivs<Real>(
const unsigned int dim,
1555 const std::vector<Point> & p,
1556 std::vector<std::vector<Real>> * comps[3],
1557 const bool add_p_level)
1559 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1561 if (elem && is_InfFE_elem(elem->type()))
1565 FEInterface::shape_derivs<Real>(fe_t, elem, i, j, p, (*comps[j])[i], add_p_level);
1570 const Order o = fe_t.order;
1575 fe_scalar_vec_error_switch(0, all_shape_derivs(elem,o,p,comps,add_p_level), , ;
return;);
1578 fe_scalar_vec_error_switch(1, all_shape_derivs(elem,o,p,comps,add_p_level), , ;
return;);
1581 fe_scalar_vec_error_switch(2, all_shape_derivs(elem,o,p,comps,add_p_level), , ;
return;);
1584 fe_scalar_vec_error_switch(3, all_shape_derivs(elem,o,p,comps,add_p_level), , ;
return;);
1587 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1596 void FEInterface::shape_derivs<RealGradient>(
const FEType & fe_t,
1598 const unsigned int i,
1599 const unsigned int j,
1600 const std::vector<Point> & p,
1601 std::vector<RealGradient> & dphi,
1602 const bool add_p_level)
1604 const Order o = fe_t.order;
1609 fe_vector_scalar_error_switch(0, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1612 fe_vector_scalar_error_switch(1, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1615 fe_vector_scalar_error_switch(2, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1618 fe_vector_scalar_error_switch(3, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ;
break;);
1621 libmesh_error_msg(
"Invalid dimension = " << elem->dim());
1629 #ifdef LIBMESH_ENABLE_DEPRECATED 1635 libmesh_deprecated();
1637 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1645 #endif // LIBMESH_ENABLE_DEPRECATED 1656 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1666 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1667 #ifdef LIBMESH_ENABLE_DEPRECATED 1671 const unsigned int i,
1672 const unsigned int j,
1675 libmesh_deprecated();
1677 libmesh_assert_greater_equal (
dim*(
dim-1),j);
1678 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1680 libmesh_not_implemented();
1696 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1705 const unsigned int i,
1706 const unsigned int j,
1709 libmesh_deprecated();
1711 libmesh_assert_greater_equal (
dim*(
dim-1),j);
1714 libmesh_not_implemented();
1729 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1733 #endif // LIBMESH_ENABLE_DEPRECATED 1739 const unsigned int i,
1740 const unsigned int j,
1745 libmesh_assert_greater_equal (
dim*(
dim-1),j);
1746 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1748 libmesh_not_implemented();
1765 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1777 const unsigned int i,
1778 const unsigned int j,
1783 libmesh_assert_greater_equal (
dim*(
dim-1),j);
1786 libmesh_not_implemented();
1790 auto total_order = fe_t.
order + extra_order;
1798 fe_family_switch (0,
shape_second_deriv(elem, total_order, i, j, p,
false),
return , ;);
1800 fe_family_switch (1,
shape_second_deriv(elem, total_order, i, j, p,
false),
return , ;);
1802 fe_family_switch (2,
shape_second_deriv(elem, total_order, i, j, p,
false),
return , ;);
1804 fe_family_switch (3,
shape_second_deriv(elem, total_order, i, j, p,
false),
return , ;);
1806 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1815 #ifdef LIBMESH_ENABLE_DEPRECATED 1821 libmesh_deprecated();
1823 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1825 libmesh_not_implemented();
1829 #endif // LIBMESH_ENABLE_DEPRECATED 1839 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1841 libmesh_not_implemented();
1846 #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES 1849 #ifdef LIBMESH_ENABLE_DEPRECATED 1851 void FEInterface::shape<RealGradient>(
const unsigned int dim,
1854 const unsigned int i,
1858 libmesh_deprecated();
1861 if (elem->infinite())
1862 libmesh_not_implemented();
1864 const Order o = fe_t.order;
1869 fe_vector_scalar_error_switch(0, shape(elem,o,i,p), phi = , ;
break;);
1872 fe_vector_scalar_error_switch(1, shape(elem,o,i,p), phi = , ;
break;);
1875 fe_vector_scalar_error_switch(2, shape(elem,o,i,p), phi = , ;
break;);
1878 fe_vector_scalar_error_switch(3, shape(elem,o,i,p), phi = , ;
break;);
1881 libmesh_error_msg(
"Invalid dimension = " <<
dim);
1886 #endif // LIBMESH_ENABLE_DEPRECATED 1894 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 1905 const unsigned int n_dof =
n_dofs (fe_t, elem);
1906 const Point & p = data.
p;
1907 data.
shape.resize(n_dof);
1911 data.
dshape.resize(n_dof);
1914 for (
unsigned int d=0; d<
dim; d++)
1918 std::vector<Point> pt = {p};
1920 fe->reinit(elem, &pt);
1943 for (
unsigned int n=0; n<n_dof; n++)
1953 for (
unsigned int j=0; j<
dim; j++)
1961 #ifdef LIBMESH_ENABLE_AMR 1965 const unsigned int variable_number,
1972 switch (elem->
dim())
2089 libmesh_error_msg(
"Invalid dimension = " << elem->
dim());
2093 #endif // #ifdef LIBMESH_ENABLE_AMR 2097 #ifdef LIBMESH_ENABLE_PERIODIC 2104 const unsigned int variable_number,
2115 constraints, dof_map, boundaries,
mesh, point_locator, variable_number, elem);
2119 constraints, dof_map, boundaries,
mesh, point_locator, variable_number, elem);
2123 "compute_periodic_constraints only set up for vector or scalar FEFieldTypes");
2127 #endif // #ifdef LIBMESH_ENABLE_PERIODIC 2138 const unsigned int unlimited = 11;
2142 const unsigned int unknown = unlimited;
2245 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES 2636 fe_family_horder_case()
2651 fe_family_vector_case()
Real(* shape_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function derivative values.
static unsigned int n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
FEFamily family
The type of finite element.
static Point map(unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
This is now deprecated; use FEMap::map instead.
bool need_derivative()
Check whether derivatives should be computed or not.
ElemType
Defines an enum for geometric element types.
Order
defines an enum for polynomial orders.
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
static shape_ptr shape_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
static unsigned int ifem_n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
static Real shape_deriv(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const unsigned int j, const Point &p)
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
class FEComputeData hides arbitrary data to be passed to and from children of FEBase through the FEIn...
We're using a class instead of a typedef to allow forward declarations and future flexibility...
const Point & p
Holds the point where the data are to be computed.
static void dofs_on_edge(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int e, std::vector< unsigned int > &di, const bool add_p_level=true)
Fills the vector di with the local degree of freedom indices associated with edge e of element elem A...
const FEType & variable_type(const unsigned int c) const
This is the base class from which all geometric element types are derived.
static FEFieldType field_type(const FEType &fe_type)
std::vector< Gradient > dshape
Storage for the computed shape derivative values.
unsigned int p_level() const
OrderWrapper order
The approximation order of the element.
std::vector< std::vector< Real > > local_transform
Storage for local to global mapping at p.
The libMesh namespace provides an interface to certain functionality in the library.
static Point ifem_inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
void init()
Inits the output data to default values, provided the fields are correctly resized.
static unsigned int max_order(const FEType &fe_t, const ElemType &el_t)
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
static void side_nodal_soln(const FEType &fe_t, const Elem *elem, const unsigned int side, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln, const bool add_p_level=true, const unsigned int vdim=1)
Build the nodal soln on one side from the (full) element soln.
This is the MeshBase class.
static bool orientation_dependent(const FEFamily &fe_family)
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
This class handles the numbering of degrees of freedom on a mesh.
FEInterface()
Empty constructor.
static void ifem_compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
static bool extra_hanging_dofs(const FEType &fe_t)
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
static bool is_hierarchic(const FEType &fe_type)
Returns whether or not the input FEType's higher-order shape functions are always hierarchic...
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
This is now deprecated; use FEMap::inverse_map instead.
static Point ifem_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
This is the base class for point locators.
static void ifem_nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
static shape_deriv_ptr shape_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
static unsigned int ifem_n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
static unsigned int n_vec_dim(const MeshBase &mesh, const FEType &fe_type)
std::vector< Number > shape
Storage for the computed shape function values.
static n_dofs_at_node_ptr n_dofs_at_node_function(const unsigned int dim, const FEType &fe_t)
static void compute_periodic_constraints(DofConstraints &constraints, DofMap &dof_map, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for meshes with periodic boundary conditions) correspon...
std::string enum_to_string(const T e)
static FEContinuity get_continuity(const FEType &fe_type)
Returns the input FEType's FEContinuity based on the underlying FEFamily and potentially the Order...
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
static Real shape_second_deriv(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const unsigned int j, const Point &p)
static unsigned int ifem_n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Real(* shape_ptr)(const FEType fe_t, const Elem *elem, const unsigned int i, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function values.
unsigned int(* n_dofs_at_node_ptr)(const ElemType, const Order, const unsigned int)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static shape_second_deriv_ptr shape_second_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
virtual unsigned short dim() const =0
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
static bool is_InfFE_elem(const ElemType et)
Real(* shape_second_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function second derivative values...
virtual bool infinite() const =0
FEFamily
defines an enum for finite element families.
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln, const bool add_p_level=true, const unsigned int vdim=1)
Build the nodal soln from the element soln.
static void compute_periodic_constraints(DofConstraints &constraints, DofMap &dof_map, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for periodic boundary conditions) corresponding to vari...
static Real ifem_shape_deriv(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const unsigned int j, const Point &p)
virtual ElemType type() const =0
static void dofs_on_side(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int s, std::vector< unsigned int > &di, const bool add_p_level=true)
Fills the vector di with the local degree of freedom indices associated with side s of element elem A...
A Point defines a location in LIBMESH_DIM dimensional Real space.
The constraint matrix storage format.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
static unsigned int ifem_n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
FEFieldType
defines an enum for finite element field types - i.e.