libMesh
Classes | Public Types | Public Member Functions | Public Attributes | List of all members
InfFERadialTest Class Reference

This class is for unit testing the "radial" basis function aspects of the InfFE. More...

Inheritance diagram for InfFERadialTest:
[legend]

Classes

struct  TabulatedGrad
 
struct  TabulatedVal
 

Public Types

typedef std::vector< TabulatedValTabulatedVals
 
typedef std::vector< TabulatedGradTabulatedGrads
 
typedef std::pair< Order, FEFamilyKeyType
 

Public Member Functions

 LIBMESH_CPPUNIT_TEST_SUITE (InfFERadialTest)
 
 CPPUNIT_TEST (testDifferentOrders)
 
 CPPUNIT_TEST (testInfQuants)
 
 CPPUNIT_TEST (testSides)
 
 CPPUNIT_TEST (testInfQuants_numericDeriv)
 
 CPPUNIT_TEST (testRefinement)
 
 CPPUNIT_TEST_SUITE_END ()
 
Point base_point (const Point physical_point, const Elem *inf_elem)
 
void testSides ()
 
void testInfQuants ()
 
void testInfQuants_numericDeriv ()
 
void testDifferentOrders ()
 
void testSingleOrder (Order radial_order, FEFamily radial_family)
 
void testRefinement ()
 
void setUp ()
 
void tearDown ()
 

Public Attributes

std::map< KeyType, TabulatedValstabulated_vals
 
std::map< KeyType, TabulatedGradstabulated_grads
 

Detailed Description

This class is for unit testing the "radial" basis function aspects of the InfFE.

There are several possible choices for these basis functions, including:

Definition at line 46 of file inf_fe_radial_test.C.

Member Typedef Documentation

◆ KeyType

Definition at line 217 of file inf_fe_radial_test.C.

◆ TabulatedGrads

Definition at line 216 of file inf_fe_radial_test.C.

◆ TabulatedVals

Definition at line 215 of file inf_fe_radial_test.C.

Member Function Documentation

◆ base_point()

Point InfFERadialTest::base_point ( const Point  physical_point,
const Elem inf_elem 
)
inline

Definition at line 78 of file inf_fe_radial_test.C.

References libMesh::TypeVector< T >::add_scaled(), libMesh::InfFEBase::build_elem(), libMesh::err, libMesh::FEMap::inverse_map(), libMesh::libmesh_assert(), libMesh::libmesh_ignore(), libMesh::libmesh_isinf(), libMesh::FE< Dim, T >::n_shape_functions(), libMesh::TypeVector< T >::norm(), libMesh::Elem::origin(), libMesh::Real, libMesh::FE< Dim, T >::shape(), libMesh::FE< Dim, T >::shape_deriv(), libMesh::TypeTensor< T >::solve(), and libMesh::TOLERANCE.

79  {
80 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
81  libmesh_assert(inf_elem);
82 
83  // The strategy is:
84  // compute the intersection of the line
85  // physical_point - origin with the base element,
86  // find its internal coordinatels using FEMap::inverse_map():
87  // The radial part can then be computed directly later on.
88 
89  // 1.)
90  // build a base element to do the map inversion in the base face
91  std::unique_ptr<const Elem> base_elem = InfFEBase::build_elem (inf_elem);
92 
93  // the origin of the infinite element
94  const Point o = inf_elem->origin();
95  // 2.)
96  // Now find the intersection of a plane represented by the base
97  // element nodes and the line given by the origin of the infinite
98  // element and the physical point.
99  Point intersection;
100 
101  const Point xi ( base_elem->point(1) - base_elem->point(0));
102  const Point eta( base_elem->point(2) - base_elem->point(0));
103  const Point zeta( physical_point - o);
104 
105  // normal vector of the base elements plane
106  Point n(xi.cross(eta));
107  Real c_factor = (base_elem->point(0) - o)*n/(zeta*n) - 1.;
108 
109  // Check whether the above system is ill-posed. It should
110  // only happen when \p physical_point is not in \p
111  // inf_elem. In this case, any point that is not in the
112  // element is a valid answer.
113  if (libmesh_isinf(c_factor))
114  return Point(0., 0., -2.);
115 
116  // Compute the intersection with
117  // {intersection} = {physical_point} + c*({physical_point}-{o}).
118  intersection.add_scaled(physical_point,1.+c_factor);
119  intersection.add_scaled(o,-c_factor);
120 
121  if (!base_elem->has_affine_map())
122  {
123  unsigned int iter_max = 20;
124 
125  // the number of shape functions needed for the base_elem
126  unsigned int n_sf = FE<2,LAGRANGE>::n_shape_functions(base_elem->type(),base_elem->default_order());
127 
128  // guess base element coordinates: p=xi,eta,0
129  // in first iteration, never run with 'secure' or
130  // 'extra_checks' to avoid false warnings.
131  Point ref_point= FEMap::inverse_map(2, base_elem.get(), intersection,
132  TOLERANCE, false, false);
133 
134  // Newton iteration
135  for(unsigned int it=0; it<iter_max; ++it)
136  {
137  // Get the shape function and derivative values at the reference coordinate
138  // phi.size() == dphi.size()
139  Point dxyz_dxi;
140  Point dxyz_deta;
141 
142  Point intersection_guess;
143  for(unsigned int i=0; i<n_sf; ++i)
144  {
145 
146  intersection_guess += base_elem->node_ref(i) * FE<2,LAGRANGE>::shape(base_elem->type(),
147  base_elem->default_order(),
148  i,
149  ref_point);
150 
151  dxyz_dxi += base_elem->node_ref(i) * FE<2,LAGRANGE>::shape_deriv(base_elem->type(),
152  base_elem->default_order(),
153  i,
154  0, // d()/dxi
155  ref_point);
156 
157  dxyz_deta+= base_elem->node_ref(i) * FE<2,LAGRANGE>::shape_deriv(base_elem->type(),
158  base_elem->default_order(),
159  i,
160  1, // d()/deta
161  ref_point);
162  } // for i
163 
164  TypeVector<Real> F(physical_point + c_factor*(physical_point-o) - intersection_guess);
165 
167  J(0,0) = (physical_point-o)(0);
168  J(0,1) = -dxyz_dxi(0);
169  J(0,2) = -dxyz_deta(0);
170  J(1,0) = (physical_point-o)(1);
171  J(1,1) = -dxyz_dxi(1);
172  J(1,2) = -dxyz_deta(1);
173  J(2,0) = (physical_point-o)(2);
174  J(2,1) = -dxyz_dxi(2);
175  J(2,2) = -dxyz_deta(2);
176 
177  // delta will be the newton step
178  TypeVector<Real> delta;
179  J.solve(F,delta);
180 
181  // check for convergence
182  Real tol = std::min( TOLERANCE*Real(0.01), TOLERANCE*base_elem->hmax() );
183  if ( delta.norm() < tol )
184  {
185  // newton solver converged, now make sure it converged to a point on the base_elem
186  if (base_elem->contains_point(intersection_guess,TOLERANCE*0.1))
187  {
188  intersection = intersection_guess;
189  }
190  else
191  {
192  err<<" No! This inverse_map failed!"<<std::endl;
193  }
194  break; // break out of 'for it'
195  }
196  else
197  {
198  c_factor -= delta(0);
199  ref_point(0) -= delta(1);
200  ref_point(1) -= delta(2);
201  }
202 
203  }
204 
205  }
206 
207  return intersection;
208 #else
209  libmesh_ignore(physical_point, inf_elem);
210  // lets make the compilers happy:
211  return Point(0.,0.,-2.);
212 #endif // LIBMESH_ENABLE_INFINITE_ELEMENTS
213  }
OStreamProxy err
void add_scaled(const TypeVector< T2 > &, const T &)
Add a scaled value to this vector without creating a temporary.
Definition: type_vector.h:651
auto norm() const -> decltype(std::norm(T()))
Definition: type_vector.h:929
virtual Point origin() const
Definition: elem.h:1767
static constexpr Real TOLERANCE
void solve(const TypeVector< T > &b, TypeVector< T > &x) const
Solve the 2x2 or 3x3 system of equations A*x = b for x by directly inverting A and multiplying it by ...
Definition: type_tensor.h:1129
bool libmesh_isinf(T x)
This class defines a tensor in LIBMESH_DIM dimensional space of type T.
Definition: tensor_tools.h:36
A specific instantiation of the FEBase class.
Definition: fe.h:127
void libmesh_ignore(const Args &...)
libmesh_assert(ctx)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39

◆ CPPUNIT_TEST() [1/5]

InfFERadialTest::CPPUNIT_TEST ( testDifferentOrders  )

◆ CPPUNIT_TEST() [2/5]

InfFERadialTest::CPPUNIT_TEST ( testInfQuants  )

◆ CPPUNIT_TEST() [3/5]

InfFERadialTest::CPPUNIT_TEST ( testSides  )

◆ CPPUNIT_TEST() [4/5]

InfFERadialTest::CPPUNIT_TEST ( testInfQuants_numericDeriv  )

◆ CPPUNIT_TEST() [5/5]

InfFERadialTest::CPPUNIT_TEST ( testRefinement  )

◆ CPPUNIT_TEST_SUITE_END()

InfFERadialTest::CPPUNIT_TEST_SUITE_END ( )

◆ LIBMESH_CPPUNIT_TEST_SUITE()

InfFERadialTest::LIBMESH_CPPUNIT_TEST_SUITE ( InfFERadialTest  )

◆ setUp()

void InfFERadialTest::setUp ( )
inline

Definition at line 970 of file inf_fe_radial_test.C.

References libMesh::FIRST, libMesh::FOURTH, libMesh::JACOBI_20_00, libMesh::JACOBI_30_00, libMesh::LEGENDRE, libMesh::SECOND, and libMesh::THIRD.

971  {
973  // Arbitrarily selected LEGENDRE reference values.
975  tabulated_vals[std::make_pair(FIRST, LEGENDRE)] = /* LEGENDRE = 14 */
976  {
977  {0,9,0.0550016}, {1,5,0.41674}, {2,8,0.0147376}, {3,7,0.111665},
978  {4,6,0.0737011}, {5,1,0.0803773}, {6,11,0.275056}, {7,15,0.021537}
979  };
980 
981  tabulated_vals[std::make_pair(SECOND, LEGENDRE)] =
982  {
983  {0,3,0.0425633}, {1,6,0.0343526}, {2,2,0.158848}, {3,7,0.128206},
984  {4,12,0.220829}, {5,5,-0.136549}, {6,0,0.0149032}, {7,10,-0.0334936},
985  {8,16,0.00399329}, {9,18,-0.00209733}, {10,2,0.0556194}, {11,17,-0.000561977}
986  };
987 
988  tabulated_vals[std::make_pair(THIRD, LEGENDRE)] =
989  {
990  {12,9,0.136657}, {13,6,0.175034}, {14,20,-0.0011016}, {15,15,0.0428935}
991  };
992 
993  tabulated_vals[std::make_pair(FOURTH, LEGENDRE)] =
994  {
995  {16,3,0.00826618}, {17,10,-0.547813}, {18,14,0.311004}, {19,27,-0.00192085}
996  };
997 
999  // Arbitrarily selected LEGENDRE reference gradients.
1001  tabulated_grads[std::make_pair(FIRST, LEGENDRE)] =
1002  {
1003  {0,9,Point(-0.0429458,-0.0115073,0.)},
1004  {1,5,Point(0.177013,-0.177013,-0.483609)},
1005  {2,8,Point(0.0115073,0.0115073,0.00842393)},
1006  {3,7,Point(-0.177013,0.0474305,0.)},
1007  {4,6,Point(-0.031305,-0.116832,0.10025)},
1008  {5,1,Point(0.0474191,-0.0474191,0.872917)},
1009  {6,11,Point(0.0575466,0.0575466,-0.11251)},
1010  {7,15,Point(-0.00353805,0.000948017,0.000111572)},
1011  };
1012 
1013  tabulated_grads[std::make_pair(SECOND, LEGENDRE)] =
1014  {
1015  {0,3,Point(-0.0959817, -0.0959817, 0.0702635)},
1016  {1,6,Point(0.0625228, -0.0625228, 0.0457699)},
1017  {2,2,Point(0.358209, 0.0959817, 0.)},
1018  {3,7,Point(-0.233338, 0.0625228, 0.)},
1019  {4,12,Point(-0.0323071, -0.0323071, -0.0729771)},
1020  {5,5,Point(0.248523, 0.0665915, -0.245097)},
1021  {6,0,Point(0.0336072, -0.00900502, 0.288589)},
1022  {7,10,Point(-0.0396234, 0.0396234, -0.0290064)},
1023  {8,16,Point(0.000443217, 0.000443217, 0.000333678)},
1024  {9,18,Point(-0.000232783, -6.23741e-05, 9.35433e-05)},
1025  {10,2,Point(-0.0336072, 0.0336072, 0.985214)},
1026  {11,17,Point(6.23741e-05, -6.23741e-05, -2.05961e-05)},
1027  };
1028 
1029  tabulated_grads[std::make_pair(THIRD, LEGENDRE)] =
1030  {
1031  {12,9,Point(0.0536552, 0.200244, -0.0849541)},
1032  {13,6,Point(-0.0921697, 0.0921697, 0.461056)},
1033  {14,20,Point(2.35811e-05, -8.80059e-05, 3.58959e-05)},
1034  {15,15,Point(-0.0386352, 0.0103523, -0.0197323)},
1035  };
1036 
1037  tabulated_grads[std::make_pair(FOURTH, LEGENDRE)] =
1038  {
1039  {16,3,Point(-0.0190603, 0.0051072, 0.308529)},
1040  {17,10,Point(0.244125, -0.244125, 0.140907)},
1041  {18,14,Point(-0.0985844, 0.0985844, -0.502591)},
1042  {19,27,Point(0.000115647, -3.09874e-05, 4.30775e-05)}
1043  };
1044 
1046  // Arbitrarily selected JACOBI_20_00 reference values.
1048  tabulated_vals[std::make_pair(FIRST, JACOBI_20_00)] = /* JACOBI_20_00 = 12 */
1049  {
1050  // These values are the same as Legendre
1051  {0,9,0.0550016}, {1,5,0.41674}, {2,8,0.0147376}, {3,7,0.111665},
1052  // These values are different from Legendre
1053  {4,6,0.147402}, {5,1,0.160755}, {6,11,0.550112}, {7,15,0.043074}
1054  };
1055 
1056  tabulated_vals[std::make_pair(SECOND, JACOBI_20_00)] =
1057  {
1058  // These values are the same as Legendre
1059  {0,3,0.0425633}, {1,6,0.0343526}, {2,2,0.158848}, {3,7,0.128206},
1060  // These values are different from Legendre
1061  {4,12,0.441658}, {5,5,-0.193445}, {6,0,0.0298063}, {7,10,-0.0279114},
1062  {8,16,0.00798659}, {9,18,0.0320146}, {10,2,0.111239}, {11,17,0.00857829}
1063  };
1064 
1065  tabulated_vals[std::make_pair(THIRD, JACOBI_20_00)] =
1066  {
1067  {12,9,0.0837876}, {13,6,0.350068}, {14,20,0.0244336}, {15,15,0.0181532}
1068  };
1069 
1070  tabulated_vals[std::make_pair(FOURTH, JACOBI_20_00)] =
1071  {
1072  {16,3,0.0165324}, {17,10,-0.720085}, {18,14,0.0777511}, {19,27,0.0453842}
1073  };
1074 
1076  // Arbitrarily selected JACOBI_20_00 reference gradients.
1078  tabulated_grads[std::make_pair(FIRST, JACOBI_20_00)] =
1079  {
1080  // These values are the same as Legendre
1081  {0,9,Point(-0.0429458, -0.0115073, 0.)},
1082  {1,5,Point(0.177013, -0.177013, -0.483609)},
1083  {2,8,Point(0.0115073, 0.0115073, 0.00842393)},
1084  {3,7,Point(-0.177013, 0.0474305, 0.)},
1085  // These values are different from Legendre
1086  {4,6,Point(-0.0626101, -0.233664, 0.2005)},
1087  {5,1,Point(0.0948382, -0.0948382, 1.74583)},
1088  {6,11,Point(0.115093, 0.115093, -0.22502)},
1089  {7,15,Point(-0.0070761, 0.00189603, 0.000223144)}
1090  };
1091 
1092  tabulated_grads[std::make_pair(SECOND, JACOBI_20_00)] =
1093  {
1094  // These values are the same as Legendre
1095  {0,3,Point(-0.0959817, -0.0959817, 0.0702635)},
1096  {1,6,Point(0.0625228, -0.0625228, 0.0457699)},
1097  {2,2,Point(0.358209, 0.0959817, 0.)},
1098  {3,7,Point(-0.233338, 0.0625228, 0.)},
1099  // These values are different from Legendre
1100  {4,12,Point(-0.0646142, -0.0646142, -0.145954)},
1101  {5,5,Point(0.352076, 0.0943384, -0.233431)},
1102  {6,0,Point(0.0672144, -0.01801, 0.577179)},
1103  {7,10,Point(-0.0330195, 0.0330195, 0.00373942)},
1104  {8,16,Point(0.000886435, 0.000886435, 0.000667355)},
1105  {9,18,Point(0.00355332, 0.000952108, 0.000319882)},
1106  {10,2,Point(-0.0672144, 0.0672144, 1.97043)},
1107  {11,17,Point(-0.000952108, 0.000952108, 0.000782704)}
1108  };
1109 
1110  tabulated_grads[std::make_pair(THIRD, JACOBI_20_00)] =
1111  {
1112  {12,9,Point(0.0328972,0.122774,-0.222472)},
1113  {13,6,Point(-0.184339,0.184339,0.922112)},
1114  {14,20,Point(-0.000523034,0.00195199,0.000121819)},
1115  {15,15,Point(-0.016351,0.00438123,0.0404817)}
1116  };
1117 
1118  tabulated_grads[std::make_pair(FOURTH, JACOBI_20_00)] =
1119  {
1120  {16,3,Point(-0.0381206,0.0102144,0.617059)},
1121  {17,10,Point(0.320895,-0.320895,0.641729)},
1122  {18,14,Point(-0.0246461,0.0246461,-0.300588)},
1123  {19,27,Point(-0.0027324,0.000732145,0.000328402)}
1124  };
1125 
1127  // Arbitrarily selected JACOBI_30_00 reference values.
1129  tabulated_vals[std::make_pair(FIRST, JACOBI_30_00)] = /* JACOBI_30_00 = 13 */
1130  {
1131  // These values are the same as Legendre
1132  {0,9,0.0550016}, {1,5,0.41674}, {2,8,0.0147376}, {3,7,0.111665},
1133  // These values are different from Legendre
1134  {4,6,0.184253}, {5,1,0.200943}, {6,11,0.68764}, {7,15,0.0538426}
1135  };
1136 
1137  tabulated_vals[std::make_pair(SECOND, JACOBI_30_00)] =
1138  {
1139  // These values are the same as Legendre
1140  {0,3,0.0425633}, {1,6,0.0343526}, {2,2,0.158848}, {3,7,0.128206},
1141  // These values are different from Legendre
1142  {4,12,0.552072}, {5,5,-0.211652}, {6,0,0.0372579}, {7,10,-0.0167468},
1143  {8,16,0.00998323}, {9,18,0.0597236}, {10,2,0.139049}, {11,17,0.0160029}
1144  };
1145 
1146  tabulated_vals[std::make_pair(THIRD, JACOBI_30_00)] =
1147  {
1148  {12,9,0.0469849}, {13,6,0.437585}, {14,20,0.0450821}, {15,15,0.0469849}
1149  };
1150 
1151  tabulated_vals[std::make_pair(FOURTH, JACOBI_30_00)] =
1152  {
1153  {16,3,0.0206655}, {17,10,-0.748341}, {18,14,0}, {19,27,0.115995}
1154  };
1155 
1157  // Arbitrarily selected JACOBI_30_00 reference gradients.
1159  tabulated_grads[std::make_pair(FIRST, JACOBI_30_00)] =
1160  {
1161  // These values are the same as Legendre
1162  {0,9,Point(-0.0429458, -0.0115073, 0.)},
1163  {1,5,Point(0.177013, -0.177013, -0.483609)},
1164  {2,8,Point(0.0115073, 0.0115073, 0.00842393)},
1165  {3,7,Point(-0.177013, 0.0474305, 0.)},
1166  // These values are different from Legendre
1167  {4,6,Point(-0.0782626,-0.29208,0.250625)},
1168  {5,1,Point(0.118548,-0.118548,2.18229)},
1169  {6,11,Point(0.143866,0.143866,-0.281275)},
1170  {7,15,Point(-0.00884512,0.00237004,0.00027893)}
1171  };
1172 
1173  tabulated_grads[std::make_pair(SECOND, JACOBI_30_00)] =
1174  {
1175  // These values are the same as Legendre
1176  {0,3,Point(-0.0959817, -0.0959817, 0.0702635)},
1177  {1,6,Point(0.0625228, -0.0625228, 0.0457699)},
1178  {2,2,Point(0.358209, 0.0959817, 0.)},
1179  {3,7,Point(-0.233338, 0.0625228, 0.)},
1180  // These values are different from Legendre
1181  {4,12,Point(-0.0807678,-0.0807678,-0.182443)},
1182  {5,5,Point(0.385213,0.103218,-0.175079)},
1183  {6,0,Point(0.0840179,-0.0225125,0.721474)},
1184  {7,10,Point(-0.0198117,0.0198117,0.0357373)},
1185  {8,16,Point(0.00110804,0.00110804,0.000834194)},
1186  {9,18,Point(0.00662875,0.00177617,0.000482244)},
1187  {10,2,Point(-0.0840179,0.0840179,2.46303)},
1188  {11,17,Point(-0.00177617,0.00177617,0.00142946)},
1189  };
1190 
1191  tabulated_grads[std::make_pair(THIRD, JACOBI_30_00)] =
1192  {
1193  {12,9,Point(0.0184475,0.0688471,-0.254748)},
1194  {13,6,Point(-0.230424,0.230424,1.15264)},
1195  {14,20,Point(-0.000965042,0.00360159,0.000183379)},
1196  {15,15,Point(-0.0423204,0.0113397,0.12514)}
1197  };
1198 
1199  tabulated_grads[std::make_pair(FOURTH, JACOBI_30_00)] =
1200  {
1201  {16,3,Point(-0.0476508,0.012768,0.771323)},
1202  {17,10,Point(0.333487,-0.333487,1.01421)},
1203  {18,14,Point(0,0,0)},
1204  {19,27,Point(-0.00698362,0.00187126,0.000667664)},
1205  };
1206  }
std::map< KeyType, TabulatedGrads > tabulated_grads
std::map< KeyType, TabulatedVals > tabulated_vals
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39

◆ tearDown()

void InfFERadialTest::tearDown ( )
inline

Definition at line 1208 of file inf_fe_radial_test.C.

1208 {}

◆ testDifferentOrders()

void InfFERadialTest::testDifferentOrders ( )
inline

◆ testInfQuants()

void InfFERadialTest::testInfQuants ( )
inline

dphase = r/|r| where r is the vector from the 'origin' to the point of interest Sob. weight = b*b/(r*r) with b being the vector of the elements base to the origin. Since Sobolev weight is 1 for finite elements, it gives a continuous function. its derivative, thus, -2*b*b*r/(r**4)

Definition at line 357 of file inf_fe_radial_test.C.

References libMesh::MeshTools::Generation::build_cube(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::FEGenericBase< OutputType >::build_InfFE(), libMesh::CARTESIAN, libMesh::FEType::default_quadrature_order(), dim, libMesh::MeshBase::elem_ptr(), libMesh::FIRST, libMesh::Elem::infinite(), libMesh::FEInterface::inverse_map(), libMesh::LAGRANGE, libMesh::FEInterface::map(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_elem(), libMesh::TypeVector< T >::norm(), libMesh::TypeVector< T >::norm_sq(), libMesh::Elem::origin(), libMesh::Real, TestCommWorld, libMesh::TET10, and libMesh::TOLERANCE.

358  {
359 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
360  LOG_UNIT_TEST;
361 
362  // std::cout << "Called testSingleOrder with radial_order = "
363  // << radial_order
364  // << std::endl;
365 
368  (mesh,
369  /*nx=*/3, /*ny=*/2, /*nz=*/5,
370  /*xmin=*/0., /*xmax=*/2.,
371  /*ymin=*/-1., /*ymax=*/1.,
372  /*zmin=*/-4., /*zmax=*/-2.,
373  TET10);
374 
378  com_x.first=true;
379  com_y.first=true;
380  com_z.first=true;
381  com_x.second=0.0;
382  com_y.second=-.5;
383  com_z.second=0.0;
384 
385  const unsigned int n_fem =mesh.n_elem();
386 
387  InfElemBuilder builder(mesh);
388  builder.build_inf_elem(com_x, com_y, com_z,
389  false, false, false,
390  true, libmesh_nullptr);
391 
392  // Get pointer to the first infinite Elem.
393  const Elem * infinite_elem = mesh.elem_ptr(n_fem);
394  if (!infinite_elem || !infinite_elem->infinite())
395  libmesh_error_msg("Error setting Elem pointer.");
396 
397  // We will construct FEs, etc. of the same dimension as the mesh elements.
398  auto dim = mesh.mesh_dimension();
399 
400  FEType fe_type(/*Order*/FIRST,
401  /*FEFamily*/LAGRANGE,
402  /*radial order*/FIRST,
403  /*radial_family*/LAGRANGE,
404  /*inf_map*/CARTESIAN);
405 
406  // Construct FE, quadrature rule, etc.
407  std::unique_ptr<FEBase> inf_fe (FEBase::build_InfFE(dim, fe_type));
408  QGauss qrule (dim, fe_type.default_quadrature_order());
409  inf_fe->attach_quadrature_rule(&qrule);
410 
411  const std::vector<Point>& _qpoint = inf_fe->get_xyz();
412  const std::vector<Real> & _weight = inf_fe->get_Sobolev_weight();
413  const std::vector<RealGradient>& dweight = inf_fe->get_Sobolev_dweight();
414  const std::vector<Real>& dxidx = inf_fe->get_dxidx();
415  const std::vector<Real>& dxidy = inf_fe->get_dxidy();
416  const std::vector<Real>& dxidz = inf_fe->get_dxidz();
417  const std::vector<Real>& detadx = inf_fe->get_detadx();
418  const std::vector<Real>& detady = inf_fe->get_detady();
419  const std::vector<Real>& detadz = inf_fe->get_detadz();
420  const std::vector<Real>& dzetadx = inf_fe->get_dzetadx();
421  const std::vector<Real>& dzetady = inf_fe->get_dzetady();
422  const std::vector<Real>& dzetadz = inf_fe->get_dzetadz();
423  const std::vector<RealGradient>& dphase = inf_fe->get_dphase();
424  // Reinit on infinite elem
425  inf_fe->reinit(infinite_elem);
426 
427  for(unsigned int qp =0 ; qp < _weight.size() ; ++qp)
428  {
437  const Point b = base_point(_qpoint[qp], infinite_elem) - infinite_elem->origin();
438  const Point p = _qpoint[qp] - infinite_elem ->origin();
439  const Point rp = FEInterface::inverse_map(3, fe_type, infinite_elem, _qpoint[qp] );
440  const Point again_qp = FEInterface::map(3, fe_type, infinite_elem, rp);
441  //const Point rb = FEInterface::inverse_map(3, fe_type, infinite_elem, b+infinite_elem->origin());
442  const Real v=rp(2);
443 
444  // check that map() does the opposite from inverse_map
445  LIBMESH_ASSERT_FP_EQUAL(again_qp(0), _qpoint[qp](0), TOLERANCE);
446  LIBMESH_ASSERT_FP_EQUAL(again_qp(1), _qpoint[qp](1), TOLERANCE);
447  LIBMESH_ASSERT_FP_EQUAL(again_qp(2), _qpoint[qp](2), TOLERANCE);
448 
449  const Point normal(dzetadx[qp],
450  dzetady[qp],
451  dzetadz[qp]);
452 
453  // check that dweight is in direction of base elements normal
454  // (holds when the base_elem has an affine map)
455  LIBMESH_ASSERT_FP_EQUAL(normal*dweight[qp], -normal.norm()*dweight[qp].norm(), TOLERANCE);
456 
457  if (rp(2) < -.85)
458  // close to the base, dphase is normal
459  LIBMESH_ASSERT_FP_EQUAL(normal*dphase[qp], normal.norm()*dphase[qp].norm(), 3e-4);
460  else if (rp(2) > .85)
461  // 'very far away' dphase goes in radial direction
462  LIBMESH_ASSERT_FP_EQUAL(p*dphase[qp]/p.norm(), dphase[qp].norm(), 1e-4);
463 
464  // check that the mapped radial coordinate corresponds to a/r
465  // - b: radius of the FEM-region (locally, i.e. distance of the projected base point)
466  // - p: distance of the point from origin
467  LIBMESH_ASSERT_FP_EQUAL(b.norm()/p.norm(), (1.-v)/2., TOLERANCE);
468 
469  // this is how weight is defined;
470  LIBMESH_ASSERT_FP_EQUAL(b.norm_sq()/p.norm_sq(), _weight[qp], TOLERANCE);
471 
472  const Point e_xi(dxidx[qp], dxidy[qp], dxidz[qp]);
473  const Point e_eta(detadx[qp], detady[qp], detadz[qp]);
474 
475  LIBMESH_ASSERT_FP_EQUAL(0., b*e_xi,TOLERANCE);
476  LIBMESH_ASSERT_FP_EQUAL(0., b*e_eta,TOLERANCE);
477 
478  }
479 
480 #endif // LIBMESH_ENABLE_INFINITE_ELEMENTS
481  }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:182
The ReplicatedMesh class is derived from the MeshBase class, and is used to store identical copies of...
auto norm() const -> decltype(std::norm(T()))
Definition: type_vector.h:929
virtual Point origin() const
Definition: elem.h:1767
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:159
static constexpr Real TOLERANCE
unsigned int dim
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
auto norm_sq() const -> decltype(std::norm(T()))
Definition: type_vector.h:948
std::pair< bool, double > InfElemOriginValue
Useful typedef.
Point base_point(const Point physical_point, const Elem *inf_elem)
This class is used to build infinite elements on top of an existing mesh.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
This class implements specific orders of Gauss quadrature.
unsigned int mesh_dimension() const
Definition: mesh_base.C:324
virtual bool infinite() const =0
virtual dof_id_type n_elem() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.

◆ testInfQuants_numericDeriv()

void InfFERadialTest::testInfQuants_numericDeriv ( )
inline

Definition at line 483 of file inf_fe_radial_test.C.

References std::abs(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::FEGenericBase< OutputType >::build_InfFE(), libMesh::MeshTools::Generation::build_sphere(), libMesh::CARTESIAN, dim, libMesh::MeshBase::elem_ptr(), libMesh::FIRST, libMesh::HEX27, libMesh::Elem::infinite(), libMesh::LAGRANGE, mesh, libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_elem(), libMesh::Elem::node_ptr(), libMesh::TensorTools::norm(), libMesh::TypeVector< T >::norm(), libMesh::TensorTools::norm_sq(), libMesh::TypeVector< T >::norm_sq(), libMesh::Elem::origin(), libMesh::Real, std::sqrt(), TestCommWorld, and libMesh::TOLERANCE.

484  {
485  // We use AMR internally to build_sphere
486 #if defined(LIBMESH_ENABLE_INFINITE_ELEMENTS) && defined(LIBMESH_ENABLE_AMR)
487  LOG_UNIT_TEST;
488 
491  (mesh, /*rad*/ 1,
492  /* nr */ 2, /*type*/ HEX27);
493 
497  com_x.first=true;
498  com_y.first=true;
499  com_z.first=true;
500  com_x.second=0.7;
501  com_y.second=-0.4;
502  com_z.second=0.0;
503 
504  const unsigned int n_fem =mesh.n_elem();
505 
506  InfElemBuilder builder(mesh);
507  builder.build_inf_elem(com_x, com_y, com_z,
508  false, false, false,
509  true, libmesh_nullptr);
510 
511  // Get pointer to the first infinite Elem.
512  Elem * infinite_elem = mesh.elem_ptr(n_fem+1);
513  if (!infinite_elem || !infinite_elem->infinite())
514  libmesh_error_msg("Error setting Elem pointer.");
515  // lets overemphasize that the base element has a non-affine map:
516  for (unsigned int n=8; n<12; ++n)
517  {
518  Node* node=infinite_elem->node_ptr(n);
519  // shift base-points on edges towards the center.
520  // This leads to strong 'non-affine effects'
521  *node -= 0.1*(*node-infinite_elem->origin()).unit();
522  }
523 
524  // We will construct FEs, etc. of the same dimension as the mesh elements.
525  auto dim = mesh.mesh_dimension();
526 
527  FEType fe_type(/*Order*/FIRST,
528  /*FEFamily*/LAGRANGE,
529  /*radial order*/FIRST,
530  /*radial_family*/LAGRANGE,
531  /*inf_map*/CARTESIAN);
532 
533  // Reinit on infinite elem
534  //
535  unsigned int num_pt=10;
536  std::vector<Point> points(2*num_pt);
537  points[0]=Point(-0.7, -0.5, -0.9);
538  points[1]=Point(-0.1, 0.9, -0.9);
539  points[2]=Point(-0.7, -0.5, -0.4);
540  points[3]=Point(-0.1, 0.9, -0.4);
541  points[4]=Point(-0.7, -0.5, -0.2);
542  points[5]=Point(-0.1, 0.9, -0.2);
543  points[6]=Point(-0.7, -0.5, 0.1);
544  points[7]=Point(-0.1, 0.9, 0.1);
545  points[8]=Point(-0.7, -0.5, 0.6);
546  points[9]=Point(-0.1, 0.9, 0.6);
547 
548  // Check derivatives along radial direction:
549  Point delta(0.,0.,1e-3);
550  for (unsigned int i=num_pt; i<2*num_pt; ++i)
551  points[i]=points[i-num_pt]+delta;
552 
553  std::unique_ptr<FEBase> inf_fe (FEBase::build_InfFE(dim, fe_type));
554  const std::vector<Point> & q_point = inf_fe->get_xyz();
555  const std::vector<Real> & sob_w = inf_fe->get_Sobolev_weightxR_sq();
556  const std::vector<RealGradient> & dsob_w = inf_fe->get_Sobolev_dweightxR_sq();
557  const std::vector<Real> & sob_now = inf_fe->get_Sobolev_weight();
558  const std::vector<RealGradient>& dsob_now = inf_fe->get_Sobolev_dweight();
559  const std::vector<RealGradient>& dphase = inf_fe->get_dphase();
560  const std::vector<std::vector<RealGradient> >& dphi = inf_fe->get_dphi();
561  const std::vector<std::vector<Real> >& phi = inf_fe->get_phi();
562  const std::vector<std::vector<RealGradient> >& dphi_w = inf_fe->get_dphi_over_decayxR();
563  const std::vector<std::vector<Real> >& phi_w = inf_fe->get_phi_over_decayxR();
564  inf_fe->reinit(infinite_elem,&points);
565 
566  // check that all quantities are sized correctly:
567  libmesh_assert_equal_to(q_point.size(), sob_w.size());
568  libmesh_assert_equal_to(q_point.size(), dsob_w.size());
569  libmesh_assert_equal_to(q_point.size(), sob_now.size());
570  libmesh_assert_equal_to(q_point.size(), dsob_now.size());
571  libmesh_assert_equal_to(q_point.size(), dphase.size());
572  libmesh_assert_equal_to(phi.size(), phi_w.size());
573  libmesh_assert_equal_to(phi.size(), dphi.size());
574  libmesh_assert_equal_to(phi.size(), dphi_w.size());
575  libmesh_assert_equal_to(q_point.size(), phi[0].size());
576  libmesh_assert_equal_to(q_point.size(), phi_w[0].size());
577  libmesh_assert_equal_to(q_point.size(), dphi[0].size());
578  libmesh_assert_equal_to(q_point.size(), dphi_w[0].size());
579 
580  for(unsigned int qp =0 ; qp < num_pt ; ++qp)
581  {
582  const Point dxyz(q_point[qp+num_pt]-q_point[qp]);
583  const Point b_i= base_point(q_point[qp], infinite_elem) - infinite_elem->origin();
584  const Point b_o= base_point(q_point[qp+num_pt], infinite_elem) - infinite_elem->origin();
585 
586  LIBMESH_ASSERT_FP_EQUAL((b_i-b_o).norm_sq(), 0, TOLERANCE);
587 
588  Real weight_i = b_i.norm_sq()/(q_point[qp]-infinite_elem->origin()).norm_sq();
589  Real weight_o = b_o.norm_sq()/(q_point[qp+num_pt]-infinite_elem->origin()).norm_sq();
590  const Real phase_i = (q_point[qp]-infinite_elem->origin()).norm() - b_i.norm();
591  const Real phase_o = (q_point[qp+num_pt]-infinite_elem->origin()).norm() - b_o.norm();
592 
593  Real tolerance = std::abs((dphase[qp+num_pt]-dphase[qp])*dxyz)+1e-10;
594  Real deriv_mean = (dphase[qp]*dxyz + dphase[qp+num_pt]*dxyz)*0.5;
595  LIBMESH_ASSERT_FP_EQUAL(phase_o - phase_i, deriv_mean, tolerance*.5);
596 
597  tolerance = std::abs((dsob_now[qp+num_pt]-dsob_now[qp])*dxyz)+1e-10;
598  deriv_mean = (dsob_now[qp]*dxyz + dsob_now[qp+num_pt]*dxyz)* 0.5;
599  LIBMESH_ASSERT_FP_EQUAL(sob_now[qp+num_pt] - sob_now[qp], deriv_mean, tolerance*.5);
600 
601  LIBMESH_ASSERT_FP_EQUAL(sob_w[qp+num_pt]*weight_o - sob_w[qp]*weight_i, dsob_w[qp]*dxyz*weight_i, tolerance);
602 
603  for (unsigned int i=0; i< phi.size(); ++i)
604  {
605  tolerance = std::abs((dphi[i][qp+num_pt]-dphi[i][qp])*dxyz)+1e-10;
606  deriv_mean = (dphi[i][qp]*dxyz + dphi[i][qp+num_pt]*dxyz)*0.5;
607  LIBMESH_ASSERT_FP_EQUAL(phi[i][qp+num_pt] - phi[i][qp], deriv_mean, tolerance*.5);
608 
609  deriv_mean = 0.5*(dphi_w[i][qp]*dxyz*sqrt(weight_i) +dphi_w[i][qp+num_pt]*dxyz*sqrt(weight_o));
610  LIBMESH_ASSERT_FP_EQUAL(phi_w[i][qp+num_pt]*sqrt(weight_o) - phi_w[i][qp]*sqrt(weight_i),
611  deriv_mean, tolerance*.5);
612  }
613 
614  }
615 
616  // Check derivatives along angular direction:
617  points[0 ]=Point(-0.7, -0.5, -0.9);
618  points[2 ]=Point(-0.1, 0.9, -0.9);
619  points[4 ]=Point(-0.7, -0.5, -0.4);
620  points[6 ]=Point(-0.1, 0.9, -0.4);
621  points[8 ]=Point(-0.7, -0.5, -0.2);
622  points[10]=Point(-0.1, 0.9, -0.2);
623  points[12]=Point(-0.7, -0.5, 0.1);
624  points[14]=Point(-0.1, 0.9, 0.1);
625  points[16]=Point(-0.7, -0.5, 0.6);
626  points[18]=Point(-0.1, 0.9, 0.6);
627 
628  delta = Point(1.2e-4,-2.7e-4,0.);
629  for (unsigned int i=0; i<2*num_pt; i+=2)
630  points[i+1]=points[i]+delta;
631 
632  const std::vector<Real>& dzetadx = inf_fe->get_dzetadx();
633  const std::vector<Real>& dzetady = inf_fe->get_dzetady();
634  const std::vector<Real>& dzetadz = inf_fe->get_dzetadz();
635 
636  inf_fe->reinit(infinite_elem,&points);
637 
638  for(unsigned int qp =0 ; qp < 2*num_pt ; qp+=2)
639  {
640  const Point dxyz(q_point[qp+1]-q_point[qp]);
641  const Point b_i= base_point(q_point[qp], infinite_elem) - infinite_elem->origin();
642  const Point b_o= base_point(q_point[qp+1], infinite_elem) - infinite_elem->origin();
643  Real weight_i = b_i.norm_sq()/(q_point[qp]-infinite_elem->origin()).norm_sq();
644  Real weight_o = b_o.norm_sq()/(q_point[qp+1]-infinite_elem->origin()).norm_sq();
645  const Real phase_i = (q_point[qp]-infinite_elem->origin()).norm() - b_i.norm();
646  const Real phase_o = (q_point[qp+1]-infinite_elem->origin()).norm() - b_o.norm();
647 
648  LIBMESH_ASSERT_FP_EQUAL(weight_o ,weight_i, TOLERANCE);
649 
650  Point normal(dzetadx[qp],
651  dzetady[qp],
652  dzetadz[qp]);
653 
654  Real a_i= (q_point[qp]-infinite_elem->origin()).norm()*0.5*(1.-points[qp](2));
655 
656  LIBMESH_ASSERT_FP_EQUAL(a_i, b_i.norm(), TOLERANCE*TOLERANCE);
657 
658  // the elements normal should be orthogonal to dxyz, but in practice deviates
659  // approx. 1 degree. So we must account for this error in direction.
660  Real err_direct = 1.5*std::abs(dxyz*normal)/(dxyz.norm()*normal.norm());
661 
662  Real tolerance = std::abs((dphase[qp+1]-dphase[qp])*dxyz)*0.5 + err_direct*dxyz.norm()*dphase[qp].norm();
663  Real deriv_mean = (dphase[qp] + dphase[qp+1])*dxyz*0.5;
664  LIBMESH_ASSERT_FP_EQUAL(phase_o - phase_i, deriv_mean, tolerance );
665 
666  LIBMESH_ASSERT_FP_EQUAL(normal.cross(dsob_now[qp]).norm_sq(), 0.0, TOLERANCE*TOLERANCE);
667 
668  tolerance = std::abs((dsob_now[qp+1]-dsob_now[qp])*dxyz)*.5+1e-10 + err_direct*dxyz.norm()*dsob_now[qp].norm();
669  deriv_mean = (dsob_now[qp]*dxyz + dsob_now[qp+1]*dxyz)*0.5;
670  LIBMESH_ASSERT_FP_EQUAL(sob_now[qp+1] - sob_now[qp], deriv_mean, tolerance);
671 
672  deriv_mean = (dsob_w[qp]*dxyz*weight_i +dsob_w[qp+1]*dxyz*weight_o)*0.5;
673  LIBMESH_ASSERT_FP_EQUAL(sob_w[qp+1]*weight_o - sob_w[qp]*weight_i, deriv_mean, tolerance);
674 
675  for (unsigned int i=0; i< phi.size(); ++i)
676  {
677  tolerance = std::abs((dphi[i][qp+1]-dphi[i][qp])*dxyz)*0.5+1e-10 + err_direct*dxyz.norm()*dphi[i][qp].norm();
678  deriv_mean = (dphi[i][qp]*dxyz + dphi[i][qp+1]*dxyz)*.5;
679  LIBMESH_ASSERT_FP_EQUAL(phi[i][qp+1] - phi[i][qp], deriv_mean, tolerance);
680  }
681  }
682 
683 #endif // LIBMESH_ENABLE_INFINITE_ELEMENTS && LIBMESH_ENABLE_AMR
684  }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:182
The ReplicatedMesh class is derived from the MeshBase class, and is used to store identical copies of...
void build_sphere(UnstructuredMesh &mesh, const Real rad=1, const unsigned int nr=2, const ElemType type=INVALID_ELEM, const unsigned int n_smooth=2, const bool flat=true)
Meshes a spherical or mapped-spherical domain.
A Node is like a Point, but with more information.
Definition: node.h:52
auto norm() const -> decltype(std::norm(T()))
Definition: type_vector.h:929
virtual Point origin() const
Definition: elem.h:1767
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:159
static constexpr Real TOLERANCE
unsigned int dim
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:53
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:57
auto norm_sq() const -> decltype(std::norm(T()))
Definition: type_vector.h:948
std::pair< bool, double > InfElemOriginValue
Useful typedef.
auto norm(const T &a) -> decltype(std::abs(a))
Definition: tensor_tools.h:74
Point base_point(const Point physical_point, const Elem *inf_elem)
This class is used to build infinite elements on top of an existing mesh.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Node * node_ptr(const unsigned int i) const
Definition: elem.h:2331
auto norm_sq(const T &a) -> decltype(std::norm(a))
Definition: tensor_tools.h:104
unsigned int mesh_dimension() const
Definition: mesh_base.C:324
virtual bool infinite() const =0
virtual dof_id_type n_elem() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39

◆ testRefinement()

void InfFERadialTest::testRefinement ( )
inline

Definition at line 813 of file inf_fe_radial_test.C.

References libMesh::MeshBase::add_elem(), libMesh::MeshBase::add_point(), libMesh::EquationSystems::add_system(), libMesh::System::add_variable(), libMesh::MeshBase::all_second_order(), libMesh::FEGenericBase< OutputType >::build(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::FEGenericBase< OutputType >::build_InfFE(), libMesh::Elem::build_with_id(), libMesh::Elem::child_ptr(), libMesh::DofMap::constrain_element_vector(), libMesh::DofMap::dof_indices(), libMesh::MeshBase::elem_ptr(), libMesh::PointLocatorTree::enable_out_of_mesh_mode(), libMesh::Elem::find_point_neighbors(), libMesh::FIRST, libMesh::System::get_dof_map(), libMesh::FEGenericBase< OutputType >::get_phi(), libMesh::HEX8, libMesh::DofObject::id(), libMesh::Elem::infinite(), libMesh::DofObject::invalid_id, libMesh::FEInterface::inverse_map(), libMesh::LAGRANGE, mesh, libMesh::NINTH, libMesh::MeshBase::node_ptr(), libMesh::MeshBase::prepare_for_use(), libMesh::FEType::radial_order, std::real(), libMesh::Elem::REFINE, libMesh::MeshRefinement::refine_elements(), libMesh::EquationSystems::reinit(), libMesh::FEAbstract::reinit(), libMesh::MeshBase::set_mesh_dimension(), libMesh::Elem::set_node(), libMesh::Elem::set_refinement_flag(), libMesh::MeshBase::set_spatial_dimension(), TestCommWorld, and libMesh::TOLERANCE.

814  {
815 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
816  LOG_UNIT_TEST;
817 
818  // 1. create a mesh with one finite element, add an infinite element at each side
820 
821  {
824 
825  // we'll create a slightly deformed hex:
826  mesh.add_point(Point(-0.7,-0.8, 0.7),0);
827  mesh.add_point(Point(-0.8, 0.8, 0.8),1);
828  mesh.add_point(Point( 0.8, 0.8, 0.7),2);
829  mesh.add_point(Point( 0.8,-0.8, 0.8),3);
830  mesh.add_point(Point(-0.8,-0.8,-0.8),4);
831  mesh.add_point(Point(-0.8, 0.8,-0.7),5);
832  mesh.add_point(Point( 0.8, 0.7,-0.8),6);
833  mesh.add_point(Point( 0.8,-0.8,-0.8),7);
834  Elem * elem = mesh.add_elem(Elem::build_with_id(HEX8, 0));
835  for (unsigned int i=0; i< 8; ++i)
836  elem->set_node(i) = mesh.node_ptr(i);
837 
838  InfElemBuilder builder(mesh);
839  builder.build_inf_elem();
840 
843  }
844 
845  EquationSystems es(mesh);
846  FEType fe_type(FIRST, LAGRANGE);
847  // set some high radial order:
848  fe_type.radial_order=NINTH;
849  LinearImplicitSystem & dummy = es.add_system<LinearImplicitSystem> ("dummy");
850  dummy.add_variable("var", fe_type);
851 
852  // 2. refine the finite element !!!
853  {
854  MeshRefinement mesh_refinement(mesh);
855  // refine one of the elements:
856  mesh.elem_ptr(2)->set_refinement_flag(Elem::REFINE);
857  mesh_refinement.refine_elements();
858  // having done this: refine its 1-st child as well to get two refinement-levels.
859  mesh.elem_ptr(2)->child_ptr(1)->set_refinement_flag(Elem::REFINE);
860  // Here, due to sanity-checking, the neighbors of the doubly-refined element
861  // are refined to at least first level.
862  mesh_refinement.refine_elements();
863  }
864 
865  // 3. go through all DOFs and check that they are continuous between elements at nodes
866  {
867  // here, the constraints are computed.
868  es.reinit();
869 
870  PointLocatorTree pt_lctr(mesh);
871  pt_lctr.enable_out_of_mesh_mode();
872 
873  std::vector<dof_id_type> dof_indices;
874 
875  std::unique_ptr<FEBase> fe (FEBase::build(3, fe_type));
876  std::unique_ptr<FEBase> inf_fe (FEBase::build_InfFE(3, fe_type));
877 
878  const DofMap& dof_map = dummy.get_dof_map();
879 
880  const auto n_end = mesh.nodes_end();
881  auto n_it = mesh.nodes_begin();
882 
883  // loop over all nodes: We will check that for all nodes each DOF has the same value on all elements.
884  for(; n_it != n_end; ++n_it)
885  {
886  const Node* node = *n_it;
887  // check that DOF-values coincide on all elements that contain this node:
888  // For this, we first search one element containing node as reference element
889  // and than loop over all other elements that contain that node
890  const Elem * some_elem = pt_lctr(*node);
891 
892  dof_map.dof_indices (some_elem, dof_indices);
893 
894  FEBase * cfe = libmesh_nullptr;
895 
896  if (some_elem->infinite())
897  cfe = inf_fe.get();
898  else
899  cfe = fe.get();
900 
901  const std::vector<std::vector<Real> >& phi = cfe->get_phi();
902 
903  const Point node_internal = FEInterface::inverse_map(3, fe_type, some_elem, *node);
904  std::vector<Point> eval_pt(1, node_internal);
905  cfe->reinit(some_elem, &eval_pt);
906 
907  // safe the values phi[i] of some_elem in a vector.
908  // Than we will test if the values are the same at the other elements.
909  DenseVector<Number> reference_values(dof_indices.size());
910  std::vector<dof_id_type> reference_indices(dof_indices.size());
911  for (unsigned int i=0; i< dof_indices.size(); ++i)
912  {
913  reference_indices[i]=dof_indices[i];
914  reference_values(i)=phi[i][0];
915  }
916  dof_map.constrain_element_vector(reference_values, reference_indices);
917 
918  std::set<const libMesh::Elem *> elements;
919  some_elem->find_point_neighbors(*node, elements);
920 
921  for (auto const_elem : elements)
922  {
923  if (some_elem->id() == const_elem->id())
924  continue;
925 
926  // we need a non-const pointer, so we just create a second pointer:
927  Elem* elem = mesh.elem_ptr(const_elem->id());
928  dof_map.dof_indices (elem, dof_indices);
929  if (elem->infinite())
930  cfe = inf_fe.get();
931  else
932  cfe = fe.get();
933 
934  const std::vector<std::vector<Real> >& phi1 = cfe->get_phi();
935  eval_pt[0]=FEInterface::inverse_map(3, fe_type, elem, *node);
936 
937  cfe->reinit(elem, &eval_pt);
938 
939  DenseVector<Number> phi_values(phi1.size());
940  for( unsigned int i=0; i < phi_values.size(); ++i)
941  phi_values(i) = phi1[i][0];
942 
943  dof_map.constrain_element_vector(phi_values, dof_indices);
944 
945  for(unsigned int i=0; i< dof_indices.size(); ++i)
946  {
947  dof_id_type glob_ind=DofObject::invalid_id;
948  for(unsigned int gi=0; gi < reference_indices.size(); ++gi)
949  {
950  if (reference_indices[gi] == dof_indices[i])
951  {
952  glob_ind = gi;
953  break;
954  }
955  }
956  // If element was found check that the values agree
957  if (glob_ind != DofObject::invalid_id)
958  LIBMESH_ASSERT_FP_EQUAL(std::real(phi_values(i)), std::real(reference_values(glob_ind)), TOLERANCE);
959  }
960 
961  }
962 
963  }
964 
965  }
966 #endif
967  }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:182
This is the EquationSystems class.
The ReplicatedMesh class is derived from the MeshBase class, and is used to store identical copies of...
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:2381
A Node is like a Point, but with more information.
Definition: node.h:52
boost::multiprecision::float128 real(const boost::multiprecision::float128 in)
This is a point locator.
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1992
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:159
static constexpr Real TOLERANCE
Manages consistently variables, degrees of freedom, coefficient vectors, matrices and linear solvers ...
void set_spatial_dimension(unsigned char d)
Sets the "spatial dimension" of the Mesh.
Definition: mesh_base.C:502
void prepare_for_use(const bool skip_renumber_nodes_and_elements, const bool skip_find_neighbors)
Prepare a newly ecreated (or read) mesh for use.
Definition: mesh_base.C:710
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
void set_refinement_flag(const RefinementState rflag)
Sets the value of the refinement flag for the element.
Definition: elem.h:3055
virtual Node * add_point(const Point &p, const dof_id_type id=DofObject::invalid_id, const processor_id_type proc_id=DofObject::invalid_processor_id)=0
Add a new Node at Point p to the end of the vertex array, with processor_id procid.
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:169
Implements (adaptive) mesh refinement algorithms for a MeshBase.
dof_id_type id() const
Definition: dof_object.h:823
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
void find_point_neighbors(const Point &p, std::set< const Elem *> &neighbor_set) const
This function finds all active elements (including this one) which are in the same manifold as this e...
Definition: elem.C:801
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
This is at the core of this class.
unsigned int add_variable(std::string_view var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=nullptr)
Adds the variable var to the list of variables for this system.
Definition: system.C:1305
void set_mesh_dimension(unsigned char d)
Resets the logical dimension of the mesh.
Definition: mesh_base.h:269
This class is used to build infinite elements on top of an existing mesh.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
void constrain_element_vector(DenseVector< Number > &rhs, std::vector< dof_id_type > &dofs, bool asymmetric_constraint_rows=true) const
Constrains the element vector.
Definition: dof_map.h:2250
void all_second_order(const bool full_ordered=true)
Calls the range-based version of this function with a range consisting of all elements in the mesh...
Definition: mesh_base.C:1535
virtual bool infinite() const =0
virtual const Node * node_ptr(const dof_id_type i) const =0
const DofMap & get_dof_map() const
Definition: system.h:2293
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
This class forms the foundation from which generic finite elements may be derived.
const Elem * child_ptr(unsigned int i) const
Definition: elem.h:3000
uint8_t dof_id_type
Definition: id_types.h:67
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:207

◆ testSides()

void InfFERadialTest::testSides ( )
inline

Definition at line 221 of file inf_fe_radial_test.C.

References libMesh::FEGenericBase< OutputType >::build(), libMesh::MeshTools::Generation::build_cube(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::FEGenericBase< OutputType >::build_InfFE(), libMesh::CARTESIAN, libMesh::FEType::default_quadrature_order(), dim, libMesh::MeshBase::elem_ptr(), libMesh::FIRST, libMesh::Elem::infinite(), libMesh::LAGRANGE, libMesh::libmesh_assert(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_elem(), libMesh::Elem::n_sides(), libMesh::Elem::neighbor_ptr(), libMesh::Elem::node_ptr(), libMesh::PRISM18, libMesh::SECOND, libMesh::Elem::set_node(), libMesh::Elem::side_ptr(), TestCommWorld, and libMesh::TOLERANCE.

222  {
223 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
224  LOG_UNIT_TEST;
225 
228  (mesh,
229  /*nx=*/2, /*ny=*/2, /*nz=*/2,
230  /*xmin=*/-2., /*xmax=*/2.,
231  /*ymin=*/-1., /*ymax=*/1.,
232  /*zmin=*/-2., /*zmax=*/2.,
233  PRISM18);
234 
235  const unsigned int n_fem =mesh.n_elem();
236 
237  InfElemBuilder builder(mesh);
238  builder.build_inf_elem();
239 
240  // Get pointer to the first infinite Elem.
241  Elem * infinite_elem = mesh.elem_ptr(n_fem+3);
242  if (!infinite_elem || !infinite_elem->infinite())
243  libmesh_error_msg("Error setting Elem pointer.");
244 
245  // We will construct FEs, etc. of the same dimension as the mesh elements.
246  auto dim = mesh.mesh_dimension();
247 
248  //FEType fe_type(/*Order*/FIRST,
249  FEType fe_type(/*Order*/SECOND,
250  /*FEFamily*/LAGRANGE,
251  /*radial order*/FIRST,
252  /*radial_family*/LAGRANGE,
253  /*inf_map*/CARTESIAN);
254 
255 
256  // Construct FE, quadrature rule, etc.
257  std::unique_ptr<FEBase> inf_fe (FEBase::build_InfFE(dim, fe_type));
258  std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
259  std::unique_ptr<FEBase> fe2 (FEBase::build(dim-1, fe_type));
260  QGauss qrule (dim-1, fe_type.default_quadrature_order());
261  inf_fe->attach_quadrature_rule(&qrule);
262  fe->attach_quadrature_rule(&qrule);
263  fe2->attach_quadrature_rule(&qrule);
264 
265  //find the FE neighbour of \p infinite_elem:
266  unsigned int side_num=7;
267  MeshBase::const_element_iterator el = mesh.elements_begin();
268  const MeshBase::const_element_iterator end_el = mesh.elements_end();
269  Point side_pt = infinite_elem->side_ptr(0)->vertex_average();
270 
271  Elem * finite_elem = infinite_elem->neighbor_ptr(0);
272  for(unsigned int i=0; i<finite_elem->n_sides(); ++i)
273  {
274  if (finite_elem->side_ptr(i)->contains_point(side_pt))
275  {
276  side_num = i;
277  break;
278  }
279  }
280  libmesh_assert(side_num < 7);
281  // in particular, side_num = 0 (used later)
282 
283  // construct a 'side' element (specific to the current geometry):
284  Tri6 side_elem;
285  side_elem.set_node(0)=finite_elem->node_ptr(0);
286  side_elem.set_node(1)=finite_elem->node_ptr(2);
287  side_elem.set_node(2)=finite_elem->node_ptr(1);
288  side_elem.set_node(3)=finite_elem->node_ptr(8);
289  side_elem.set_node(4)=finite_elem->node_ptr(7);
290  side_elem.set_node(5)=finite_elem->node_ptr(6);
291 
292  const std::vector<Point>& i_qpoint = inf_fe->get_xyz();
293  const std::vector<Real> & i_weight = inf_fe->get_Sobolev_weight();
294  const std::vector<Real> & i_JxW = inf_fe->get_JxWxdecay_sq();
295  const std::vector<std::vector<Real> >& i_phi = inf_fe->get_phi();
296  const std::vector<Point> & i_normal = inf_fe->get_normals();
297  const std::vector<std::vector<Point> >& i_tangents = inf_fe->get_tangents();
298 
299  const std::vector<Point>& f_qpoint = fe->get_xyz();
300  const std::vector<Real> & f_weight = fe->get_Sobolev_weight();
301  const std::vector<Real> & f_JxW = fe->get_JxWxdecay_sq();
302  const std::vector<std::vector<Real> >& f_phi = fe->get_phi();
303  const std::vector<Point> & f_normal = fe->get_normals();
304  const std::vector<std::vector<Point> >& f_tangents = fe->get_tangents();
305 
306  const std::vector<Point>& s_qpoint = fe2->get_xyz();
307  const std::vector<Real> & s_JxW = fe2->get_JxWxdecay_sq();
308  const std::vector<std::vector<Real> >& s_phi = fe2->get_phi();
309 
310  // Reinit on infinite elem
311  inf_fe->reinit(infinite_elem, (unsigned)0);
312  fe->reinit(finite_elem, side_num);
313  fe2->reinit(&side_elem);
314 
315  LIBMESH_ASSERT_FP_EQUAL(i_weight.size(), f_weight.size(), TOLERANCE);
316  for(unsigned int qp =0 ; qp < i_weight.size() ; ++qp)
317  {
318  LIBMESH_ASSERT_FP_EQUAL(i_qpoint[qp](0), f_qpoint[qp](0), TOLERANCE);
319  LIBMESH_ASSERT_FP_EQUAL(i_qpoint[qp](1), f_qpoint[qp](1), TOLERANCE);
320  LIBMESH_ASSERT_FP_EQUAL(i_qpoint[qp](2), f_qpoint[qp](2), TOLERANCE);
321 
322  LIBMESH_ASSERT_FP_EQUAL(s_qpoint[qp](0), f_qpoint[qp](0), TOLERANCE);
323  LIBMESH_ASSERT_FP_EQUAL(s_qpoint[qp](1), f_qpoint[qp](1), TOLERANCE);
324  LIBMESH_ASSERT_FP_EQUAL(s_qpoint[qp](2), f_qpoint[qp](2), TOLERANCE);
325 
326  LIBMESH_ASSERT_FP_EQUAL(i_weight[qp], f_weight[qp], TOLERANCE); // this is both 1
327  LIBMESH_ASSERT_FP_EQUAL(i_JxW[qp] , f_JxW[qp] , TOLERANCE);
328  LIBMESH_ASSERT_FP_EQUAL(s_JxW[qp] , f_JxW[qp] , TOLERANCE);
329 
330  // initialize index maps: (specific to current geometry:
331  unsigned int inf_index[] ={0, 1, 2, 6, 7, 8};
332  unsigned int fe_index[] ={0, 2, 1, 8, 7, 6};
333  unsigned int s_index[] ={0, 1, 2, 3, 4, 5};
334  for (unsigned int i=0; i<6; ++i)
335  {
336  LIBMESH_ASSERT_FP_EQUAL(i_phi[inf_index[i]][qp], f_phi[fe_index[i]][qp], TOLERANCE);
337  LIBMESH_ASSERT_FP_EQUAL(i_phi[inf_index[i]][qp], s_phi[s_index[i]][qp], TOLERANCE);
338  }
339 
340  // note that infinite element normals point OUTWARD of the element
341  LIBMESH_ASSERT_FP_EQUAL(i_normal[qp](0), f_normal[qp](0), TOLERANCE);
342  LIBMESH_ASSERT_FP_EQUAL(i_normal[qp](1), f_normal[qp](1), TOLERANCE);
343  LIBMESH_ASSERT_FP_EQUAL(i_normal[qp](2), f_normal[qp](2), TOLERANCE);
344 
345  for (unsigned int i=0; i< i_tangents[0].size(); ++i)
346  {
347  LIBMESH_ASSERT_FP_EQUAL(i_tangents[qp][i](0), f_tangents[qp][i](0), TOLERANCE);
348  LIBMESH_ASSERT_FP_EQUAL(i_tangents[qp][i](1), f_tangents[qp][i](1), TOLERANCE);
349  LIBMESH_ASSERT_FP_EQUAL(i_tangents[qp][i](2), f_tangents[qp][i](2), TOLERANCE);
350  }
351 
352  }
353 
354 #endif
355  }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:182
The ReplicatedMesh class is derived from the MeshBase class, and is used to store identical copies of...
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:2381
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:159
The definition of the const_element_iterator struct.
Definition: mesh_base.h:2121
static constexpr Real TOLERANCE
unsigned int dim
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
The Tri6 is an element in 2D composed of 6 nodes.
Definition: face_tri6.h:61
libmesh_assert(ctx)
This class is used to build infinite elements on top of an existing mesh.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
virtual unsigned int n_sides() const =0
const Elem * neighbor_ptr(unsigned int i) const
Definition: elem.h:2407
virtual std::unique_ptr< Elem > side_ptr(unsigned int i)=0
const Node * node_ptr(const unsigned int i) const
Definition: elem.h:2331
This class implements specific orders of Gauss quadrature.
unsigned int mesh_dimension() const
Definition: mesh_base.C:324
virtual bool infinite() const =0
virtual dof_id_type n_elem() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.

◆ testSingleOrder()

void InfFERadialTest::testSingleOrder ( Order  radial_order,
FEFamily  radial_family 
)
inline

Definition at line 706 of file inf_fe_radial_test.C.

References libMesh::MeshTools::Generation::build_cube(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::FEGenericBase< OutputType >::build_InfFE(), libMesh::CARTESIAN, libMesh::FEType::default_quadrature_order(), dim, libMesh::MeshBase::elem_ptr(), libMesh::FIRST, libMesh::HEX8, libMesh::Elem::infinite(), libMesh::LAGRANGE, libMesh::libmesh_ignore(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_elem(), TestCommWorld, and libMesh::TOLERANCE.

708  {
709  // Avoid warnings when not compiled with infinite element support.
710  libmesh_ignore(radial_order, radial_family);
711 
712 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
713  // std::cout << "Called testSingleOrder with radial_order = "
714  // << radial_order
715  // << " radial_family: "<<radial_family
716  // << std::endl;
717 
720  (mesh,
721  /*nx=*/1, /*ny=*/1, /*nz=*/1,
722  /*xmin=*/-1., /*xmax=*/1.,
723  /*ymin=*/-1., /*ymax=*/1.,
724  /*zmin=*/-1., /*zmax=*/1.,
725  HEX8);
726 
727  // Add infinite elements to the mesh. The optional verbose flag only
728  // prints extra information in non-optimized mode.
729  InfElemBuilder builder(mesh);
730  builder.build_inf_elem(/*verbose=true*/);
731 
732  // Get pointer to the last infinite Elem. I originally intended to
733  // get a pointer to the *first* infinite Elem but then I computed
734  // all the reference values on the last Elem, so that's the one we
735  // are using now...
736  const Elem * infinite_elem = mesh.elem_ptr(mesh.n_elem() - 1);
737  libmesh_error_msg_if(!infinite_elem || !infinite_elem->infinite(),
738  "Error setting Elem pointer.");
739 
740  // We will construct FEs, etc. of the same dimension as the mesh elements.
741  auto dim = mesh.mesh_dimension();
742 
743  FEType fe_type(/*Order*/FIRST,
744  /*FEFamily*/LAGRANGE,
745  radial_order,
746  radial_family,
747  /*inf_map*/CARTESIAN);
748 
749  // Construct FE, quadrature rule, etc.
750  std::unique_ptr<FEBase> inf_fe (FEBase::build_InfFE(dim, fe_type));
751  QGauss qrule (dim, fe_type.default_quadrature_order());
752  inf_fe->attach_quadrature_rule(&qrule);
753  const std::vector<std::vector<Real>> & phi = inf_fe->get_phi();
754  const std::vector<std::vector<RealGradient>> & dphi = inf_fe->get_dphi();
755 
756  // Reinit on infinite elem
757  inf_fe->reinit(infinite_elem);
758 
759  // Check phi values against reference values. The number of
760  // quadrature points and shape functions seem to both increase by
761  // 4 as the radial_order is increased.
762  //
763  // radial_order | n_qp | n_sf
764  // --------------------------
765  // FIRST | 16 | 8
766  // SECOND | 20 | 12
767  // THIRD | 24 | 16
768  // FOURTH | 28 | 20
769  // FIFTH | 32 | 24
770 
771  // Debugging print values so we can tabulate them
772  // auto n_qp = inf_fe->n_quadrature_points();
773  // auto n_sf = inf_fe->n_shape_functions();
774  // for (unsigned int i=0; i<n_sf; ++i)
775  // for (unsigned qp=0; qp<n_qp; ++qp)
776  // {
777  // libMesh::out << "{" << i << "," << qp << "," << phi[i][qp] << "}," << std::endl;
778  // }
779  // for (unsigned int i=0; i<n_sf; ++i)
780  // for (unsigned qp=0; qp<n_qp; ++qp)
781  // {
782  // libMesh::out << "{" << i << "," << qp << ",Point("
783  // << dphi[i][qp](0)
784  // << ","
785  // << dphi[i][qp](1)
786  // << ","
787  // << dphi[i][qp](2)
788  // << ")},"
789  // << std::endl;
790  // }
791 
792  // Get reference vals for this Order/Family combination.
793  const TabulatedVals & val_table =
794  tabulated_vals[std::make_pair(radial_order, radial_family)];
795  const TabulatedGrads & grad_table =
796  tabulated_grads[std::make_pair(radial_order, radial_family)];
797 
798  // Test whether the computed values match the tabulated values to
799  // the specified accuracy.
800  for (const auto & t : val_table)
801  LIBMESH_ASSERT_FP_EQUAL(t.val, phi[t.i][t.qp], TOLERANCE);
802  for (const auto & t : grad_table)
803  LIBMESH_ASSERT_FP_EQUAL(0., (dphi[t.i][t.qp] - t.grad).norm_sq(), 1e-10);
804 
805  // Make sure there are actually reference values
806  libmesh_error_msg_if(val_table.empty(), "No tabulated values found!");
807  libmesh_error_msg_if(grad_table.empty(), "No tabulated gradients found!");
808 
809 #endif // LIBMESH_ENABLE_INFINITE_ELEMENTS
810  }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:182
The ReplicatedMesh class is derived from the MeshBase class, and is used to store identical copies of...
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:159
static constexpr Real TOLERANCE
unsigned int dim
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
std::vector< TabulatedGrad > TabulatedGrads
std::vector< TabulatedVal > TabulatedVals
void libmesh_ignore(const Args &...)
std::map< KeyType, TabulatedGrads > tabulated_grads
This class is used to build infinite elements on top of an existing mesh.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
std::map< KeyType, TabulatedVals > tabulated_vals
This class implements specific orders of Gauss quadrature.
unsigned int mesh_dimension() const
Definition: mesh_base.C:324
virtual bool infinite() const =0
virtual dof_id_type n_elem() const =0
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.

Member Data Documentation

◆ tabulated_grads

std::map<KeyType, TabulatedGrads> InfFERadialTest::tabulated_grads

Definition at line 219 of file inf_fe_radial_test.C.

◆ tabulated_vals

std::map<KeyType, TabulatedVals> InfFERadialTest::tabulated_vals

Definition at line 218 of file inf_fe_radial_test.C.


The documentation for this class was generated from the following file: