43 #ifndef PANZER_INTREPID_BASIS_FACTORY_H 44 #define PANZER_INTREPID_BASIS_FACTORY_H 49 #include "Teuchos_RCP.hpp" 50 #include "Intrepid2_Basis.hpp" 52 #include "Shards_CellTopology.hpp" 54 #include "Intrepid2_HGRAD_QUAD_C1_FEM.hpp" 55 #include "Intrepid2_HGRAD_QUAD_C2_FEM.hpp" 57 #include "Intrepid2_HGRAD_HEX_C1_FEM.hpp" 58 #include "Intrepid2_HGRAD_HEX_C2_FEM.hpp" 60 #include "Intrepid2_HGRAD_TET_C1_FEM.hpp" 61 #include "Intrepid2_HGRAD_TET_C2_FEM.hpp" 63 #include "Intrepid2_HGRAD_TRI_C1_FEM.hpp" 64 #include "Intrepid2_HGRAD_TRI_C2_FEM.hpp" 66 #include "Intrepid2_HGRAD_LINE_C1_FEM.hpp" 68 #include "Intrepid2_HCURL_TRI_I1_FEM.hpp" 69 #include "Intrepid2_HCURL_TET_I1_FEM.hpp" 70 #include "Intrepid2_HCURL_QUAD_I1_FEM.hpp" 71 #include "Intrepid2_HCURL_HEX_I1_FEM.hpp" 73 #include "Intrepid2_HDIV_TRI_I1_FEM.hpp" 74 #include "Intrepid2_HDIV_QUAD_I1_FEM.hpp" 75 #include "Intrepid2_HDIV_TET_I1_FEM.hpp" 76 #include "Intrepid2_HDIV_HEX_I1_FEM.hpp" 96 template <
typename ScalarT,
typename ArrayT>
97 Teuchos::RCP<Intrepid2::Basis<ScalarT,ArrayT> >
99 const shards::CellTopology & cell_topology)
104 std::string cell_topology_type = cell_topology.getName();
105 std::size_t end_position = 0;
106 end_position = cell_topology_type.find(
"_");
107 std::string cell_type = cell_topology_type.substr(0,end_position);
109 Teuchos::RCP<Intrepid2::Basis<ScalarT,ArrayT> >
basis;
111 if ( (basis_type ==
"Const") && (basis_order == 0) )
114 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Hexahedron") && (basis_order == 1) )
115 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_HEX_C1_FEM<ScalarT,ArrayT> );
117 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Hexahedron") && (basis_order == 2) )
118 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_HEX_C2_FEM<ScalarT,ArrayT> );
120 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Hexahedron") && (basis_order == 1) )
121 basis = Teuchos::rcp(
new Intrepid2::Basis_HCURL_HEX_I1_FEM<ScalarT,ArrayT> );
123 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Hexahedron") && (basis_order == 1) )
124 basis = Teuchos::rcp(
new Intrepid2::Basis_HDIV_HEX_I1_FEM<ScalarT,ArrayT> );
126 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Tetrahedron") && (basis_order == 1) )
127 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TET_C1_FEM<ScalarT,ArrayT> );
129 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Tetrahedron") && (basis_order == 2) )
130 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TET_C1_FEM<ScalarT,ArrayT> );
132 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Tetrahedron") && (basis_order == 1) )
133 basis = Teuchos::rcp(
new Intrepid2::Basis_HCURL_TET_I1_FEM<ScalarT,ArrayT> );
135 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Tetrahedron") && (basis_order == 1) )
136 {
basis = Teuchos::rcp(
new Intrepid2::Basis_HDIV_TET_I1_FEM<ScalarT,ArrayT> ); }
138 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Quadrilateral") && (basis_order == 1) )
139 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_QUAD_C1_FEM<ScalarT,ArrayT> );
141 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Quadrilateral") && (basis_order == 2) )
142 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_QUAD_C2_FEM<ScalarT,ArrayT> );
144 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Quadrilateral") && (basis_order == 1) )
145 basis = Teuchos::rcp(
new Intrepid2::Basis_HCURL_QUAD_I1_FEM<ScalarT,ArrayT> );
147 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Quadrilateral") && (basis_order == 1) )
148 basis = Teuchos::rcp(
new Intrepid2::Basis_HDIV_QUAD_I1_FEM<ScalarT,ArrayT> );
150 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Triangle") && (basis_order == 1) )
151 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TRI_C1_FEM<ScalarT,ArrayT> );
153 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Triangle") && (basis_order == 2) )
154 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_TRI_C2_FEM<ScalarT,ArrayT> );
156 else if ( (basis_type ==
"HCurl") && (cell_type ==
"Triangle") && (basis_order == 1) )
157 basis = Teuchos::rcp(
new Intrepid2::Basis_HCURL_TRI_I1_FEM<ScalarT,ArrayT> );
159 else if ( (basis_type ==
"HDiv") && (cell_type ==
"Triangle") && (basis_order == 1) )
160 basis = Teuchos::rcp(
new Intrepid2::Basis_HDIV_TRI_I1_FEM<ScalarT,ArrayT> );
162 else if ( (basis_type ==
"HGrad") && (cell_type ==
"Line") && (basis_order == 1) )
163 basis = Teuchos::rcp(
new Intrepid2::Basis_HGRAD_LINE_C1_FEM<ScalarT,ArrayT> );
165 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::is_null(
basis), std::runtime_error,
166 "Failed to create the requestedbasis with basis_type=\"" << basis_type <<
167 "\", basis_order=\"" << basis_order <<
"\", and cell_type=\"" << cell_type <<
"\"!\n");
169 TEUCHOS_TEST_FOR_EXCEPTION(cell_topology!=
basis->getBaseCellTopology(),
171 "Failed to create basis. Intrepid2 basis topology does not match mesh cell topology!");
189 template <
typename ScalarT,
typename ArrayT>
190 Teuchos::RCP<Intrepid2::Basis<ScalarT,ArrayT> >
192 const Teuchos::RCP<const shards::CellTopology> & cell_topology)
194 return createIntrepid2Basis<ScalarT,ArrayT>(basis_type,basis_order,*cell_topology);
Teuchos::RCP< Intrepid2::Basis< ScalarT, ArrayT > > 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.
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.