43 #ifndef PANZER_INTREPID_BASIS_FACTORY_H
44 #define PANZER_INTREPID_BASIS_FACTORY_H
50 #include "Intrepid2_HVOL_C0_FEM.hpp"
51 #include "Intrepid2_HVOL_TRI_Cn_FEM.hpp"
52 #include "Intrepid2_HVOL_QUAD_Cn_FEM.hpp"
53 #include "Intrepid2_HVOL_HEX_Cn_FEM.hpp"
54 #include "Intrepid2_HVOL_TET_Cn_FEM.hpp"
57 #include "Intrepid2_Basis.hpp"
59 #include "Shards_CellTopology.hpp"
61 #include "Intrepid2_HGRAD_QUAD_C1_FEM.hpp"
62 #include "Intrepid2_HGRAD_QUAD_C2_FEM.hpp"
63 #include "Intrepid2_HGRAD_QUAD_Cn_FEM.hpp"
65 #include "Intrepid2_HGRAD_HEX_C1_FEM.hpp"
66 #include "Intrepid2_HGRAD_HEX_C2_FEM.hpp"
67 #include "Intrepid2_HGRAD_HEX_Cn_FEM.hpp"
69 #include "Intrepid2_HGRAD_TET_C1_FEM.hpp"
70 #include "Intrepid2_HGRAD_TET_C2_FEM.hpp"
71 #include "Intrepid2_HGRAD_TET_Cn_FEM.hpp"
73 #include "Intrepid2_HGRAD_TRI_C1_FEM.hpp"
74 #include "Intrepid2_HGRAD_TRI_C2_FEM.hpp"
75 #include "Intrepid2_HGRAD_TRI_Cn_FEM.hpp"
77 #include "Intrepid2_HGRAD_LINE_C1_FEM.hpp"
78 #include "Intrepid2_HGRAD_LINE_Cn_FEM.hpp"
80 #include "Intrepid2_HCURL_TRI_I1_FEM.hpp"
81 #include "Intrepid2_HCURL_TRI_In_FEM.hpp"
83 #include "Intrepid2_HCURL_TET_I1_FEM.hpp"
84 #include "Intrepid2_HCURL_TET_In_FEM.hpp"
86 #include "Intrepid2_HCURL_QUAD_I1_FEM.hpp"
87 #include "Intrepid2_HCURL_QUAD_In_FEM.hpp"
89 #include "Intrepid2_HCURL_HEX_I1_FEM.hpp"
90 #include "Intrepid2_HCURL_HEX_In_FEM.hpp"
92 #include "Intrepid2_HDIV_TRI_I1_FEM.hpp"
93 #include "Intrepid2_HDIV_TRI_In_FEM.hpp"
95 #include "Intrepid2_HDIV_QUAD_I1_FEM.hpp"
96 #include "Intrepid2_HDIV_QUAD_In_FEM.hpp"
98 #include "Intrepid2_HDIV_TET_I1_FEM.hpp"
99 #include "Intrepid2_HDIV_TET_In_FEM.hpp"
101 #include "Intrepid2_HDIV_HEX_I1_FEM.hpp"
102 #include "Intrepid2_HDIV_HEX_In_FEM.hpp"
121 template <
typename ExecutionSpace,
typename OutputValueType,
typename Po
intValueType>
124 const shards::CellTopology & cell_topology)
129 std::string cell_topology_type = cell_topology.getName();
130 std::size_t end_position = 0;
131 end_position = cell_topology_type.find(
"_");
132 std::string cell_type = cell_topology_type.substr(0,end_position);
139 const Intrepid2::EPointType point_type = Intrepid2::POINTTYPE_EQUISPACED;
141 if ( (basis_type ==
"Const") && (basis_order == 0) )
142 basis =
Teuchos::rcp(
new Intrepid2::Basis_HVOL_C0_FEM<ExecutionSpace,OutputValueType,PointValueType>(cell_topology) );
144 else if ( (basis_type ==
"HVol") && (basis_order == 0) )
145 basis =
Teuchos::rcp(
new Intrepid2::Basis_HVOL_C0_FEM<ExecutionSpace,OutputValueType,PointValueType>(cell_topology) );
147 else if ( (basis_type ==
"HVol") && (cell_type ==
"Quadrilateral") && (basis_order > 0) )
148 basis =
Teuchos::rcp(
new Intrepid2::Basis_HVOL_QUAD_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
150 else if ( (basis_type ==
"HVol") && (cell_type ==
"Triangle") && (basis_order > 0) )
151 basis =
Teuchos::rcp(
new Intrepid2::Basis_HVOL_TRI_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
153 else if ( (basis_type ==
"HVol") && (cell_type ==
"Hexahedron") )
154 basis =
Teuchos::rcp(
new Intrepid2::Basis_HVOL_HEX_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
156 else if ( (basis_type ==
"HVol") && (cell_type ==
"Tetrahedron") )
157 basis =
Teuchos::rcp(
new Intrepid2::Basis_HVOL_TET_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
159 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Hexahedron") && (basis_order == 1) )
160 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_HEX_C1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
162 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Hexahedron") && (basis_order == 2) )
163 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_HEX_C2_FEM<ExecutionSpace,OutputValueType,PointValueType> );
165 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Hexahedron") && (basis_order > 2) )
166 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_HEX_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
168 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Hexahedron") && (basis_order == 1) )
169 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_HEX_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
171 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Hexahedron") && (basis_order > 1) )
172 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_HEX_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
174 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Hexahedron") && (basis_order == 1) )
175 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_HEX_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
177 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Hexahedron") && (basis_order > 1) )
178 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_HEX_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
180 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Tetrahedron") && (basis_order == 1) )
181 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TET_C1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
183 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Tetrahedron") && (basis_order == 2) )
184 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TET_C2_FEM<ExecutionSpace,OutputValueType,PointValueType> );
186 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Tetrahedron") && (basis_order > 2) )
187 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TET_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
189 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Tetrahedron") && (basis_order == 1) )
190 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_TET_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
192 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Tetrahedron") && (basis_order > 1) )
193 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_TET_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
195 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Tetrahedron") && (basis_order == 1) )
196 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_TET_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
198 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Tetrahedron") && (basis_order > 1) )
199 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_TET_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
201 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Quadrilateral") && (basis_order == 1) )
202 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_QUAD_C1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
204 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Quadrilateral") && (basis_order == 2) )
205 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_QUAD_C2_FEM<ExecutionSpace,OutputValueType,PointValueType> );
207 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Quadrilateral") && (basis_order > 2) )
208 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_QUAD_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
210 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Quadrilateral") && (basis_order == 1) )
211 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_QUAD_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
213 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Quadrilateral") && (basis_order > 1) )
214 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_QUAD_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
216 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Quadrilateral") && (basis_order == 1) )
217 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_QUAD_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
219 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Quadrilateral") && (basis_order > 1) )
220 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_QUAD_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
222 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Triangle") && (basis_order == 1) )
223 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TRI_C1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
225 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Triangle") && (basis_order == 2) )
226 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TRI_C2_FEM<ExecutionSpace,OutputValueType,PointValueType> );
228 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Triangle") && (basis_order > 2) )
229 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TRI_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
231 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Triangle") && (basis_order == 1) )
232 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_TRI_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
234 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Triangle") && (basis_order > 1) )
235 basis =
Teuchos::rcp(
new Intrepid2::Basis_HCURL_TRI_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
237 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Triangle") && (basis_order == 1) )
238 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_TRI_I1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
240 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Triangle") && (basis_order > 1) )
241 basis =
Teuchos::rcp(
new Intrepid2::Basis_HDIV_TRI_In_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
243 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Line") && (basis_order == 1) )
244 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_LINE_C1_FEM<ExecutionSpace,OutputValueType,PointValueType> );
246 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Line") && (basis_order >= 2) )
247 basis =
Teuchos::rcp(
new Intrepid2::Basis_HGRAD_LINE_Cn_FEM<ExecutionSpace,OutputValueType,PointValueType>(basis_order, point_type) );
250 "Failed to create the requestedbasis with basis_type=\"" << basis_type <<
251 "\", basis_order=\"" << basis_order <<
"\", and cell_type=\"" << cell_type <<
"\"!\n");
259 "Failed to create basis. Intrepid2 basis base topology does not match mesh cell base topology!");
278 template <
typename ExecutionSpace,
typename OutputValueType,
typename Po
intValueType>
283 return createIntrepid2Basis<ExecutionSpace,OutputValueType,PointValueType>(basis_type,basis_order,*cell_topology);
bool is_null(const std::shared_ptr< T > &p)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< Intrepid2::Basis< ExecutionSpace, OutputValueType, PointValueType > > createIntrepid2Basis(const std::string basis_type, int basis_order, const shards::CellTopology &cell_topology)
Creates an Intrepid2::Basis object given the basis, order and cell topology.