48 #ifndef __INTREPID2_ORIENTATIONTOOLS_DEF_MATRIX_DATA_HPP__ 49 #define __INTREPID2_ORIENTATIONTOOLS_DEF_MATRIX_DATA_HPP__ 52 #if defined (__clang__) && !defined (__INTEL_COMPILER) 53 #pragma clang system_header 59 template<
typename BasisHostType>
61 OrientationTools<DT>::createCoeffMatrixInternal(
const BasisHostType* basis) {
62 const std::string name(basis->getName());
63 CoeffMatrixDataViewType matData;
68 const auto cellTopo = basis->getBaseCellTopology();
69 const ordinal_type numEdges = cellTopo.getEdgeCount();
70 const ordinal_type numFaces = cellTopo.getFaceCount();
71 ordinal_type matDim = 0, matDim1 = 0, matDim2 = 0, numOrts = 0, numSubCells;
72 for(ordinal_type i=0; i<numEdges; ++i) {
73 matDim1 = std::max(matDim1, basis->getDofCount(1,i));
74 numOrts = std::max(numOrts,2);
76 for(ordinal_type i=0; i<numFaces; ++i) {
77 matDim2 = std::max(matDim2, basis->getDofCount(2,i));
78 numOrts = std::max(numOrts,2*ordinal_type(cellTopo.getSideCount(2,i)));
80 matDim = std::max(matDim1,matDim2);
81 numSubCells = (matDim1>0)*numEdges + (matDim2>0)*numFaces;
84 matData = CoeffMatrixDataViewType(
"Orientation::CoeffMatrix::"+name,
90 if(basis->getFunctionSpace() == FUNCTION_SPACE_HGRAD) {
91 init_HGRAD(matData, basis);
92 }
else if (basis->getFunctionSpace() == FUNCTION_SPACE_HCURL) {
93 init_HCURL(matData, basis);
94 }
else if (basis->getFunctionSpace() == FUNCTION_SPACE_HDIV) {
95 init_HDIV(matData, basis);
103 template<
typename DT>
104 template<
typename BasisHostType>
108 BasisHostType
const *cellBasis) {
110 const auto cellTopo = cellBasis->getBaseCellTopology();
111 const ordinal_type numEdges = cellTopo.getEdgeCount();
112 const ordinal_type numFaces = cellTopo.getFaceCount();
115 const ordinal_type numOrt = 2;
116 for (ordinal_type edgeId=0;edgeId<numEdges;++edgeId) {
117 if(cellBasis->getDofCount(1, edgeId) < 2)
continue;
118 for (ordinal_type edgeOrt=0;edgeOrt<numOrt;++edgeOrt) {
119 auto mat = Kokkos::subview(matData,
121 Kokkos::ALL(), Kokkos::ALL());
127 *cellBasis->getSubCellRefBasis(1,edgeId), *cellBasis,
134 for (ordinal_type faceId=0;faceId<numFaces;++faceId) {
136 const ordinal_type numOrt = 2*cellTopo.getSideCount(2,faceId);
137 if(cellBasis->getDofCount(2, faceId) < 1)
continue;
138 for (ordinal_type faceOrt=0;faceOrt<numOrt;++faceOrt) {
139 auto mat = Kokkos::subview(matData,
140 numEdges+faceId, faceOrt,
141 Kokkos::ALL(), Kokkos::ALL());
147 *cellBasis->getSubCellRefBasis(2,faceId), *cellBasis,
157 template<
typename DT>
158 template<
typename BasisHostType>
162 BasisHostType
const *cellBasis) {
163 const auto cellTopo = cellBasis->getBaseCellTopology();
164 const ordinal_type numEdges = cellTopo.getEdgeCount();
165 const ordinal_type numFaces = cellTopo.getFaceCount();
168 const ordinal_type numOrt = 2;
169 for (ordinal_type edgeId=0;edgeId<numEdges;++edgeId) {
170 if(cellBasis->getDofCount(1, edgeId) < 1)
continue;
171 for (ordinal_type edgeOrt=0;edgeOrt<numOrt;++edgeOrt) {
172 auto mat = Kokkos::subview(matData,
174 Kokkos::ALL(), Kokkos::ALL());
176 *cellBasis->getSubCellRefBasis(1,edgeId), *cellBasis,
182 for (ordinal_type faceId=0;faceId<numFaces;++faceId) {
184 const ordinal_type numOrt = 2*cellTopo.getSideCount(2,faceId);
185 if(cellBasis->getDofCount(2, faceId) < 1)
continue;
186 for (ordinal_type faceOrt=0;faceOrt<numOrt;++faceOrt) {
187 auto mat = Kokkos::subview(matData,
188 numEdges+faceId, faceOrt,
189 Kokkos::ALL(), Kokkos::ALL());
192 *cellBasis->getSubCellRefBasis(2,faceId), *cellBasis,
202 template<
typename DT>
203 template<
typename BasisHostType>
207 BasisHostType
const *cellBasis) {
208 const auto cellTopo = cellBasis->getBaseCellTopology();
209 const ordinal_type numSides = cellTopo.getSideCount();
210 const ordinal_type sideDim = cellTopo.getDimension()-1;
213 for (ordinal_type sideId=0;sideId<numSides;++sideId) {
214 if(cellBasis->getDofCount(sideDim, sideId) < 1)
continue;
215 const ordinal_type numOrt = (sideDim == 1) ? 2 : 2*cellTopo.getSideCount(sideDim,sideId);
216 for (ordinal_type faceOrt=0;faceOrt<numOrt;++faceOrt) {
217 auto mat = Kokkos::subview(matData,
219 Kokkos::ALL(), Kokkos::ALL());
221 *cellBasis->getSubCellRefBasis(sideDim,sideId), *cellBasis,
228 template<
typename DT>
229 template<
typename BasisType>
232 #ifdef HAVE_INTREPID2_DEBUG 233 INTREPID2_TEST_FOR_EXCEPTION( !basis->requireOrientation(), std::invalid_argument,
234 ">>> ERROR (OrientationTools::createCoeffMatrix): basis does not require orientations." );
236 Kokkos::push_finalize_hook( [=] {
240 const std::pair<std::string,ordinal_type> key(basis->getName(), basis->getDegree());
241 const auto found = ortCoeffData.find(key);
244 if (found == ortCoeffData.end()) {
246 auto basis_host = basis->getHostBasis();
247 matData = createCoeffMatrixInternal(basis_host.getRawPtr());
248 ortCoeffData.insert(std::make_pair(key, matData));
251 matData = found->second;
257 template<
typename DT>
259 ortCoeffData.clear();