Intrepid2
Intrepid2_ProjectionStructDef.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid2 Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38 // Mauro Perego (mperego@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 
47 #ifndef __INTREPID2_PROJECTIONSTRUCTDEF_HPP__
48 #define __INTREPID2_PROJECTIONSTRUCTDEF_HPP__
49 
50 #include "Intrepid2_Utils.hpp"
51 #include "Shards_CellTopology.hpp"
52 #include "Shards_BasicTopologies.hpp"
53 
58 
60 
61 
62 namespace Intrepid2 {
63 
64 namespace Experimental {
65 template<typename DeviceType, typename ValueType>
66 template<typename BasisPtrType>
68  const ordinal_type targetCubDegree) {
69  const auto cellTopo = cellBasis->getBaseCellTopology();
70  std::string name = cellBasis->getName();
71  ordinal_type dim = cellTopo.getDimension();
72  numBasisEvalPoints = 0;
73  numBasisDerivEvalPoints = 0;
74  numTargetEvalPoints = 0;
75  numTargetDerivEvalPoints = 0;
76  const ordinal_type edgeDim = 1;
77  const ordinal_type faceDim = 2;
78 
79  ordinal_type basisCubDegree = cellBasis->getDegree();
80  ordinal_type edgeBasisCubDegree, faceBasisCubDegree;
81 
82  ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
83  ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
84  ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
85 
86  ordinal_type hasCellDofs = (cellBasis->getDofCount(dim, 0) > 0);
87 
88  INTREPID2_TEST_FOR_ABORT( (numVertices > maxSubCellsCount) || (numFaces > maxSubCellsCount) || (numEdges > maxSubCellsCount),
89  ">>> ERROR (Intrepid2::ProjectionStruct:createHDivProjectionStruct, Projections do not support a cell topology with this cub cells count");
90 
91 
92  numBasisEvalPoints += numVertices;
93  numTargetEvalPoints += numVertices;
94  view_type coord("vertex_coord", dim);
95 
96  basisPointsRange = range_tag("basisPointsRange", 4,maxSubCellsCount);
97  targetPointsRange = range_tag("targetPointsRange", 4,maxSubCellsCount);
98  basisDerivPointsRange = range_tag("basisDerivPointsRange", 4,maxSubCellsCount);
99  targetDerivPointsRange = range_tag("targetDerivPointsRange", numberSubCellDims,maxSubCellsCount);
100  subCellTopologyKey = key_tag("subCellTopologyKey",numberSubCellDims,maxSubCellsCount);
101 
102  maxNumBasisEvalPoints = numVertices; maxNumTargetEvalPoints = numVertices;
103  for(ordinal_type iv=0; iv<numVertices; ++iv) {
104  basisPointsRange(0,iv) = range_type(iv, iv+1);
105  basisCubPoints[0][iv] = view_type("basisCubPoints",1,dim);
106  targetPointsRange(0,iv) = range_type(iv, iv+1);
107  targetCubPoints[0][iv] = view_type("targetCubPoints",1,dim);
109  for(ordinal_type d=0; d<dim; d++) {
110  basisCubPoints[0][iv](0,d) = coord(d);
111  targetCubPoints[0][iv](0,d) = coord(d);
112  }
113  }
114 
115  if (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HCURL) {
116  edgeBasisCubDegree = basisCubDegree - 1;
117  faceBasisCubDegree = basisCubDegree;
118  }
119  else if (cellBasis->getFunctionSpace() == FUNCTION_SPACE_HDIV) {
120  edgeBasisCubDegree = basisCubDegree - 1;
121  faceBasisCubDegree = basisCubDegree -1;
122  }
123  else {
124  edgeBasisCubDegree = basisCubDegree;
125  faceBasisCubDegree = basisCubDegree;
126  }
127 
128  DefaultCubatureFactory cub_factory;
129  for(ordinal_type ie=0; ie<numEdges; ++ie) {
130  ordinal_type cub_degree = 2*edgeBasisCubDegree;
131  subCellTopologyKey(edgeDim,ie) = cellBasis->getBaseCellTopology().getKey(edgeDim, ie);
132  auto edgeBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(cellBasis->getBaseCellTopology().getKey(edgeDim, ie), cub_degree);
133  basisPointsRange(edgeDim,ie) = range_type(numBasisEvalPoints, numBasisEvalPoints+edgeBasisCub->getNumPoints());
134  numBasisEvalPoints += edgeBasisCub->getNumPoints();
135  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, edgeBasisCub->getNumPoints());
136  basisCubPoints[edgeDim][ie] = view_type("basisCubPoints",edgeBasisCub->getNumPoints(),edgeDim);
137  basisCubWeights[edgeDim][ie] = view_type("basisCubWeights",edgeBasisCub->getNumPoints());
138  edgeBasisCub->getCubature(basisCubPoints[edgeDim][ie], basisCubWeights[edgeDim][ie]);
139 
140  cub_degree = edgeBasisCubDegree + targetCubDegree;
141  auto edgeTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(cellBasis->getBaseCellTopology().getKey(edgeDim, ie), cub_degree);
142  targetPointsRange(edgeDim,ie) = range_type(numTargetEvalPoints, numTargetEvalPoints+edgeTargetCub->getNumPoints());
143  numTargetEvalPoints += edgeTargetCub->getNumPoints();
144  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, edgeTargetCub->getNumPoints());
145  targetCubPoints[edgeDim][ie] = view_type("targetCubPoints",edgeTargetCub->getNumPoints(),edgeDim);
146  targetCubWeights[edgeDim][ie] = view_type("targetCubWeights",edgeTargetCub->getNumPoints());
147  edgeTargetCub->getCubature(targetCubPoints[edgeDim][ie], targetCubWeights[edgeDim][ie]);
148  }
149 
150  for(ordinal_type iface=0; iface<numFaces; ++iface) {
151  ordinal_type cub_degree = 2*faceBasisCubDegree;
152  subCellTopologyKey(faceDim,iface) = cellBasis->getBaseCellTopology().getKey(faceDim, iface);
153  auto faceBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
154  basisPointsRange(faceDim,iface) = range_type(numBasisEvalPoints, numBasisEvalPoints+faceBasisCub->getNumPoints());
155  numBasisEvalPoints += faceBasisCub->getNumPoints();
156  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, faceBasisCub->getNumPoints());
157  basisCubPoints[faceDim][iface] = view_type("basisCubPoints",faceBasisCub->getNumPoints(),faceDim);
158  basisCubWeights[faceDim][iface] = view_type("basisCubWeights",faceBasisCub->getNumPoints());
159  faceBasisCub->getCubature(basisCubPoints[faceDim][iface], basisCubWeights[faceDim][iface]);
160 
161  cub_degree = faceBasisCubDegree + targetCubDegree;
162  auto faceTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
163  targetPointsRange(faceDim,iface) = range_type(numTargetEvalPoints, numTargetEvalPoints+faceTargetCub->getNumPoints());
164  numTargetEvalPoints += faceTargetCub->getNumPoints();
165  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, faceTargetCub->getNumPoints());
166  targetCubPoints[faceDim][iface] = view_type("targetCubPoints",faceTargetCub->getNumPoints(),faceDim);
167  targetCubWeights[faceDim][iface] = view_type("targetCubWeights",faceTargetCub->getNumPoints());
168  faceTargetCub->getCubature(targetCubPoints[faceDim][iface], targetCubWeights[faceDim][iface]);
169  }
170  subCellTopologyKey(dim,0) = cellBasis->getBaseCellTopology().getBaseKey();
171  if(hasCellDofs) {
172  ordinal_type cub_degree = 2*basisCubDegree;
173  auto elemBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
174  basisPointsRange(dim,0) = range_type(numBasisEvalPoints, numBasisEvalPoints+elemBasisCub->getNumPoints());
175  numBasisEvalPoints += elemBasisCub->getNumPoints();
176  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, elemBasisCub->getNumPoints());
177  basisCubPoints[dim][0] = view_type("basisCubPoints",elemBasisCub->getNumPoints(),dim);
178  basisCubWeights[dim][0] = view_type("basisCubWeights",elemBasisCub->getNumPoints());
179  elemBasisCub->getCubature(basisCubPoints[dim][0], basisCubWeights[dim][0]);
180 
181  cub_degree = basisCubDegree + targetCubDegree;
182  auto elemTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
183  targetPointsRange(dim,0) = range_type(numTargetEvalPoints, numTargetEvalPoints+elemTargetCub->getNumPoints());
184  numTargetEvalPoints += elemTargetCub->getNumPoints();
185  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, elemTargetCub->getNumPoints());
186  targetCubPoints[dim][0] = view_type("targetCubPoints",elemTargetCub->getNumPoints(),dim);
187  targetCubWeights[dim][0] = view_type("targetCubWeights",elemTargetCub->getNumPoints());
188  elemTargetCub->getCubature(targetCubPoints[dim][0], targetCubWeights[dim][0]);
189  }
190 }
191 
192 template<typename DeviceType, typename ValueType>
193 template<typename BasisPtrType>
195  const ordinal_type targetCubDegree,
196  const ordinal_type targetGradCubDegree) {
197  const auto cellTopo = cellBasis->getBaseCellTopology();
198  std::string name = cellBasis->getName();
199  ordinal_type dim = cellTopo.getDimension();
200  numBasisEvalPoints = 0;
201  numBasisDerivEvalPoints = 0;
202  numTargetEvalPoints = 0;
203  numTargetDerivEvalPoints = 0;
204  const ordinal_type edgeDim = 1;
205  const ordinal_type faceDim = 2;
206 
207  ordinal_type basisCubDegree = cellBasis->getDegree();
208  ordinal_type edgeBasisCubDegree = basisCubDegree;
209  ordinal_type faceBasisCubDegree = basisCubDegree;
210 
211  ordinal_type numVertices = (cellBasis->getDofCount(0, 0) > 0) ? cellTopo.getVertexCount() : 0;
212  ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount(): 0;
213  ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
214 
215  INTREPID2_TEST_FOR_ABORT( (numFaces > maxSubCellsCount) || (numEdges > maxSubCellsCount),
216  ">>> ERROR (Intrepid2::ProjectionStruct:createHDivProjectionStruct, Projections do not support a cell topology with this cub cells count");
217 
218 
219  ordinal_type hasCellDofs = (cellBasis->getDofCount(dim, 0) > 0);
220 
221  maxNumBasisEvalPoints = numVertices; maxNumTargetEvalPoints = numVertices;
222  maxNumBasisDerivEvalPoints = 0; maxNumTargetDerivEvalPoints = 0;
223 
224  basisPointsRange = range_tag("basisPointsRange", numberSubCellDims,maxSubCellsCount);
225  targetPointsRange = range_tag("targetPointsRange", numberSubCellDims,maxSubCellsCount);
226  basisDerivPointsRange = range_tag("basisDerivPointsRange", numberSubCellDims,maxSubCellsCount);
227  targetDerivPointsRange = range_tag("targetDerivPointsRange", numberSubCellDims,maxSubCellsCount);
228  subCellTopologyKey = key_tag("subCellTopologyKey",numberSubCellDims,maxSubCellsCount);
229 
230  numBasisEvalPoints += numVertices;
231  numTargetEvalPoints += numVertices;
232  view_type coord("vertex_coord", dim);
233  for(ordinal_type iv=0; iv<numVertices; ++iv) {
234  basisPointsRange(0,iv) = range_type(iv, iv+1);
235  basisCubPoints[0][iv] = view_type("basisCubPoints",1,dim);
236  targetPointsRange(0,iv) = range_type(iv, iv+1);
237  targetCubPoints[0][iv] = view_type("targetCubPoints",1,dim);
239  for(ordinal_type d=0; d<dim; d++) {
240  basisCubPoints[0][iv](0,d) = coord(d);
241  targetCubPoints[0][iv](0,d) = coord(d);
242  }
243  }
244 
245  DefaultCubatureFactory cub_factory;
246  for(ordinal_type ie=0; ie<numEdges; ++ie) {
247  ordinal_type cub_degree = 2*edgeBasisCubDegree;
248  subCellTopologyKey(edgeDim,ie) = cellBasis->getBaseCellTopology().getKey(edgeDim, ie);
249  auto edgeBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(cellBasis->getBaseCellTopology().getKey(edgeDim, ie), cub_degree);
250  basisDerivPointsRange(edgeDim,ie) = range_type(numBasisDerivEvalPoints, numBasisDerivEvalPoints+edgeBasisCub->getNumPoints());
251  numBasisDerivEvalPoints += edgeBasisCub->getNumPoints();
252  maxNumBasisDerivEvalPoints = std::max(maxNumBasisDerivEvalPoints, edgeBasisCub->getNumPoints());
253  basisDerivCubPoints[edgeDim][ie] = view_type("basisDerivCubPoints",edgeBasisCub->getNumPoints(),edgeDim);
254  basisDerivCubWeights[edgeDim][ie] = view_type("basisDerivCubWeights",edgeBasisCub->getNumPoints());
255  edgeBasisCub->getCubature(basisDerivCubPoints[edgeDim][ie], basisDerivCubWeights[edgeDim][ie]);
256 
257  cub_degree = edgeBasisCubDegree + targetGradCubDegree;
258  auto edgeTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(cellBasis->getBaseCellTopology().getKey(edgeDim, ie), cub_degree);
259  targetDerivPointsRange(edgeDim,ie) = range_type(numTargetDerivEvalPoints, numTargetDerivEvalPoints+edgeTargetCub->getNumPoints());
260  numTargetDerivEvalPoints += edgeTargetCub->getNumPoints();
261  maxNumTargetDerivEvalPoints = std::max(maxNumTargetDerivEvalPoints, edgeTargetCub->getNumPoints());
262  targetDerivCubPoints[edgeDim][ie] = view_type("targetDerivCubPoints",edgeTargetCub->getNumPoints(),edgeDim);
263  targetDerivCubWeights[edgeDim][ie] = view_type("targetDerivCubWeights",edgeTargetCub->getNumPoints());
264  edgeTargetCub->getCubature(targetDerivCubPoints[edgeDim][ie], targetDerivCubWeights[edgeDim][ie]);
265  }
266 
267  for(ordinal_type iface=0; iface<numFaces; ++iface) {
268  ordinal_type cub_degree = 2*faceBasisCubDegree;
269  subCellTopologyKey(faceDim,iface) = cellBasis->getBaseCellTopology().getKey(faceDim, iface);
270  auto faceBasisGradCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
271  basisDerivPointsRange(faceDim,iface) = range_type(numBasisDerivEvalPoints, numBasisDerivEvalPoints+faceBasisGradCub->getNumPoints());
272  numBasisDerivEvalPoints += faceBasisGradCub->getNumPoints();
273  maxNumBasisDerivEvalPoints = std::max(maxNumBasisDerivEvalPoints, faceBasisGradCub->getNumPoints());
274  basisDerivCubPoints[faceDim][iface] = view_type("basisDerivCubPoints",faceBasisGradCub->getNumPoints(),faceDim);
275  basisDerivCubWeights[faceDim][iface] = view_type("basisDerivCubWeights",faceBasisGradCub->getNumPoints());
276  faceBasisGradCub->getCubature(basisDerivCubPoints[faceDim][iface], basisDerivCubWeights[faceDim][iface]);
277 
278  cub_degree = faceBasisCubDegree + targetGradCubDegree;
279  auto faceTargetDerivCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
280  targetDerivPointsRange(faceDim,iface) = range_type(numTargetDerivEvalPoints, numTargetDerivEvalPoints+faceTargetDerivCub->getNumPoints());
281  numTargetDerivEvalPoints += faceTargetDerivCub->getNumPoints();
282  maxNumTargetDerivEvalPoints = std::max(maxNumTargetDerivEvalPoints, faceTargetDerivCub->getNumPoints());
283  targetDerivCubPoints[faceDim][iface] = view_type("targetDerivCubPoints",faceTargetDerivCub->getNumPoints(),faceDim);
284  targetDerivCubWeights[faceDim][iface] = view_type("targetDerivCubWeights",faceTargetDerivCub->getNumPoints());
285  faceTargetDerivCub->getCubature(targetDerivCubPoints[faceDim][iface], targetDerivCubWeights[faceDim][iface]);
286  }
287  subCellTopologyKey(dim,0) = cellBasis->getBaseCellTopology().getBaseKey();
288  if(hasCellDofs) {
289  ordinal_type cub_degree = 2*basisCubDegree;
290  auto elemBasisGradCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
291  basisDerivPointsRange(dim,0) = range_type(numBasisDerivEvalPoints, numBasisDerivEvalPoints+elemBasisGradCub->getNumPoints());
292  numBasisDerivEvalPoints += elemBasisGradCub->getNumPoints();
293  maxNumBasisDerivEvalPoints = std::max(maxNumBasisDerivEvalPoints, elemBasisGradCub->getNumPoints());
294  basisDerivCubPoints[dim][0] = view_type("basisDerivCubPoints",elemBasisGradCub->getNumPoints(),dim);
295  basisDerivCubWeights[dim][0] = view_type("basisDerivCubWeights",elemBasisGradCub->getNumPoints());
296  elemBasisGradCub->getCubature(basisDerivCubPoints[dim][0], basisDerivCubWeights[dim][0]);
297 
298  cub_degree = basisCubDegree + targetGradCubDegree;
299  auto elemTargetGradCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
300  targetDerivPointsRange(dim,0) = range_type(numTargetDerivEvalPoints, numTargetDerivEvalPoints+elemTargetGradCub->getNumPoints());
301  numTargetDerivEvalPoints += elemTargetGradCub->getNumPoints();
302  maxNumTargetDerivEvalPoints = std::max(maxNumTargetDerivEvalPoints, elemTargetGradCub->getNumPoints());
303  targetDerivCubPoints[dim][0] = view_type("targetDerivCubPoints",elemTargetGradCub->getNumPoints(),dim);
304  targetDerivCubWeights[dim][0] = view_type("targetDerivCubWeights",elemTargetGradCub->getNumPoints());
305  elemTargetGradCub->getCubature(targetDerivCubPoints[dim][0], targetDerivCubWeights[dim][0]);
306  }
307 }
308 
309 template<typename DeviceType, typename ValueType>
310 template<typename BasisPtrType>
312  const ordinal_type targetCubDegree,
313  const ordinal_type targetCurlCubDegre) {
314  const auto cellTopo = cellBasis->getBaseCellTopology();
315  std::string name = cellBasis->getName();
316  ordinal_type dim = cellTopo.getDimension();
317  numBasisEvalPoints = 0;
318  numBasisDerivEvalPoints = 0;
319  numTargetEvalPoints = 0;
320  numTargetDerivEvalPoints = 0;
321  const ordinal_type edgeDim = 1;
322  const ordinal_type faceDim = 2;
323 
324  ordinal_type basisCubDegree = cellBasis->getDegree();
325  ordinal_type edgeBasisCubDegree = basisCubDegree - 1;
326  ordinal_type faceBasisCubDegree = basisCubDegree;
327  ordinal_type numFaces = (cellBasis->getDofCount(2, 0) > 0) ? cellTopo.getFaceCount() : 0;
328  ordinal_type numEdges = (cellBasis->getDofCount(1, 0) > 0) ? cellTopo.getEdgeCount() : 0;
329  ordinal_type hasCellDofs = (cellBasis->getDofCount(dim, 0) > 0);
330 
331  maxNumBasisEvalPoints = 0; maxNumTargetEvalPoints = 0;
332  maxNumBasisDerivEvalPoints = 0; maxNumTargetDerivEvalPoints = 0;
333 
334  basisPointsRange = range_tag("basisPointsRange", numberSubCellDims,maxSubCellsCount);
335  targetPointsRange = range_tag("targetPointsRange", numberSubCellDims,maxSubCellsCount);
336  basisDerivPointsRange = range_tag("basisDerivPointsRange", numberSubCellDims,maxSubCellsCount);
337  targetDerivPointsRange = range_tag("targetDerivPointsRange", numberSubCellDims,maxSubCellsCount);
338  subCellTopologyKey = key_tag("subCellTopologyKey",numberSubCellDims,maxSubCellsCount);
339 
340  DefaultCubatureFactory cub_factory;
341  for(ordinal_type ie=0; ie<numEdges; ++ie) {
342  ordinal_type cub_degree = 2*edgeBasisCubDegree;
343  subCellTopologyKey(edgeDim,ie) = cellBasis->getBaseCellTopology().getKey(edgeDim, ie);
344  auto edgeBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(edgeDim,ie), cub_degree);
345  basisPointsRange(edgeDim,ie) = range_type(numBasisEvalPoints, numBasisEvalPoints+edgeBasisCub->getNumPoints());
346  numBasisEvalPoints += edgeBasisCub->getNumPoints();
347  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, edgeBasisCub->getNumPoints());
348  basisCubPoints[edgeDim][ie] = view_type("basisCubPoints",edgeBasisCub->getNumPoints(),edgeDim);
349  basisCubWeights[edgeDim][ie] = view_type("basisCubWeights",edgeBasisCub->getNumPoints());
350  edgeBasisCub->getCubature(basisCubPoints[edgeDim][ie], basisCubWeights[edgeDim][ie]);
351 
352  cub_degree = edgeBasisCubDegree + targetCubDegree;
353  auto edgeTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(edgeDim,ie), cub_degree);
354  targetPointsRange(edgeDim,ie) = range_type(numTargetEvalPoints, numTargetEvalPoints+edgeTargetCub->getNumPoints());
355  numTargetEvalPoints += edgeTargetCub->getNumPoints();
356  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, edgeTargetCub->getNumPoints());
357  targetCubPoints[edgeDim][ie] = view_type("targetCubPoints",edgeTargetCub->getNumPoints(),edgeDim);
358  targetCubWeights[edgeDim][ie] = view_type("targetCubWeights",edgeTargetCub->getNumPoints());
359  edgeTargetCub->getCubature(targetCubPoints[edgeDim][ie], targetCubWeights[edgeDim][ie]);
360  }
361 
362  for(ordinal_type iface=0; iface<numFaces; ++iface) {
363  ordinal_type cub_degree = 2*faceBasisCubDegree;
364  subCellTopologyKey(faceDim,iface) = cellBasis->getBaseCellTopology().getKey(faceDim, iface);
365  auto faceBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
366  basisPointsRange(faceDim,iface) = range_type(numBasisEvalPoints, numBasisEvalPoints+faceBasisCub->getNumPoints());
367  numBasisEvalPoints += faceBasisCub->getNumPoints();
368  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, faceBasisCub->getNumPoints());
369  basisCubPoints[faceDim][iface] = view_type("basisCubPoints",faceBasisCub->getNumPoints(),faceDim);
370  basisCubWeights[faceDim][iface] = view_type("basisCubWeights",faceBasisCub->getNumPoints());
371  faceBasisCub->getCubature(basisCubPoints[faceDim][iface], basisCubWeights[faceDim][iface]);
372 
373  auto faceBasisDerivCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
374  basisDerivPointsRange(faceDim,iface) = range_type(numBasisDerivEvalPoints, numBasisDerivEvalPoints+faceBasisCub->getNumPoints());
375  numBasisDerivEvalPoints += faceBasisCub->getNumPoints();
376  maxNumBasisDerivEvalPoints = std::max(maxNumBasisDerivEvalPoints, faceBasisCub->getNumPoints());
377  basisDerivCubPoints[faceDim][iface] = view_type("basisDerivCubPoints",faceBasisCub->getNumPoints(),faceDim);
378  basisDerivCubWeights[faceDim][iface] = view_type("basisDerivCubWeights",faceBasisCub->getNumPoints());
379  faceBasisCub->getCubature(basisDerivCubPoints[faceDim][iface], basisDerivCubWeights[faceDim][iface]);
380 
381  cub_degree = faceBasisCubDegree + targetCubDegree;
382  auto faceTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
383  targetPointsRange(faceDim,iface) = range_type(numTargetEvalPoints, numTargetEvalPoints+faceTargetCub->getNumPoints());
384  numTargetEvalPoints += faceTargetCub->getNumPoints();
385  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, faceTargetCub->getNumPoints());
386  targetCubPoints[faceDim][iface] = view_type("targetCubPoints",faceTargetCub->getNumPoints(),faceDim);
387  targetCubWeights[faceDim][iface] = view_type("targetCubWeights",faceTargetCub->getNumPoints());
388  faceTargetCub->getCubature(targetCubPoints[faceDim][iface], targetCubWeights[faceDim][iface]);
389 
390  cub_degree = faceBasisCubDegree + targetCurlCubDegre;
391  auto faceTargetDerivCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(faceDim,iface), cub_degree);
392  targetDerivPointsRange(faceDim,iface) = range_type(numTargetDerivEvalPoints, numTargetDerivEvalPoints+faceTargetDerivCub->getNumPoints());
393  numTargetDerivEvalPoints += faceTargetDerivCub->getNumPoints();
394  maxNumTargetDerivEvalPoints = std::max(maxNumTargetDerivEvalPoints, faceTargetDerivCub->getNumPoints());
395  targetDerivCubPoints[faceDim][iface] = view_type("targetDerivCubPoints",faceTargetDerivCub->getNumPoints(),faceDim);
396  targetDerivCubWeights[faceDim][iface] = view_type("targetDerivCubWeights",faceTargetDerivCub->getNumPoints());
397  faceTargetDerivCub->getCubature(targetDerivCubPoints[faceDim][iface], targetDerivCubWeights[faceDim][iface]);
398  }
399 
400  subCellTopologyKey(dim,0) = cellBasis->getBaseCellTopology().getBaseKey();
401  if(hasCellDofs) {
402  ordinal_type cub_degree = 2*basisCubDegree;
403  auto elemBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
404  basisPointsRange(dim,0) = range_type(numBasisEvalPoints, numBasisEvalPoints+elemBasisCub->getNumPoints());
405  numBasisEvalPoints += elemBasisCub->getNumPoints();
406  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, elemBasisCub->getNumPoints());
407  basisCubPoints[dim][0] = view_type("basisCubPoints",elemBasisCub->getNumPoints(),dim);
408  basisCubWeights[dim][0] = view_type("basisCubWeights",elemBasisCub->getNumPoints());
409  elemBasisCub->getCubature(basisCubPoints[dim][0], basisCubWeights[dim][0]);
410 
411  basisDerivPointsRange(dim,0) = range_type(numBasisDerivEvalPoints, numBasisDerivEvalPoints+elemBasisCub->getNumPoints());
412  numBasisDerivEvalPoints += elemBasisCub->getNumPoints();
413  maxNumBasisDerivEvalPoints = std::max(maxNumBasisDerivEvalPoints, elemBasisCub->getNumPoints());
414  basisDerivCubPoints[dim][0] = view_type("basisDerivCubPoints",elemBasisCub->getNumPoints(),dim);
415  basisDerivCubWeights[dim][0] = view_type("basisDerivCubWeights",elemBasisCub->getNumPoints());
416  elemBasisCub->getCubature(basisDerivCubPoints[dim][0], basisDerivCubWeights[dim][0]);
417 
418  cub_degree = basisCubDegree + targetCubDegree;
419  auto elemTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
420  targetPointsRange(dim,0) = range_type(numTargetEvalPoints, numTargetEvalPoints+elemTargetCub->getNumPoints());
421  numTargetEvalPoints += elemTargetCub->getNumPoints();
422  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, elemTargetCub->getNumPoints());
423  targetCubPoints[dim][0] = view_type("targetCubPoints",elemTargetCub->getNumPoints(),dim);
424  targetCubWeights[dim][0] = view_type("targetCubWeights",elemTargetCub->getNumPoints());
425  elemTargetCub->getCubature(targetCubPoints[dim][0], targetCubWeights[dim][0]);
426 
427  cub_degree = basisCubDegree + targetCurlCubDegre;
428  auto elemTargetCurlCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
429  targetDerivPointsRange(dim,0) = range_type(numTargetDerivEvalPoints, numTargetDerivEvalPoints+elemTargetCurlCub->getNumPoints());
430  numTargetDerivEvalPoints += elemTargetCurlCub->getNumPoints();
431  maxNumTargetDerivEvalPoints = std::max(maxNumTargetDerivEvalPoints, elemTargetCurlCub->getNumPoints());
432  targetDerivCubPoints[dim][0] = view_type("targetDerivCubPoints",elemTargetCurlCub->getNumPoints(),dim);
433  targetDerivCubWeights[dim][0] = view_type("targetDerivCubWeights",elemTargetCurlCub->getNumPoints());
434  elemTargetCurlCub->getCubature(targetDerivCubPoints[dim][0], targetDerivCubWeights[dim][0]);
435  }
436 }
437 
438 template<typename DeviceType, typename ValueType>
439 template<typename BasisPtrType>
441  const ordinal_type targetCubDegree,
442  const ordinal_type targetDivCubDegre) {
443  const auto cellTopo = cellBasis->getBaseCellTopology();
444  std::string name = cellBasis->getName();
445  ordinal_type dim = cellTopo.getDimension();
446  numBasisEvalPoints = 0;
447  numBasisDerivEvalPoints = 0;
448  numTargetEvalPoints = 0;
449  numTargetDerivEvalPoints = 0;
450  const ordinal_type sideDim = dim - 1;
451 
452  ordinal_type basisCubDegree = cellBasis->getDegree();
453  ordinal_type sideBasisCubDegree = basisCubDegree - 1;
454  ordinal_type numSides = cellTopo.getSideCount()*ordinal_type(cellBasis->getDofCount(sideDim, 0) > 0);
455  ordinal_type hasCellDofs = (cellBasis->getDofCount(dim, 0) > 0);
456 
457  INTREPID2_TEST_FOR_ABORT( numSides > maxSubCellsCount,
458  ">>> ERROR (Intrepid2::ProjectionStruct:createHDivProjectionStruct, Projections do not support a cell topology with so many sides");
459 
460  maxNumBasisEvalPoints = 0; maxNumTargetEvalPoints = 0;
461  maxNumBasisDerivEvalPoints = 0; maxNumTargetDerivEvalPoints = 0;
462 
463  basisPointsRange = range_tag("basisPointsRange", numberSubCellDims,maxSubCellsCount);
464  targetPointsRange = range_tag("targetPointsRange", numberSubCellDims,maxSubCellsCount);
465  basisDerivPointsRange = range_tag("basisDerivPointsRange", numberSubCellDims,maxSubCellsCount);
466  targetDerivPointsRange = range_tag("targetDerivPointsRange", numberSubCellDims,maxSubCellsCount);
467  subCellTopologyKey = key_tag("subCellTopologyKey",numberSubCellDims,maxSubCellsCount);
468 
469  Basis<HostSpaceType,ValueType,ValueType> *hcurlBasis = NULL;
470  if(cellTopo.getKey() == shards::getCellTopologyData<shards::Hexahedron<8> >()->key)
471  hcurlBasis = new Basis_HCURL_HEX_In_FEM<HostSpaceType,ValueType,ValueType>(cellBasis->getDegree());
472  else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Tetrahedron<4> >()->key)
473  hcurlBasis = new Basis_HCURL_TET_In_FEM<HostSpaceType,ValueType,ValueType>(cellBasis->getDegree());
474  else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Quadrilateral<4> >()->key)
475  hcurlBasis = new Basis_HGRAD_QUAD_Cn_FEM<HostSpaceType,ValueType,ValueType>(cellBasis->getDegree());
476  else if(cellTopo.getKey() == shards::getCellTopologyData<shards::Triangle<3> >()->key)
477  hcurlBasis = new Basis_HGRAD_TRI_Cn_FEM<HostSpaceType,ValueType,ValueType>(cellBasis->getDegree());
478  else {
479  std::stringstream ss;
480  ss << ">>> ERROR (Intrepid2::ProjectionTools::createHDivProjectionStruct): "
481  << "Method not implemented for basis " << name;
482  INTREPID2_TEST_FOR_EXCEPTION( true, std::runtime_error, ss.str().c_str() );
483  }
484 
485  bool haveHCurlConstraint = (hcurlBasis != NULL) && (hcurlBasis->getDofCount(dim,0) > 0);
486  if(hcurlBasis != NULL) delete hcurlBasis;
487 
488  DefaultCubatureFactory cub_factory;
489 
490  for(ordinal_type is=0; is<numSides; ++is) {
491  ordinal_type cub_degree = 2*sideBasisCubDegree;
492  subCellTopologyKey(sideDim,is) = cellBasis->getBaseCellTopology().getKey(sideDim, is);
493  auto sideBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(sideDim,is), cub_degree);
494  basisPointsRange(sideDim,is) = range_type(numBasisEvalPoints, numBasisEvalPoints+sideBasisCub->getNumPoints());
495  numBasisEvalPoints += sideBasisCub->getNumPoints();
496  basisCubPoints[sideDim][is] = view_type("basisCubPoints",sideBasisCub->getNumPoints(),sideDim);
497  basisCubWeights[sideDim][is] = view_type("basisCubWeights",sideBasisCub->getNumPoints());
498  sideBasisCub->getCubature(basisCubPoints[sideDim][is], basisCubWeights[sideDim][is]);
499  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, sideBasisCub->getNumPoints());
500 
501  cub_degree = sideBasisCubDegree + targetCubDegree;
502  auto sideTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(sideDim,is), cub_degree);
503  targetPointsRange(sideDim,is) = range_type(numTargetEvalPoints, numTargetEvalPoints+sideTargetCub->getNumPoints());
504  numTargetEvalPoints += sideTargetCub->getNumPoints();
505  targetCubPoints[sideDim][is] = view_type("targetCubPoints",sideTargetCub->getNumPoints(),sideDim);
506  targetCubWeights[sideDim][is] = view_type("targetCubWeights",sideTargetCub->getNumPoints());
507  sideTargetCub->getCubature(targetCubPoints[sideDim][is], targetCubWeights[sideDim][is]);
508  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, sideTargetCub->getNumPoints());
509  }
510 
511  subCellTopologyKey(dim,0) = cellBasis->getBaseCellTopology().getBaseKey();
512  if(hasCellDofs) {
513  ordinal_type cub_degree = 2*basisCubDegree - 1;
514  auto elemBasisDivCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
515  basisDerivPointsRange(dim,0) = range_type(numBasisDerivEvalPoints, numBasisDerivEvalPoints+elemBasisDivCub->getNumPoints());
516  numBasisDerivEvalPoints += elemBasisDivCub->getNumPoints();
517  basisDerivCubPoints[dim][0] = view_type("basisDerivCubPoints",elemBasisDivCub->getNumPoints(),dim);
518  basisDerivCubWeights[dim][0] = view_type("basisDerivCubWeights",elemBasisDivCub->getNumPoints());
519  elemBasisDivCub->getCubature(basisDerivCubPoints[dim][0], basisDerivCubWeights[dim][0]);
520  maxNumBasisDerivEvalPoints = std::max(maxNumBasisDerivEvalPoints, elemBasisDivCub->getNumPoints());
521 
522  cub_degree = basisCubDegree - 1 + targetDivCubDegre;
523  auto elemTargetDivCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
524  targetDerivPointsRange(dim,0) = range_type(numTargetDerivEvalPoints, numTargetDerivEvalPoints+elemTargetDivCub->getNumPoints());
525  numTargetDerivEvalPoints += elemTargetDivCub->getNumPoints();
526  targetDerivCubPoints[dim][0] = view_type("targetDerivCubPoints",elemTargetDivCub->getNumPoints(),dim);
527  targetDerivCubWeights[dim][0] = view_type("targetDerivCubWeights",elemTargetDivCub->getNumPoints());
528  elemTargetDivCub->getCubature(targetDerivCubPoints[dim][0], targetDerivCubWeights[dim][0]);
529  maxNumTargetDerivEvalPoints = std::max(maxNumTargetDerivEvalPoints, elemTargetDivCub->getNumPoints());
530 
531  if(haveHCurlConstraint)
532  {
533  cub_degree = 2*basisCubDegree;
534  auto elemBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
535  basisPointsRange(dim,0) = range_type(numBasisEvalPoints, numBasisEvalPoints + elemBasisCub->getNumPoints());
536  numBasisEvalPoints += elemBasisCub->getNumPoints();
537  basisCubPoints[dim][0] = view_type("basisCubPoints",elemBasisCub->getNumPoints(),dim);
538  basisCubWeights[dim][0] = view_type("basisCubWeights",elemBasisCub->getNumPoints());
539  elemBasisCub->getCubature(basisCubPoints[dim][0], basisCubWeights[dim][0]);
540  maxNumBasisEvalPoints = std::max(maxNumBasisEvalPoints, elemBasisCub->getNumPoints());
541 
542  cub_degree = basisCubDegree + targetCubDegree;
543  auto elemTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(subCellTopologyKey(dim,0), cub_degree);
544  targetPointsRange(dim,0) = range_type(numTargetEvalPoints, numTargetEvalPoints + elemTargetCub->getNumPoints());
545  numTargetEvalPoints += elemTargetCub->getNumPoints();
546  targetCubPoints[dim][0] = view_type("targetCubPoints",elemTargetCub->getNumPoints(),dim);
547  targetCubWeights[dim][0] = view_type("targetCubWeights",elemTargetCub->getNumPoints());
548  elemTargetCub->getCubature(targetCubPoints[dim][0], targetCubWeights[dim][0]);
549  maxNumTargetEvalPoints = std::max(maxNumTargetEvalPoints, elemTargetCub->getNumPoints());
550  }
551  }
552 }
553 
554 template<typename DeviceType, typename ValueType>
555 template<typename BasisPtrType>
557  const ordinal_type targetCubDegree) {
558  const auto cellTopo = cellBasis->getBaseCellTopology();
559  ordinal_type dim = cellTopo.getDimension();
560  numBasisEvalPoints = 0;
561  numBasisDerivEvalPoints = 0;
562  numTargetEvalPoints = 0;
563  numTargetDerivEvalPoints = 0;
564 
565  basisPointsRange = range_tag("basisPointsRange", 4,maxSubCellsCount);
566  targetPointsRange = range_tag("targetPointsRange", 4,maxSubCellsCount);
567  basisDerivPointsRange = range_tag("basisDerivPointsRange", 4,maxSubCellsCount);
568  targetDerivPointsRange = range_tag("targetDerivPointsRange", 4,maxSubCellsCount);
569  subCellTopologyKey = key_tag("subCellTopologyKey",4,maxSubCellsCount);
570 
571  ordinal_type basisCubDegree = cellBasis->getDegree();
572  DefaultCubatureFactory cub_factory;
573  subCellTopologyKey(dim,0) = cellBasis->getBaseCellTopology().getBaseKey();
574 
575  maxNumBasisEvalPoints = 0; maxNumTargetEvalPoints =0;
576 
577  ordinal_type cub_degree = 2*basisCubDegree;
578  auto elemBasisCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(cellTopo.getBaseKey(), cub_degree);
579  basisPointsRange(dim,0) = range_type(0, elemBasisCub->getNumPoints());
580  numBasisEvalPoints += elemBasisCub->getNumPoints();
581  maxNumBasisEvalPoints = elemBasisCub->getNumPoints();
582  basisCubPoints[dim][0] = view_type("basisCubPoints",elemBasisCub->getNumPoints(),dim);
583  basisCubWeights[dim][0] = view_type("basisCubWeights",elemBasisCub->getNumPoints());
584  elemBasisCub->getCubature(basisCubPoints[dim][0], basisCubWeights[dim][0]);
585 
586  cub_degree = basisCubDegree + targetCubDegree;
587  auto elemTargetCub = cub_factory.create<HostSpaceType, ValueType, ValueType>(cellTopo.getBaseKey(), cub_degree);
588  targetPointsRange(dim,0) = range_type(0, elemTargetCub->getNumPoints());
589  numTargetEvalPoints += elemTargetCub->getNumPoints();
590  maxNumTargetEvalPoints = elemTargetCub->getNumPoints();
591  targetCubPoints[dim][0] = view_type("targetCubPoints",elemTargetCub->getNumPoints(),dim);
592  targetCubWeights[dim][0] = view_type("targetCubWeights",elemTargetCub->getNumPoints());
593  elemTargetCub->getCubature(targetCubPoints[dim][0], targetCubWeights[dim][0]);
594 }
595 
596 }
597 }
598 #endif
599 
600 
601 
602 
603 
Header file for the abstract base class Intrepid2::DefaultCubatureFactory.
Header file for the Intrepid2::Basis_HCURL_HEX_In_FEM class.
A factory class that generates specific instances of cubatures.
void createHDivProjectionStruct(const BasisPtrType cellBasis, const ordinal_type targetCubDegree, const ordinal_type targetDivCubDegre)
Initialize the ProjectionStruct for HDIV projections.
Implementation of the default H(curl)-compatible Nedelec (first kind) basis of arbitrary degree on Te...
void createL2ProjectionStruct(const BasisPtrType cellBasis, const ordinal_type targetCubDegree)
Initialize the ProjectionStruct for L2 projections.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
Implementation of the default H(curl)-compatible FEM basis on Hexahedron cell.
Header function for Intrepid2::Util class and other utility functions.
Header file for the Intrepid2::Basis_HGRAD_TRI_Cn_FEM class.
Header file for the Intrepid2::Basis_HCURL_TET_In_FEM class.
void createHCurlProjectionStruct(const BasisPtrType cellBasis, const ordinal_type targetCubDegree, const ordinal_type targetCurlCubDegre)
Initialize the ProjectionStruct for HCURL projections.
static Teuchos::RCP< Cubature< DeviceType, pointValueType, weightValueType > > create(unsigned topologyKey, const std::vector< ordinal_type > &degree, const EPolyType polytype=POLYTYPE_MAX)
Factory method.
Header file for the Intrepid2::Basis_HGRAD_QUAD_Cn_FEM class.
ordinal_type getDofCount(const ordinal_type subcDim, const ordinal_type subcOrd) const
DoF count for specified subcell.
static void getReferenceVertex(Kokkos::DynRankView< cellVertexValueType, cellVertexProperties... > cellVertex, const shards::CellTopology cell, const ordinal_type vertexOrd)
Retrieves the Cartesian coordinates of a reference cell vertex.
Implementation of the default H(grad)-compatible Lagrange basis of arbitrary degree on Triangle cell...
void createHVolProjectionStruct(const BasisPtrType cellBasis, const ordinal_type targetCubDegree)
Initialize the ProjectionStruct for HVOL (local-L2) projection.
Implementation of the default H(grad)-compatible FEM basis of degree n on Quadrilateral cell Implemen...
void createHGradProjectionStruct(const BasisPtrType cellBasis, const ordinal_type targetCubDegree, const ordinal_type targetGradCubDegre)
Initialize the ProjectionStruct for HGRAD projections.