Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_BasisValues2_impl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Panzer: A partial differential equation assembly
5 // engine for strongly coupled complex multiphysics systems
6 // Copyright (2011) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Roger P. Pawlowski (rppawlo@sandia.gov) and
39 // Eric C. Cyr (eccyr@sandia.gov)
40 // ***********************************************************************
41 // @HEADER
42 
43 #include "PanzerDiscFE_config.hpp"
44 #include "Panzer_Traits.hpp"
45 
47 #include "Kokkos_ViewFactory.hpp"
48 
49 #include "Intrepid2_Utils.hpp"
50 #include "Intrepid2_FunctionSpaceTools.hpp"
51 #include "Intrepid2_Orientation.hpp"
52 #include "Intrepid2_OrientationTools.hpp"
53 
54 
55 namespace panzer {
56 
57 template <typename Scalar>
59 evaluateValues(const PHX::MDField<Scalar,IP,Dim> & cub_points,
60  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
61  const PHX::MDField<Scalar,Cell,IP> & jac_det,
62  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
63  const int in_num_cells)
64 {
65  PHX::MDField<Scalar,Cell,IP> weighted_measure;
66  PHX::MDField<Scalar,Cell,NODE,Dim> vertex_coordinates;
67  build_weighted = false;
68  evaluateValues(cub_points,jac,jac_det,jac_inv,weighted_measure,vertex_coordinates,false,in_num_cells);
69 }
70 
71 template <typename Scalar>
73 evaluateValues(const PHX::MDField<Scalar,IP,Dim> & cub_points,
74  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
75  const PHX::MDField<Scalar,Cell,IP> & jac_det,
76  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
77  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
78  const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
79  bool use_vertex_coordinates,
80  const int in_num_cells)
81 {
82  MDFieldArrayFactory af("",ddims_,true);
83 
84  int num_dim = basis_layout->dimension();
85 
86  // currently this just copies on the basis objects are converted
87  // in intrepid
88  evaluateReferenceValues(cub_points,compute_derivatives,use_vertex_coordinates);
89 
90  const int num_cells = in_num_cells < 0 ? jac.extent(0) : in_num_cells;
91  const std::pair<int,int> cell_range(0,num_cells);
92 
93  PureBasis::EElementSpace elmtspace = getElementSpace();
94  if(elmtspace==PureBasis::CONST ||
95  elmtspace==PureBasis::HGRAD) {
96  auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
97  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
98  HGRADtransformVALUE(s_basis_scalar,
99  basis_ref_scalar.get_view());
100 
101  if(build_weighted) {
102  auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
103  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
104  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
105  multiplyMeasure(s_weighted_basis_scalar,
106  s_weighted_measure,
107  s_basis_scalar);
108  }
109  }
110  else if(elmtspace==PureBasis::HCURL) {
111  auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
112  auto s_jac_inv = Kokkos::subview(jac_inv.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
113  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
114  HCURLtransformVALUE(s_basis_vector,
115  s_jac_inv,
116  basis_ref_vector.get_view());
117 
118  if(build_weighted) {
119  auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
120  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
121  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
122  multiplyMeasure(s_weighted_basis_vector,
123  s_weighted_measure,
124  s_basis_vector);
125  }
126  }
127  else if(elmtspace==PureBasis::HDIV)
128  {
129  auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
130  auto s_jac = Kokkos::subview(jac.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
131  auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
132  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
133  HDIVtransformVALUE(s_basis_vector,
134  s_jac,
135  s_jac_det,
136  basis_ref_vector.get_view());
137 
138  if(build_weighted) {
139  auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
140  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
141  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
142  multiplyMeasure(s_weighted_basis_vector,
143  s_weighted_measure,
144  s_basis_vector);
145  }
146  }
147  else if(elmtspace==PureBasis::HVOL)
148  {
149  auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
150  auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
151  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
152  HVOLtransformVALUE(s_basis_scalar,
153  s_jac_det,
154  basis_ref_scalar.get_view());
155 
156  if(build_weighted) {
157  auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
158  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
159  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
160  multiplyMeasure(s_weighted_basis_scalar,
161  s_weighted_measure,
162  s_basis_scalar);
163  }
164  }
165  else { TEUCHOS_ASSERT(false); }
166 
167  if(elmtspace==PureBasis::HGRAD && compute_derivatives) {
168  auto s_grad_basis = Kokkos::subview(grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
169  auto s_jac_inv = Kokkos::subview(jac_inv.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
170  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
171  HGRADtransformGRAD(s_grad_basis,
172  s_jac_inv,
173  grad_basis_ref.get_view());
174 
175  if(build_weighted) {
176  auto s_weighted_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
177  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
178  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
179  multiplyMeasure(s_weighted_grad_basis,
180  s_weighted_measure,
181  s_grad_basis);
182  }
183  }
184  else if(elmtspace==PureBasis::HCURL && num_dim==2 && compute_derivatives) {
185  auto s_curl_basis_scalar = Kokkos::subview(curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
186  auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
187  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
188  HDIVtransformDIV(s_curl_basis_scalar,
189  s_jac_det, // note only volume deformation is needed!
190  // this relates directly to this being in
191  // the divergence space in 2D!
192  curl_basis_ref_scalar.get_view());
193 
194  if(build_weighted) {
195  auto s_weighted_curl_basis_scalar = Kokkos::subview(weighted_curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
196  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
197  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
198  multiplyMeasure(s_weighted_curl_basis_scalar,
199  s_weighted_measure,
200  s_curl_basis_scalar);
201  }
202  }
203  else if(elmtspace==PureBasis::HCURL && num_dim==3 && compute_derivatives) {
204  auto s_curl_basis_vector = Kokkos::subview(curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
205  auto s_jac = Kokkos::subview(jac.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
206  auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
207  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
208  HCURLtransformCURL(s_curl_basis_vector,
209  s_jac,
210  s_jac_det,
211  curl_basis_ref_vector.get_view());
212 
213  if(build_weighted) {
214  auto s_weighted_curl_basis_vector = Kokkos::subview(weighted_curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
215  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
216  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
217  multiplyMeasure(s_weighted_curl_basis_vector,
218  s_weighted_measure,
219  s_curl_basis_vector);
220  }
221  }
222  else if(elmtspace==PureBasis::HDIV && compute_derivatives) {
223  auto s_div_basis = Kokkos::subview(div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
224  auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
225  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
226  HDIVtransformDIV(s_div_basis,
227  s_jac_det,
228  div_basis_ref.get_view());
229 
230  if(build_weighted) {
231  auto s_weighted_div_basis = Kokkos::subview(weighted_div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
232  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
233  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::
234  multiplyMeasure(s_weighted_div_basis,
235  s_weighted_measure,
236  s_div_basis);
237  }
238  }
239 
240  // If basis supports coordinate values at basis points, then
241  // compute these values
242  if(use_vertex_coordinates) {
243  // Teuchos::RCP<Intrepid2::DofCoordsInterface<ArrayDynamic> > coords
244  // = Teuchos::rcp_dynamic_cast<Intrepid2::DofCoordsInterface<ArrayDynamic> >(intrepid_basis);
245  // if (!Teuchos::is_null(coords)) {
246 /*
247  ArrayDynamic dyn_basis_coordinates_ref = af.buildArray<Scalar,BASIS,Dim>("basis_coordinates_ref",basis_coordinates_ref.extent(0),basis_coordinates_ref.extent(1));
248  coords->getDofCoords(dyn_basis_coordinates_ref);
249 
250  // fill in basis coordinates
251  for (size_type i = 0; i < basis_coordinates_ref.extent(0); ++i)
252  for (size_type j = 0; j < basis_coordinates_ref.extent(1); ++j)
253  basis_coordinates_ref(i,j) = dyn_basis_coordinates_ref(i,j);
254 */
255 
256  auto s_basis_coordinates = Kokkos::subview(basis_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
257  auto s_vertex_coordinates = Kokkos::subview(vertex_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
258  Intrepid2::CellTools<PHX::Device::execution_space> cell_tools;
259  cell_tools.mapToPhysicalFrame(s_basis_coordinates,
260  basis_coordinates_ref.get_view(),
261  s_vertex_coordinates,
262  intrepid_basis->getBaseCellTopology());
263  }
264 }
265 
266 
267 
268 
269 
270 
271 
272 template <typename Scalar>
274 evaluateBasisCoordinates(const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
275  const int in_num_cells)
276 {
277  MDFieldArrayFactory af("",ddims_,true);
278 
279  // intrepid_basis->getDofCoords requires DynRankView, but basis_coordinates_ref is more of a static View
280  // We use an auxiliary 'dyn' array to get around this
281  using coordsScalarType = typename Intrepid2::Basis<PHX::Device::execution_space,Scalar,Scalar>::scalarType;
282  auto dyn_basis_coordinates_ref = af.buildArray<coordsScalarType,BASIS,Dim>("basis_coordinates_ref",
283  basis_coordinates_ref.extent(0),
284  basis_coordinates_ref.extent(1));
285  intrepid_basis->getDofCoords(dyn_basis_coordinates_ref.get_view());
286 
287  // fill in basis coordinates
288  for (int i = 0; i < basis_coordinates_ref.extent_int(0); ++i)
289  for (int j = 0; j < basis_coordinates_ref.extent_int(1); ++j)
290  basis_coordinates_ref(i,j) = dyn_basis_coordinates_ref(i,j);
291 
292  const int num_cells = in_num_cells < 0 ? vertex_coordinates.extent(0) : in_num_cells;
293  const std::pair<int,int> cell_range(0,num_cells);
294  const auto s_basis_coordinates = Kokkos::subview(basis_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
295  const auto s_vertex_coordinates = Kokkos::subview(vertex_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
296 
297  Intrepid2::CellTools<PHX::Device::execution_space> cell_tools;
298  cell_tools.mapToPhysicalFrame(s_basis_coordinates,
299  basis_coordinates_ref.get_view(),
300  s_vertex_coordinates,
301  intrepid_basis->getBaseCellTopology());
302 }
303 
304 
305 
306 
307 
308 template <typename Scalar>
310 evaluateValues(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
311  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
312  const PHX::MDField<Scalar,Cell,IP> & jac_det,
313  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
314  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
315  const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
316  bool use_vertex_coordinates,
317  const int in_num_cells)
318 {
319 
320  PureBasis::EElementSpace elmtspace = getElementSpace();
321 
322  if(elmtspace == PureBasis::CONST){
323  evaluateValues_Const(cub_points,jac_inv,weighted_measure,in_num_cells);
324  } else if(elmtspace == PureBasis::HVOL){
325  evaluateValues_HVol(cub_points,jac_det,jac_inv,weighted_measure,in_num_cells);
326  } else if(elmtspace == PureBasis::HGRAD){
327  evaluateValues_HGrad(cub_points,jac_inv,weighted_measure,in_num_cells);
328  } else if(elmtspace == PureBasis::HCURL){
329  evaluateValues_HCurl(cub_points,jac,jac_det,jac_inv,weighted_measure,in_num_cells);
330  } else if(elmtspace == PureBasis::HDIV){
331  evaluateValues_HDiv(cub_points,jac,jac_det,weighted_measure,in_num_cells);
332  } else {
333  TEUCHOS_TEST_FOR_EXCEPT_MSG(true,"panzer::BasisValues2::evaluateValues : Element space not recognized.");
334  }
335 
336  if(use_vertex_coordinates) {
337  TEUCHOS_TEST_FOR_EXCEPT_MSG(elmtspace == PureBasis::CONST,"panzer::BasisValues2::evaluateValues : Const basis cannot have basis coordinates.");
338  evaluateBasisCoordinates(vertex_coordinates);
339  }
340 
341 }
342 
343 template <typename Scalar>
345 evaluateValues_Const(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
346  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
347  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
348  const int in_num_cells)
349 {
350 
351  TEUCHOS_ASSERT(getElementSpace() == PureBasis::CONST);
352 
353  typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
354  MDFieldArrayFactory af("",ddims_,true);
355 
356  const panzer::PureBasis & basis = *(basis_layout->getBasis());
357 
358  const int num_points = basis_layout->numPoints();
359  const int num_basis = basis.cardinality();
360  const int num_dim = basis_layout->dimension();
361  const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
362 
363  auto cell_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("cell_basis_scalar",1,num_basis,num_points);
364  auto cell_cub_points = af.buildStaticArray<Scalar,IP,Dim>("cell_cub_points",num_points,num_dim);
365  auto cell_grad_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("cell_grad_basis",1,num_basis,num_points,num_dim);
366  auto cell_jac_inv = af.buildStaticArray<Scalar,Cell,IP,Dim,Dim>("cell_jac_inv",1,num_points,num_dim,num_dim);
367 
368  auto cell_basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("cell_basis_ref_scalar",num_basis,num_points);
369  auto cell_grad_basis_ref = af.buildStaticArray<Scalar,BASIS,IP,Dim>("cell_grad_basis_ref",num_basis,num_points,num_dim);
370 
371  for(int cell=0;cell<num_cells;++cell){
372 
373  // =============================================
374  // Load external into cell-local arrays
375 
376  for(int p=0;p<num_points;++p)
377  for(int d=0;d<num_dim;++d)
378  for(int d2=0;d2<num_dim;++d2)
379  cell_jac_inv(0,p,d,d2)=jac_inv(cell,p,d,d2);
380  for(int p=0;p<num_points;++p)
381  for(int d=0;d<num_dim;++d)
382  cell_cub_points(p,d)=cub_points(cell,p,d);
383 
384  // =============================================
385  // Load Reference Values
386 
387  intrepid_basis->getValues(cell_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
388 
389  if(compute_derivatives){
390  Kokkos::deep_copy(cell_grad_basis_ref.get_view(),0.0);
391  }
392 
393  // =============================================
394  // Transform reference values to physical values
395 
396  fst::HGRADtransformVALUE(cell_basis_scalar.get_view(),cell_basis_ref_scalar.get_view());
397  for(int b=0;b<num_basis;++b)
398  for(int p=0;p<num_points;++p)
399  basis_scalar(cell,b,p)=cell_basis_scalar(0,b,p);
400 
401  if(compute_derivatives){
402  fst::HGRADtransformGRAD(cell_grad_basis.get_view(),cell_jac_inv.get_view(),cell_grad_basis_ref.get_view());
403  for(int b=0;b<num_basis;++b)
404  for(int p=0;p<num_points;++p)
405  for(int d=0;d<num_dim;++d)
406  grad_basis(cell,b,p,d)=cell_grad_basis(0,b,p,d);
407  }
408  // =============================================
409  }
410 
411 
412  if(build_weighted){
413  const std::pair<int,int> cell_range(0,num_cells);
414  auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
415  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
416  auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
417  fst::multiplyMeasure(s_weighted_basis_scalar,s_weighted_measure,s_basis_scalar);
418  if(compute_derivatives){
419  auto s_weighted_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
420  auto s_grad_basis = Kokkos::subview(grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
421  fst::multiplyMeasure(s_weighted_grad_basis,s_weighted_measure,s_grad_basis);
422  }
423  }
424 
425 
426 }
427 
428 template <typename Scalar>
430 evaluateValues_HVol(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
431  const PHX::MDField<Scalar,Cell,IP> & jac_det,
432  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & /* jac_inv */,
433  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
434  const int in_num_cells)
435 {
436 
437  TEUCHOS_ASSERT(getElementSpace() == PureBasis::HVOL);
438 
439  typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
440  MDFieldArrayFactory af("",ddims_,true);
441 
442  const panzer::PureBasis & basis = *(basis_layout->getBasis());
443 
444  const int num_points = basis_layout->numPoints();
445  const int num_basis = basis.cardinality();
446  const int num_dim = basis_layout->dimension();
447  const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
448 
449  auto cell_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("cell_basis_scalar",1,num_basis,num_points);
450  auto cell_cub_points = af.buildStaticArray<Scalar,IP,Dim>("cell_cub_points",num_points,num_dim);
451  auto cell_jac_det = af.buildStaticArray<Scalar,Cell,IP>("cell_jac_det",1,num_points);
452 
453  auto cell_basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("cell_basis_ref_scalar",num_basis,num_points);
454 
455  for(int cell=0;cell<num_cells;++cell){
456 
457  // =============================================
458  // Load external into cell-local arrays
459 
460  for(int p=0;p<num_points;++p)
461  cell_jac_det(0,p)=jac_det(cell,p);
462  for(int p=0;p<num_points;++p)
463  for(int d=0;d<num_dim;++d)
464  cell_cub_points(p,d)=cub_points(cell,p,d);
465 
466  // =============================================
467  // Load Reference Values
468 
469  intrepid_basis->getValues(cell_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
470 
471  // =============================================
472  // Transform reference values to physical values
473 
474  fst::HVOLtransformVALUE(cell_basis_scalar.get_view(),cell_jac_det.get_view(),cell_basis_ref_scalar.get_view());
475  for(int b=0;b<num_basis;++b)
476  for(int p=0;p<num_points;++p)
477  basis_scalar(cell,b,p)=cell_basis_scalar(0,b,p);
478 
479  // =============================================
480  }
481 
482 
483  if(build_weighted) {
484  const std::pair<int,int> cell_range(0,num_cells);
485  auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
486  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
487  auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
488  fst::multiplyMeasure(s_weighted_basis_scalar,s_weighted_measure,s_basis_scalar);
489  }
490 
491 }
492 
493 template <typename Scalar>
495 evaluateValues_HGrad(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
496  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
497  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
498  const int in_num_cells)
499 {
500 
501  TEUCHOS_ASSERT(getElementSpace() == PureBasis::HGRAD);
502 
503  typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
504  MDFieldArrayFactory af("",ddims_,true);
505 
506  const panzer::PureBasis & basis = *(basis_layout->getBasis());
507 
508  const int num_points = basis_layout->numPoints();
509  const int num_basis = basis.cardinality();
510  const int num_dim = basis_layout->dimension();
511  const int num_cells = in_num_cells < 0 ? cub_points.extent(0) : in_num_cells;
512 
513  auto cell_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("cell_basis_scalar",1,num_basis,num_points);
514  auto cell_cub_points = af.buildStaticArray<Scalar,IP,Dim>("cell_cub_points",num_points,num_dim);
515  auto cell_grad_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("cell_grad_basis",1,num_basis,num_points,num_dim);
516  auto cell_jac_inv = af.buildStaticArray<Scalar,Cell,IP,Dim,Dim>("cell_jac_inv",1,num_points,num_dim,num_dim);
517 
518  auto cell_basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("cell_basis_ref_scalar",num_basis,num_points);
519  auto cell_grad_basis_ref = af.buildStaticArray<Scalar,BASIS,IP,Dim>("cell_grad_basis_ref",num_basis,num_points,num_dim);
520 
521  for(int cell=0;cell<num_cells;++cell){
522 
523  // =============================================
524  // Load external into cell-local arrays
525 
526  for(int p=0;p<num_points;++p)
527  for(int d=0;d<num_dim;++d)
528  for(int d2=0;d2<num_dim;++d2)
529  cell_jac_inv(0,p,d,d2)=jac_inv(cell,p,d,d2);
530  for(int p=0;p<num_points;++p)
531  for(int d=0;d<num_dim;++d)
532  cell_cub_points(p,d)=cub_points(cell,p,d);
533 
534  // =============================================
535  // Load Reference Values
536 
537  intrepid_basis->getValues(cell_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
538 
539  if(compute_derivatives){
540  intrepid_basis->getValues(cell_grad_basis_ref.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_GRAD);
541  }
542 
543  // =============================================
544  // Transform reference values to physical values
545 
546  fst::HGRADtransformVALUE(cell_basis_scalar.get_view(),cell_basis_ref_scalar.get_view());
547  for(int b=0;b<num_basis;++b)
548  for(int p=0;p<num_points;++p)
549  basis_scalar(cell,b,p)=cell_basis_scalar(0,b,p);
550 
551  if(compute_derivatives){
552  fst::HGRADtransformGRAD(cell_grad_basis.get_view(),cell_jac_inv.get_view(),cell_grad_basis_ref.get_view());
553  for(int b=0;b<num_basis;++b)
554  for(int p=0;p<num_points;++p)
555  for(int d=0;d<num_dim;++d)
556  grad_basis(cell,b,p,d)=cell_grad_basis(0,b,p,d);
557  }
558  // =============================================
559  }
560 
561  if(build_weighted){
562  const std::pair<int,int> cell_range(0,num_cells);
563  auto s_weighted_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
564  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
565  auto s_basis_scalar = Kokkos::subview(basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
566  fst::multiplyMeasure(s_weighted_basis_scalar,s_weighted_measure,s_basis_scalar);
567  if(compute_derivatives){
568  auto s_weighted_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
569  auto s_grad_basis = Kokkos::subview(grad_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
570  fst::multiplyMeasure(s_weighted_grad_basis,s_weighted_measure,s_grad_basis);
571  }
572  }
573 
574 }
575 
576 
577 template <typename Scalar>
579 evaluateValues_HCurl(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
580  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
581  const PHX::MDField<Scalar,Cell,IP> & jac_det,
582  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
583  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
584  const int in_num_cells)
585 {
586 
587  TEUCHOS_ASSERT(getElementSpace() == PureBasis::HCURL);
588 
589 
590  typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
591  MDFieldArrayFactory af("",ddims_,true);
592 
593  const panzer::PureBasis & basis = *(basis_layout->getBasis());
594 
595  const int num_points = basis_layout->numPoints();
596  const int num_basis = basis.cardinality();
597  const int num_dim = basis_layout->dimension();
598  const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
599 
600  auto cell_cub_points = af.buildStaticArray<Scalar,IP,Dim>("cell_cub_points",num_points,num_dim);
601  auto cell_jac = af.buildStaticArray<Scalar,Cell,IP,Dim,Dim>("cell_jac",1,num_points,num_dim,num_dim);
602  auto cell_jac_inv = af.buildStaticArray<Scalar,Cell,IP,Dim,Dim>("cell_jac_inv",1,num_points,num_dim,num_dim);
603  auto cell_jac_det = af.buildStaticArray<Scalar,Cell,IP>("cell_jac_det",1,num_points);
604 
605  auto cell_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("cell_basis_vector",1,num_basis,num_points,num_dim);
606  auto cell_curl_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("cell_curl_basis_scalar",1,num_basis,num_points);
607  auto cell_curl_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("cell_curl_basis_vector",1,num_basis,num_points,num_dim);
608 
609  auto cell_curl_basis_ref = af.buildArray<Scalar,BASIS,IP,Dim>("cell_curl_basis_ref",num_basis,num_points,num_dim);
610  auto cell_curl_basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("cell_curl_basis_ref_scalar",num_basis,num_points);
611  auto cell_basis_ref_vector = af.buildArray<Scalar,BASIS,IP,Dim>("cell_basis_ref_vector",num_basis,num_points,num_dim);
612 
613  for(int cell=0;cell<num_cells;++cell){
614 
615  // =============================================
616  // Load external into cell-local arrays
617 
618  for(int p=0;p<num_points;++p)
619  for(int d=0;d<num_dim;++d)
620  for(int d2=0;d2<num_dim;++d2)
621  cell_jac(0,p,d,d2)=jac(cell,p,d,d2);
622  for(int p=0;p<num_points;++p)
623  for(int d=0;d<num_dim;++d)
624  for(int d2=0;d2<num_dim;++d2)
625  cell_jac_inv(0,p,d,d2)=jac_inv(cell,p,d,d2);
626  for(int p=0;p<num_points;++p)
627  cell_jac_det(0,p)=jac_det(cell,p);
628  for(int p=0;p<num_points;++p)
629  for(int d=0;d<num_dim;++d)
630  cell_cub_points(p,d)=cub_points(cell,p,d);
631 
632  // =============================================
633  // Load Reference Values
634 
635  intrepid_basis->getValues(cell_basis_ref_vector.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
636 
637  if(compute_derivatives){
638  if(num_dim==2){
639  intrepid_basis->getValues(cell_curl_basis_ref_scalar.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_CURL);
640  } else if(num_dim==3){
641  intrepid_basis->getValues(cell_curl_basis_ref.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_CURL);
642  }
643  }
644 
645  // =============================================
646  // Transform reference values to physical values
647 
648  fst::HCURLtransformVALUE(cell_basis_vector.get_view(),cell_jac_inv.get_view(),cell_basis_ref_vector.get_view());
649  for(int b=0;b<num_basis;++b)
650  for(int p=0;p<num_points;++p)
651  for(int d=0;d<num_dim;++d)
652  basis_vector(cell,b,p,d)=cell_basis_vector(0,b,p,d);
653 
654  if(compute_derivatives){
655  if(num_dim==2){
656  // note only volume deformation is needed!
657  // this relates directly to this being in
658  // the divergence space in 2D!
659  fst::HDIVtransformDIV(cell_curl_basis_scalar.get_view(),cell_jac_det.get_view(),cell_curl_basis_ref_scalar.get_view());
660  for(int b=0;b<num_basis;++b)
661  for(int p=0;p<num_points;++p)
662  curl_basis_scalar(cell,b,p)=cell_curl_basis_scalar(0,b,p);
663  } else if(num_dim==3) {
664  fst::HCURLtransformCURL(cell_curl_basis_vector.get_view(),cell_jac.get_view(),cell_jac_det.get_view(),cell_curl_basis_ref.get_view());
665  for(int b=0;b<num_basis;++b)
666  for(int p=0;p<num_points;++p)
667  for(int d=0;d<num_dim;++d)
668  curl_basis_vector(cell,b,p,d)=cell_curl_basis_vector(0,b,p,d);
669  } else {
670  TEUCHOS_TEST_FOR_EXCEPT_MSG(true,"panzer::BasisValues2::evaluateValues_HCurl : HCurl only setup for 2D and 3D.");
671  }
672  }
673  }
674 
675  if(build_weighted){
676  const std::pair<int,int> cell_range(0,num_cells);
677  auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
678  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
679  auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
680  fst::multiplyMeasure(s_weighted_basis_vector,s_weighted_measure,s_basis_vector);
681  if(compute_derivatives){
682  if(num_dim==2){
683  auto s_weighted_curl_basis_scalar = Kokkos::subview(weighted_curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
684  auto s_curl_basis_scalar = Kokkos::subview(curl_basis_scalar.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
685  fst::multiplyMeasure(s_weighted_curl_basis_scalar,s_weighted_measure,s_curl_basis_scalar);
686  } else if(num_dim==3){
687  auto s_weighted_curl_basis_vector = Kokkos::subview(weighted_curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
688  auto s_curl_basis_vector = Kokkos::subview(curl_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
689  fst::multiplyMeasure(s_weighted_curl_basis_vector,s_weighted_measure,s_curl_basis_vector);
690  }
691  }
692  }
693 
694 }
695 
696 template <typename Scalar>
698 evaluateValues_HDiv(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
699  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
700  const PHX::MDField<Scalar,Cell,IP> & jac_det,
701  const PHX::MDField<Scalar,Cell,IP> & weighted_measure,
702  const int in_num_cells)
703 {
704 
705  TEUCHOS_ASSERT(getElementSpace() == PureBasis::HDIV);
706 
707  typedef Intrepid2::FunctionSpaceTools<PHX::Device::execution_space> fst;
708  MDFieldArrayFactory af("",ddims_,true);
709 
710  const panzer::PureBasis & basis = *(basis_layout->getBasis());
711 
712  const int num_points = basis_layout->numPoints();
713  const int num_basis = basis.cardinality();
714  const int num_dim = basis_layout->dimension();
715  const int num_cells = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
716 
717  auto cell_cub_points = af.buildStaticArray<Scalar,IP,Dim>("cell_cub_points",num_points,num_dim);
718  auto cell_jac = af.buildStaticArray<Scalar,Cell,IP,Dim,Dim>("cell_jac",1,num_points,num_dim,num_dim);
719  auto cell_jac_det = af.buildStaticArray<Scalar,Cell,IP>("cell_jac_det",1,num_points);
720 
721  auto cell_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("cell_basis_vector",1,num_basis,num_points,num_dim);
722  auto cell_div_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP>("cell_div_basis",1,num_basis,num_points);
723 
724  auto cell_basis_ref_vector = af.buildArray<Scalar,BASIS,IP,Dim>("cell_basis_ref_vector",num_basis,num_points,num_dim);
725  auto cell_div_basis_ref = af.buildStaticArray<Scalar,BASIS,IP>("cell_div_basis_ref",num_basis,num_points);
726 
727  for(int cell=0;cell<num_cells;++cell){
728 
729  // =============================================
730  // Load external into cell-local arrays
731 
732  for(int p=0;p<num_points;++p)
733  for(int d=0;d<num_dim;++d)
734  for(int d2=0;d2<num_dim;++d2)
735  cell_jac(0,p,d,d2)=jac(cell,p,d,d2);
736  for(int p=0;p<num_points;++p)
737  cell_jac_det(0,p)=jac_det(cell,p);
738  for(int p=0;p<num_points;++p)
739  for(int d=0;d<num_dim;++d)
740  cell_cub_points(p,d)=cub_points(cell,p,d);
741  // =============================================
742  // Load Reference Values
743 
744  intrepid_basis->getValues(cell_basis_ref_vector.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_VALUE);
745 
746  if(compute_derivatives){
747  intrepid_basis->getValues(cell_div_basis_ref.get_view(),cell_cub_points.get_view(),Intrepid2::OPERATOR_DIV);
748  }
749 
750  // =============================================
751  // Transform reference values to physical values
752 
753  fst::HDIVtransformVALUE(cell_basis_vector.get_view(),cell_jac.get_view(),cell_jac_det.get_view(),cell_basis_ref_vector.get_view());
754  for(int b=0;b<num_basis;++b)
755  for(int p=0;p<num_points;++p)
756  for(int d=0;d<num_dim;++d)
757  basis_vector(cell,b,p,d)=cell_basis_vector(0,b,p,d);
758 
759  if(compute_derivatives){
760  fst::HDIVtransformDIV(cell_div_basis.get_view(),cell_jac_det.get_view(),cell_div_basis_ref.get_view());
761  for(int b=0;b<num_basis;++b)
762  for(int p=0;p<num_points;++p)
763  div_basis(cell,b,p)=cell_div_basis(0,b,p);
764  }
765  }
766 
767  if(build_weighted){
768  const std::pair<int,int> cell_range(0,num_cells);
769  auto s_weighted_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
770  auto s_weighted_measure = Kokkos::subview(weighted_measure.get_view(),cell_range,Kokkos::ALL());
771  auto s_basis_vector = Kokkos::subview(basis_vector.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
772  fst::multiplyMeasure(s_weighted_basis_vector,s_weighted_measure,s_basis_vector);
773  if(compute_derivatives){
774  auto s_weighted_div_basis = Kokkos::subview(weighted_div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
775  auto s_div_basis = Kokkos::subview(div_basis.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
776  fst::multiplyMeasure(s_weighted_div_basis,s_weighted_measure,s_div_basis);
777  }
778  }
779 
780 }
781 
782 template <typename Scalar>
784 evaluateValuesCV(const PHX::MDField<Scalar,Cell,IP,Dim> & cub_points,
785  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
786  const PHX::MDField<Scalar,Cell,IP> & jac_det,
787  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv)
788 {
789 
790  PHX::MDField<Scalar,Cell,NODE,Dim> vertex_coordinates;
791  const int in_num_cells = jac.extent(0);
792  evaluateValuesCV(cub_points,jac,jac_det,jac_inv,vertex_coordinates,false,in_num_cells);
793 
794 }
795 
796 template <typename Scalar>
798 evaluateValuesCV(const PHX::MDField<Scalar,Cell,IP,Dim> & cell_cub_points,
799  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac,
800  const PHX::MDField<Scalar,Cell,IP> & jac_det,
801  const PHX::MDField<Scalar,Cell,IP,Dim,Dim> & jac_inv,
802  const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
803  bool use_vertex_coordinates,
804  const int in_num_cells)
805 {
806  MDFieldArrayFactory af("",ddims_,true);
807 
808  int num_ip = basis_layout->numPoints();
809  int num_card = basis_layout->cardinality();
810  int num_dim = basis_layout->dimension();
811 
812  size_type num_cells = jac.extent(0);
813 
814  PureBasis::EElementSpace elmtspace = getElementSpace();
815  ArrayDynamic dyn_cub_points = af.buildArray<Scalar,IP,Dim>("dyn_cub_points", num_ip, num_dim);
816 
817  // Integration points are located on physical cells rather than reference cells,
818  // so we evaluate the basis in a loop over cells.
819  for (size_type icell = 0; icell < num_cells; ++icell)
820  {
821  for (int ip = 0; ip < num_ip; ++ip)
822  for (int d = 0; d < num_dim; ++d)
823  dyn_cub_points(ip,d) = cell_cub_points(icell,ip,d);
824 
825  if(elmtspace==PureBasis::CONST) {
826  ArrayDynamic dyn_basis_ref_scalar = af.buildArray<Scalar,BASIS,IP>("dyn_basis_ref_scalar",num_card,num_ip);
827 
828  intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
829  dyn_cub_points.get_view(),
830  Intrepid2::OPERATOR_VALUE);
831 
832  // transform values method just transfers values to array with cell index - no need to call
833  for (int b = 0; b < num_card; ++b)
834  for (int ip = 0; ip < num_ip; ++ip)
835  basis_scalar(icell,b,ip) = dyn_basis_ref_scalar(b,ip);
836 
837  }
838  if(elmtspace==PureBasis::HVOL) {
839  ArrayDynamic dyn_basis_ref_scalar = af.buildArray<Scalar,BASIS,IP>("dyn_basis_ref_scalar",num_card,num_ip);
840 
841  intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
842  dyn_cub_points.get_view(),
843  Intrepid2::OPERATOR_VALUE);
844 
845  int one_cell= 1;
846  ArrayDynamic dyn_basis_scalar = af.buildArray<Scalar,Cell,BASIS,IP>("dyn_basis_vector",one_cell,num_card,num_ip);
847  ArrayDynamic dyn_jac = af.buildArray<Scalar,Cell,IP,Dim,Dim>("dyn_jac",one_cell,num_ip,num_dim,num_dim);
848  ArrayDynamic dyn_jac_det = af.buildArray<Scalar,Cell,IP>("dyn_jac_det",one_cell,num_ip);
849 
850  int cellInd = 0;
851  for (int ip = 0; ip < num_ip; ++ip)
852  dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
853 
854  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HVOLtransformVALUE(dyn_basis_scalar.get_view(),
855  dyn_jac_det.get_view(),
856  dyn_basis_ref_scalar.get_view());
857 
858  for (int b = 0; b < num_card; ++b)
859  for (int ip = 0; ip < num_ip; ++ip)
860  basis_scalar(icell,b,ip) = dyn_basis_scalar(0,b,ip);
861 
862  }
863  if(elmtspace==PureBasis::HGRAD) {
864  ArrayDynamic dyn_basis_ref_scalar = af.buildArray<Scalar,BASIS,IP>("dyn_basis_ref_scalar",num_card,num_ip);
865 
866  intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
867  dyn_cub_points.get_view(),
868  Intrepid2::OPERATOR_VALUE);
869 
870  // transform values method just transfers values to array with cell index - no need to call
871  for (int b = 0; b < num_card; ++b)
872  for (int ip = 0; ip < num_ip; ++ip)
873  basis_scalar(icell,b,ip) = dyn_basis_ref_scalar(b,ip);
874 
875  if(compute_derivatives) {
876 
877  int one_cell = 1;
878  ArrayDynamic dyn_grad_basis_ref = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_grad_basis_ref",num_card,num_ip,num_dim);
879  ArrayDynamic dyn_grad_basis = af.buildArray<Scalar,Cell,BASIS,IP,Dim>("dyn_grad_basis",one_cell,num_card,num_ip,num_dim);
880  ArrayDynamic dyn_jac_inv = af.buildArray<Scalar,Cell,IP,Dim,Dim>("dyn_jac_inv",one_cell,num_ip,num_dim,num_dim);
881 
882  intrepid_basis->getValues(dyn_grad_basis_ref.get_view(),
883  dyn_cub_points.get_view(),
884  Intrepid2::OPERATOR_GRAD);
885 
886  int cellInd = 0;
887  for (int ip = 0; ip < num_ip; ++ip)
888  for (int d1 = 0; d1 < num_dim; ++d1)
889  for (int d2 = 0; d2 < num_dim; ++d2)
890  dyn_jac_inv(cellInd,ip,d1,d2) = jac_inv(icell,ip,d1,d2);
891 
892  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HGRADtransformGRAD<Scalar>(dyn_grad_basis.get_view(),
893  dyn_jac_inv.get_view(),
894  dyn_grad_basis_ref.get_view());
895 
896  for (int b = 0; b < num_card; ++b)
897  for (int ip = 0; ip < num_ip; ++ip)
898  for (int d = 0; d < num_dim; ++d)
899  grad_basis(icell,b,ip,d) = dyn_grad_basis(0,b,ip,d);
900 
901  }
902  }
903  else if(elmtspace==PureBasis::HCURL) {
904  ArrayDynamic dyn_basis_ref_vector = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_basis_ref_vector",num_card,num_ip,num_dim);
905 
906  intrepid_basis->getValues(dyn_basis_ref_vector.get_view(),
907  dyn_cub_points.get_view(),
908  Intrepid2::OPERATOR_VALUE);
909 
910  const int one_cell = 1;
911  ArrayDynamic dyn_basis_vector = af.buildArray<Scalar,Cell,BASIS,IP,Dim>("dyn_basis_vector",one_cell,num_card,num_ip,num_dim);
912  ArrayDynamic dyn_jac_inv = af.buildArray<Scalar,Cell,IP,Dim,Dim>("dyn_jac_inv",one_cell,num_ip,num_dim,num_dim);
913 
914  const int cellInd = 0;
915  for (int ip = 0; ip < num_ip; ++ip)
916  for (int d1 = 0; d1 < num_dim; ++d1)
917  for (int d2 = 0; d2 < num_dim; ++d2)
918  dyn_jac_inv(cellInd,ip,d1,d2) = jac_inv(icell,ip,d1,d2);
919 
920  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HCURLtransformVALUE(dyn_basis_vector.get_view(),
921  dyn_jac_inv.get_view(),
922  dyn_basis_ref_vector.get_view());
923 
924  for (int b = 0; b < num_card; ++b)
925  for (int ip = 0; ip < num_ip; ++ip)
926  for (int d = 0; d < num_dim; ++d)
927  basis_vector(icell,b,ip,d) = dyn_basis_vector(0,b,ip,d);
928 
929  if(compute_derivatives && num_dim ==2) {
930 
931  ArrayDynamic dyn_curl_basis_ref_scalar = af.buildArray<Scalar,BASIS,IP>("dyn_curl_basis_ref_scalar",num_card,num_ip);
932  ArrayDynamic dyn_curl_basis_scalar = af.buildArray<Scalar,Cell,BASIS,IP>("dyn_curl_basis_scalar",one_cell,num_card,num_ip);
933  ArrayDynamic dyn_jac_det = af.buildArray<Scalar,Cell,IP>("dyn_jac_det",one_cell,num_ip);
934 
935  intrepid_basis->getValues(dyn_curl_basis_ref_scalar.get_view(),
936  dyn_cub_points.get_view(),
937  Intrepid2::OPERATOR_CURL);
938 
939  for (int ip = 0; ip < num_ip; ++ip)
940  dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
941 
942  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HDIVtransformDIV(dyn_curl_basis_scalar.get_view(),
943  dyn_jac_det.get_view(),
944  dyn_curl_basis_ref_scalar.get_view());
945 
946  for (int b = 0; b < num_card; ++b)
947  for (int ip = 0; ip < num_ip; ++ip)
948  curl_basis_scalar(icell,b,ip) = dyn_curl_basis_scalar(0,b,ip);
949 
950  }
951  if(compute_derivatives && num_dim ==3) {
952 
953  ArrayDynamic dyn_curl_basis_ref = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_curl_basis_ref_vector",num_card,num_ip,num_dim);
954  ArrayDynamic dyn_curl_basis = af.buildArray<Scalar,Cell,BASIS,IP,Dim>("dyn_curl_basis_vector",one_cell,num_card,num_ip,num_dim);
955  ArrayDynamic dyn_jac_det = af.buildArray<Scalar,Cell,IP>("dyn_jac_det",one_cell,num_ip);
956  ArrayDynamic dyn_jac = af.buildArray<Scalar,Cell,IP,Dim,Dim>("dyn_jac",one_cell,num_ip,num_dim,num_dim);
957 
958  intrepid_basis->getValues(dyn_curl_basis_ref.get_view(),
959  dyn_cub_points.get_view(),
960  Intrepid2::OPERATOR_CURL);
961 
962  for (int ip = 0; ip < num_ip; ++ip)
963  {
964  dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
965  for (int d1 = 0; d1 < num_dim; ++d1)
966  for (int d2 = 0; d2 < num_dim; ++d2)
967  dyn_jac(cellInd,ip,d1,d2) = jac(icell,ip,d1,d2);
968  }
969 
970  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HCURLtransformCURL(dyn_curl_basis.get_view(),
971  dyn_jac.get_view(),
972  dyn_jac_det.get_view(),
973  dyn_curl_basis_ref.get_view());
974 
975  for (int b = 0; b < num_card; ++b)
976  for (int ip = 0; ip < num_ip; ++ip)
977  for (int d = 0; d < num_dim; ++d)
978  curl_basis_vector(icell,b,ip,d) = dyn_curl_basis(0,b,ip,d);
979 
980  }
981 
982  }
983  else if(elmtspace==PureBasis::HDIV) {
984 
985  ArrayDynamic dyn_basis_ref_vector = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_basis_ref_vector",num_card,num_ip,num_dim);
986 
987  intrepid_basis->getValues(dyn_basis_ref_vector.get_view(),
988  dyn_cub_points.get_view(),
989  Intrepid2::OPERATOR_VALUE);
990 
991  int one_cell= 1;
992  ArrayDynamic dyn_basis_vector = af.buildArray<Scalar,Cell,BASIS,IP,Dim>("dyn_basis_vector",one_cell,num_card,num_ip,num_dim);
993  ArrayDynamic dyn_jac = af.buildArray<Scalar,Cell,IP,Dim,Dim>("dyn_jac",one_cell,num_ip,num_dim,num_dim);
994  ArrayDynamic dyn_jac_det = af.buildArray<Scalar,Cell,IP>("dyn_jac_det",one_cell,num_ip);
995 
996  int cellInd = 0;
997  for (int ip = 0; ip < num_ip; ++ip)
998  {
999  dyn_jac_det(cellInd,ip) = jac_det(icell,ip);
1000  for (int d1 = 0; d1 < num_dim; ++d1)
1001  for (int d2 = 0; d2 < num_dim; ++d2)
1002  dyn_jac(cellInd,ip,d1,d2) = jac(icell,ip,d1,d2);
1003  }
1004 
1005  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HDIVtransformVALUE(dyn_basis_vector.get_view(),
1006  dyn_jac.get_view(),
1007  dyn_jac_det.get_view(),
1008  dyn_basis_ref_vector.get_view());
1009 
1010  for (int b = 0; b < num_card; ++b)
1011  for (int ip = 0; ip < num_ip; ++ip)
1012  for (int d = 0; d < num_dim; ++d)
1013  basis_vector(icell,b,ip,d) = dyn_basis_vector(0,b,ip,d);
1014 
1015  if(compute_derivatives) {
1016 
1017  ArrayDynamic dyn_div_basis_ref = af.buildArray<Scalar,BASIS,IP>("dyn_div_basis_ref_scalar",num_card,num_ip);
1018  ArrayDynamic dyn_div_basis = af.buildArray<Scalar,Cell,BASIS,IP>("dyn_div_basis_scalar",one_cell,num_card,num_ip);
1019 
1020  intrepid_basis->getValues(dyn_div_basis_ref.get_view(),
1021  dyn_cub_points.get_view(),
1022  Intrepid2::OPERATOR_DIV);
1023 
1024  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::HDIVtransformDIV<Scalar>(dyn_div_basis.get_view(),
1025  dyn_jac_det.get_view(),
1026  dyn_div_basis_ref.get_view());
1027 
1028  for (int b = 0; b < num_card; ++b)
1029  for (int ip = 0; ip < num_ip; ++ip)
1030  div_basis(icell,b,ip) = dyn_div_basis(0,b,ip);
1031 
1032  }
1033 
1034  }
1035  else { TEUCHOS_ASSERT(false); }
1036 
1037  } // cell loop
1038 
1039  if(use_vertex_coordinates) {
1040  TEUCHOS_TEST_FOR_EXCEPT_MSG(elmtspace == PureBasis::CONST,"panzer::BasisValues2::evaluateValues : Const basis cannot have basis coordinates.");
1041  evaluateBasisCoordinates(vertex_coordinates);
1042  }
1043 
1044 }
1045 
1046 template <typename Scalar>
1048 evaluateReferenceValues(const PHX::MDField<Scalar,IP,Dim> & cub_points,bool in_compute_derivatives,bool use_vertex_coordinates)
1049 {
1050  MDFieldArrayFactory af("",ddims_,true);
1051 
1052  int num_quad = basis_layout->numPoints();
1053  int num_dim = basis_layout->dimension();
1054  int num_card = basis_layout->cardinality();
1055 
1056  ArrayDynamic dyn_cub_points = af.buildArray<Scalar,IP,Dim>("dyn_cub_points", num_quad,num_dim);
1057 
1058  for (int ip = 0; ip < num_quad; ++ip)
1059  for (int d = 0; d < num_dim; ++d)
1060  dyn_cub_points(ip,d) = cub_points(ip,d);
1061 
1062  PureBasis::EElementSpace elmtspace = getElementSpace();
1063  if(elmtspace==PureBasis::HGRAD || elmtspace==PureBasis::CONST || elmtspace==PureBasis::HVOL) {
1064  ArrayDynamic dyn_basis_ref_scalar = af.buildArray<Scalar,BASIS,IP>("dyn_basis_ref_scalar",num_card,num_quad);
1065 
1066  intrepid_basis->getValues(dyn_basis_ref_scalar.get_view(),
1067  dyn_cub_points.get_view(),
1068  Intrepid2::OPERATOR_VALUE);
1069 
1070  for (int b = 0; b < num_card; ++b)
1071  for (int ip = 0; ip < num_quad; ++ip)
1072  basis_ref_scalar(b,ip) = dyn_basis_ref_scalar(b,ip);
1073  }
1074  else if(elmtspace==PureBasis::HDIV || elmtspace==PureBasis::HCURL) {
1075  ArrayDynamic dyn_basis_ref_vector = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_basis_ref_vector",num_card,num_quad,num_dim);
1076 
1077  intrepid_basis->getValues(dyn_basis_ref_vector.get_view(),
1078  dyn_cub_points.get_view(),
1079  Intrepid2::OPERATOR_VALUE);
1080 
1081  for (int b = 0; b < num_card; ++b)
1082  for (int ip = 0; ip < num_quad; ++ip)
1083  for (int d = 0; d < num_dim; ++d)
1084  basis_ref_vector(b,ip,d) = dyn_basis_ref_vector(b,ip,d);
1085  }
1086  else { TEUCHOS_ASSERT(false); }
1087 
1088  if(elmtspace==PureBasis::HGRAD && in_compute_derivatives) {
1089  ArrayDynamic dyn_grad_basis_ref = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_basis_ref_vector",num_card,num_quad,num_dim);
1090 
1091  intrepid_basis->getValues(dyn_grad_basis_ref.get_view(),
1092  dyn_cub_points.get_view(),
1093  Intrepid2::OPERATOR_GRAD);
1094 
1095  for (int b = 0; b < num_card; ++b)
1096  for (int ip = 0; ip < num_quad; ++ip)
1097  for (int d = 0; d < num_dim; ++d)
1098  grad_basis_ref(b,ip,d) = dyn_grad_basis_ref(b,ip,d);
1099  }
1100  else if(elmtspace==PureBasis::HCURL && in_compute_derivatives && num_dim==2) {
1101  ArrayDynamic dyn_curl_basis_ref = af.buildArray<Scalar,BASIS,IP>("dyn_curl_basis_ref_scalar",num_card,num_quad);
1102 
1103  intrepid_basis->getValues(dyn_curl_basis_ref.get_view(),
1104  dyn_cub_points.get_view(),
1105  Intrepid2::OPERATOR_CURL);
1106 
1107  for (int b = 0; b < num_card; ++b)
1108  for (int ip = 0; ip < num_quad; ++ip)
1109  curl_basis_ref_scalar(b,ip) = dyn_curl_basis_ref(b,ip);
1110  }
1111  else if(elmtspace==PureBasis::HCURL && in_compute_derivatives && num_dim==3) {
1112  ArrayDynamic dyn_curl_basis_ref = af.buildArray<Scalar,BASIS,IP,Dim>("dyn_curl_basis_ref_vector",num_card,num_quad,num_dim);
1113 
1114  intrepid_basis->getValues(dyn_curl_basis_ref.get_view(),
1115  dyn_cub_points.get_view(),
1116  Intrepid2::OPERATOR_CURL);
1117 
1118  for (int b = 0; b < num_card; ++b)
1119  for (int ip = 0; ip < num_quad; ++ip)
1120  for (int d = 0; d < num_dim; ++d)
1121  curl_basis_ref_vector(b,ip,d) = dyn_curl_basis_ref(b,ip,d);
1122  }
1123  else if(elmtspace==PureBasis::HDIV && in_compute_derivatives) {
1124  ArrayDynamic dyn_div_basis_ref = af.buildArray<Scalar,BASIS,IP>("dyn_div_basis_ref_scalar",num_card,num_quad);
1125 
1126  intrepid_basis->getValues(dyn_div_basis_ref.get_view(),
1127  dyn_cub_points.get_view(),
1128  Intrepid2::OPERATOR_DIV);
1129 
1130  for (int b = 0; b < num_card; ++b)
1131  for (int ip = 0; ip < num_quad; ++ip)
1132  div_basis_ref(b,ip) = dyn_div_basis_ref(b,ip);
1133  }
1134 
1135 
1136  if(use_vertex_coordinates) {
1137  // Intrepid removes fad types from the coordinate scalar type. We
1138  // pull the actual field scalar type from the basis object to be
1139  // consistent.
1140  if (elmtspace != PureBasis::CONST) {
1141  using coordsScalarType = typename Intrepid2::Basis<PHX::Device::execution_space,Scalar,Scalar>::scalarType;
1142  auto dyn_basis_coordinates_ref = af.buildArray<coordsScalarType,BASIS,Dim>("basis_coordinates_ref",
1143  basis_coordinates_ref.extent(0),
1144  basis_coordinates_ref.extent(1));
1145  intrepid_basis->getDofCoords(dyn_basis_coordinates_ref.get_view());
1146 
1147  // fill in basis coordinates
1148  for (int i = 0; i < basis_coordinates_ref.extent_int(0); ++i)
1149  for (int j = 0; j < basis_coordinates_ref.extent_int(1); ++j)
1150  basis_coordinates_ref(i,j) = dyn_basis_coordinates_ref(i,j);
1151  }
1152  }
1153 
1154  references_evaluated = true;
1155 }
1156 
1157 // method for applying orientations
1158 template <typename Scalar>
1160 applyOrientations(const std::vector<Intrepid2::Orientation> & orientations,
1161  const int in_num_cells)
1162 {
1163  if (!intrepid_basis->requireOrientation())
1164  return;
1165 
1166  typedef Intrepid2::OrientationTools<PHX::Device> ots;
1167  const PureBasis::EElementSpace elmtspace = getElementSpace();
1168 
1169  // orientation (right now std vector) is created using push_back method.
1170  // thus, its size is the actual number of elements to be applied.
1171  // on the other hand, basis_layout num cell indicates workset size.
1172  // to get right size of cells, use minimum of them.
1173  const int num_cell_basis_layout = in_num_cells < 0 ? basis_layout->numCells() : in_num_cells;
1174  const int num_cell_orientation = orientations.size();
1175  const int num_cell = num_cell_basis_layout < num_cell_orientation ? num_cell_basis_layout : num_cell_orientation;
1176  const int num_dim = basis_layout->dimension();
1177  const Kokkos::pair<int,int> range_cell(0, num_cell);
1178 
1179  // Kokkos::DynRankView<Intrepid2::Orientation,PHX::Device>
1180  // drv_orts((Intrepid2::Orientation*)orientations.data(), num_cell);
1181  Kokkos::DynRankView<Intrepid2::Orientation,PHX::Device> drv_orts("drv_orts", num_cell);
1182  auto host_drv_orts = Kokkos::create_mirror_view(drv_orts);
1183  for (size_t i=0; i < drv_orts.size(); ++i)
1184  host_drv_orts(i) = orientations[i];
1185  Kokkos::deep_copy(drv_orts,host_drv_orts);
1186  typename PHX::Device().fence();
1187 
1191  if (elmtspace==PureBasis::HGRAD) {
1192  {
1193  {
1194  auto drv_basis_scalar = Kokkos::subview(basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1195  auto drv_basis_scalar_tmp = Kokkos::createDynRankView(basis_scalar.get_view(),
1196  "drv_basis_scalar_tmp",
1197  drv_basis_scalar.extent(0), // C
1198  drv_basis_scalar.extent(1), // F
1199  drv_basis_scalar.extent(2)); // P
1200  Kokkos::deep_copy(drv_basis_scalar_tmp, drv_basis_scalar);
1201  ots::modifyBasisByOrientation(drv_basis_scalar,
1202  drv_basis_scalar_tmp,
1203  drv_orts,
1204  intrepid_basis.getRawPtr());
1205  }
1206  if(build_weighted) {
1207  auto drv_basis_scalar = Kokkos::subview(weighted_basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1208  auto drv_basis_scalar_tmp = Kokkos::createDynRankView(weighted_basis_scalar.get_view(),
1209  "drv_basis_scalar_tmp",
1210  drv_basis_scalar.extent(0), // C
1211  drv_basis_scalar.extent(1), // F
1212  drv_basis_scalar.extent(2)); // P
1213  Kokkos::deep_copy(drv_basis_scalar_tmp, drv_basis_scalar);
1214  ots::modifyBasisByOrientation(drv_basis_scalar,
1215  drv_basis_scalar_tmp,
1216  drv_orts,
1217  intrepid_basis.getRawPtr());
1218  }
1219 
1220  }
1221 
1222  if (compute_derivatives) {
1223  {
1224  auto drv_grad_basis = Kokkos::subview(grad_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1225  auto drv_grad_basis_tmp = Kokkos::createDynRankView(grad_basis.get_view(),
1226  "drv_grad_basis_tmp",
1227  drv_grad_basis.extent(0), // C
1228  drv_grad_basis.extent(1), // F
1229  drv_grad_basis.extent(2), // P
1230  drv_grad_basis.extent(3)); // D
1231  Kokkos::deep_copy(drv_grad_basis_tmp, drv_grad_basis);
1232  ots::modifyBasisByOrientation(drv_grad_basis,
1233  drv_grad_basis_tmp,
1234  drv_orts,
1235  intrepid_basis.getRawPtr());
1236  }
1237  if(build_weighted) {
1238  auto drv_grad_basis = Kokkos::subview(weighted_grad_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1239  auto drv_grad_basis_tmp = Kokkos::createDynRankView(weighted_grad_basis.get_view(),
1240  "drv_grad_basis_tmp",
1241  drv_grad_basis.extent(0), // C
1242  drv_grad_basis.extent(1), // F
1243  drv_grad_basis.extent(2), // P
1244  drv_grad_basis.extent(3)); // D
1245  Kokkos::deep_copy(drv_grad_basis_tmp, drv_grad_basis);
1246  ots::modifyBasisByOrientation(drv_grad_basis,
1247  drv_grad_basis_tmp,
1248  drv_orts,
1249  intrepid_basis.getRawPtr());
1250  }
1251  }
1252  }
1253 
1257  else if (elmtspace==PureBasis::HCURL && num_dim==2) {
1258  {
1259  {
1260  auto drv_basis_vector = Kokkos::subview(basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1261  auto drv_basis_vector_tmp = Kokkos::createDynRankView(basis_vector.get_view(),
1262  "drv_basis_vector_tmp",
1263  drv_basis_vector.extent(0), // C
1264  drv_basis_vector.extent(1), // F
1265  drv_basis_vector.extent(2), // P
1266  drv_basis_vector.extent(3)); // D
1267  Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1268  ots::modifyBasisByOrientation(drv_basis_vector,
1269  drv_basis_vector_tmp,
1270  drv_orts,
1271  intrepid_basis.getRawPtr());
1272  }
1273  if(build_weighted) {
1274  auto drv_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1275  auto drv_basis_vector_tmp = Kokkos::createDynRankView(weighted_basis_vector.get_view(),
1276  "drv_basis_vector_tmp",
1277  drv_basis_vector.extent(0), // C
1278  drv_basis_vector.extent(1), // F
1279  drv_basis_vector.extent(2), // P
1280  drv_basis_vector.extent(3)); // D
1281  Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1282  ots::modifyBasisByOrientation(drv_basis_vector,
1283  drv_basis_vector_tmp,
1284  drv_orts,
1285  intrepid_basis.getRawPtr());
1286  }
1287  }
1288 
1289  if (compute_derivatives) {
1290  {
1291  auto drv_curl_basis_scalar = Kokkos::subview(curl_basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1292  auto drv_curl_basis_scalar_tmp = Kokkos::createDynRankView(curl_basis_scalar.get_view(),
1293  "drv_curl_basis_scalar_tmp",
1294  drv_curl_basis_scalar.extent(0), // C
1295  drv_curl_basis_scalar.extent(1), // F
1296  drv_curl_basis_scalar.extent(2)); // P
1297  Kokkos::deep_copy(drv_curl_basis_scalar_tmp, drv_curl_basis_scalar);
1298  ots::modifyBasisByOrientation(drv_curl_basis_scalar,
1299  drv_curl_basis_scalar_tmp,
1300  drv_orts,
1301  intrepid_basis.getRawPtr());
1302  }
1303 
1304  if(build_weighted) {
1305  auto drv_curl_basis_scalar = Kokkos::subview(weighted_curl_basis_scalar.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1306  auto drv_curl_basis_scalar_tmp = Kokkos::createDynRankView(weighted_curl_basis_scalar.get_view(),
1307  "drv_curl_basis_scalar_tmp",
1308  drv_curl_basis_scalar.extent(0), // C
1309  drv_curl_basis_scalar.extent(1), // F
1310  drv_curl_basis_scalar.extent(2)); // P
1311  Kokkos::deep_copy(drv_curl_basis_scalar_tmp, drv_curl_basis_scalar);
1312  ots::modifyBasisByOrientation(drv_curl_basis_scalar,
1313  drv_curl_basis_scalar_tmp,
1314  drv_orts,
1315  intrepid_basis.getRawPtr());
1316  }
1317  }
1318  }
1319 
1323  else if (elmtspace==PureBasis::HCURL && num_dim==3) {
1324  {
1325  {
1326  auto drv_basis_vector = Kokkos::subview(basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1327  auto drv_basis_vector_tmp = Kokkos::createDynRankView(basis_vector.get_view(),
1328  "drv_basis_vector_tmp",
1329  drv_basis_vector.extent(0), // C
1330  drv_basis_vector.extent(1), // F
1331  drv_basis_vector.extent(2), // P
1332  drv_basis_vector.extent(3)); // D
1333  Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1334  ots::modifyBasisByOrientation(drv_basis_vector,
1335  drv_basis_vector_tmp,
1336  drv_orts,
1337  intrepid_basis.getRawPtr());
1338  }
1339  if(build_weighted) {
1340  auto drv_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1341  auto drv_basis_vector_tmp = Kokkos::createDynRankView(weighted_basis_vector.get_view(),
1342  "drv_basis_vector_tmp",
1343  drv_basis_vector.extent(0), // C
1344  drv_basis_vector.extent(1), // F
1345  drv_basis_vector.extent(2), // P
1346  drv_basis_vector.extent(3)); // D
1347  Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1348  ots::modifyBasisByOrientation(drv_basis_vector,
1349  drv_basis_vector_tmp,
1350  drv_orts,
1351  intrepid_basis.getRawPtr());
1352  }
1353  }
1354 
1355  if (compute_derivatives) {
1356  {
1357  auto drv_curl_basis_vector = Kokkos::subview(curl_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1358  auto drv_curl_basis_vector_tmp = Kokkos::createDynRankView(curl_basis_vector.get_view(),
1359  "drv_curl_basis_vector_tmp",
1360  drv_curl_basis_vector.extent(0), // C
1361  drv_curl_basis_vector.extent(1), // F
1362  drv_curl_basis_vector.extent(2), // P
1363  drv_curl_basis_vector.extent(3)); // D
1364  Kokkos::deep_copy(drv_curl_basis_vector_tmp, drv_curl_basis_vector);
1365  ots::modifyBasisByOrientation(drv_curl_basis_vector,
1366  drv_curl_basis_vector_tmp,
1367  drv_orts,
1368  intrepid_basis.getRawPtr());
1369  }
1370  if(build_weighted) {
1371  auto drv_curl_basis_vector = Kokkos::subview(weighted_curl_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1372  auto drv_curl_basis_vector_tmp = Kokkos::createDynRankView(weighted_curl_basis_vector.get_view(),
1373  "drv_curl_basis_vector_tmp",
1374  drv_curl_basis_vector.extent(0), // C
1375  drv_curl_basis_vector.extent(1), // F
1376  drv_curl_basis_vector.extent(2), // P
1377  drv_curl_basis_vector.extent(3)); // D
1378  Kokkos::deep_copy(drv_curl_basis_vector_tmp, drv_curl_basis_vector);
1379  ots::modifyBasisByOrientation(drv_curl_basis_vector,
1380  drv_curl_basis_vector_tmp,
1381  drv_orts,
1382  intrepid_basis.getRawPtr());
1383  }
1384  }
1385  }
1389  else if (elmtspace==PureBasis::HDIV) {
1390  {
1391  {
1392  auto drv_basis_vector = Kokkos::subview(basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1393  auto drv_basis_vector_tmp = Kokkos::createDynRankView(basis_vector.get_view(),
1394  "drv_basis_vector_tmp",
1395  drv_basis_vector.extent(0), // C
1396  drv_basis_vector.extent(1), // F
1397  drv_basis_vector.extent(2), // P
1398  drv_basis_vector.extent(3)); // D
1399  Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1400  ots::modifyBasisByOrientation(drv_basis_vector,
1401  drv_basis_vector_tmp,
1402  drv_orts,
1403  intrepid_basis.getRawPtr());
1404  }
1405  if(build_weighted) {
1406  auto drv_basis_vector = Kokkos::subview(weighted_basis_vector.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL(), Kokkos::ALL());
1407  auto drv_basis_vector_tmp = Kokkos::createDynRankView(weighted_basis_vector.get_view(),
1408  "drv_basis_vector_tmp",
1409  drv_basis_vector.extent(0), // C
1410  drv_basis_vector.extent(1), // F
1411  drv_basis_vector.extent(2), // P
1412  drv_basis_vector.extent(3)); // D
1413  Kokkos::deep_copy(drv_basis_vector_tmp, drv_basis_vector);
1414  ots::modifyBasisByOrientation(drv_basis_vector,
1415  drv_basis_vector_tmp,
1416  drv_orts,
1417  intrepid_basis.getRawPtr());
1418  }
1419  }
1420  if (compute_derivatives) {
1421  {
1422  auto drv_div_basis = Kokkos::subview(div_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1423  auto drv_div_basis_tmp = Kokkos::createDynRankView(div_basis.get_view(),
1424  "drv_div_basis_tmp",
1425  drv_div_basis.extent(0), // C
1426  drv_div_basis.extent(1), // F
1427  drv_div_basis.extent(2)); // P
1428  Kokkos::deep_copy(drv_div_basis_tmp, drv_div_basis);
1429  ots::modifyBasisByOrientation(drv_div_basis,
1430  drv_div_basis_tmp,
1431  drv_orts,
1432  intrepid_basis.getRawPtr());
1433  }
1434  if(build_weighted) {
1435  auto drv_div_basis = Kokkos::subview(weighted_div_basis.get_view(), range_cell, Kokkos::ALL(), Kokkos::ALL());
1436  auto drv_div_basis_tmp = Kokkos::createDynRankView(weighted_div_basis.get_view(),
1437  "drv_div_basis_tmp",
1438  drv_div_basis.extent(0), // C
1439  drv_div_basis.extent(1), // F
1440  drv_div_basis.extent(2)); // P
1441  Kokkos::deep_copy(drv_div_basis_tmp, drv_div_basis);
1442  ots::modifyBasisByOrientation(drv_div_basis,
1443  drv_div_basis_tmp,
1444  drv_orts,
1445  intrepid_basis.getRawPtr());
1446  }
1447  }
1448  }
1449 }
1450 
1451 // method for applying orientations
1452 template <typename Scalar>
1454 applyOrientations(const PHX::MDField<const Scalar,Cell,BASIS> & orientations)
1455 {
1456 
1457  TEUCHOS_TEST_FOR_EXCEPT_MSG(true,"panzer::BasisValues2::applyOrientations : this should not be called.");
1458 
1459  int num_cell = orientations.extent(0);
1460  int num_basis = orientations.extent(1);
1461  int num_dim = basis_layout->dimension();
1462  int num_ip = basis_layout->numPoints();
1463  PureBasis::EElementSpace elmtspace = getElementSpace();
1464 
1465  if(elmtspace==PureBasis::HCURL && num_dim==2) {
1466 
1467  // setup the orientations for the trial space
1468  // Intrepid2::FunctionSpaceTools::applyFieldSigns<Scalar>(basis_vector,orientations);
1469 
1470  for (int c=0; c<num_cell; c++)
1471  for (int b=0; b<num_basis; b++)
1472  for (int p=0; p<num_ip; p++)
1473  for (int d=0; d<num_dim; d++)
1474  basis_vector(c, b, p, d) *= orientations(c, b);
1475 
1476  if(compute_derivatives) {
1477  // Intrepid2::FunctionSpaceTools::applyFieldSigns<Scalar>(curl_basis_scalar,orientations);
1478  for (int c=0; c<num_cell; c++)
1479  for (int b=0; b<num_basis; b++)
1480  for (int p=0; p<num_ip; p++)
1481  curl_basis_scalar(c, b, p) *= orientations(c, b);
1482  }
1483 
1484  // setup the orientations for the test space
1485  if(build_weighted) {
1486  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_basis_vector.get_view(),orientations.get_view());
1487 
1488  if(compute_derivatives)
1489  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_curl_basis_scalar.get_view(),orientations.get_view());
1490  }
1491  }
1492  else if(elmtspace==PureBasis::HCURL && num_dim==3) {
1493 
1494  // setup the orientations for the trial space
1495  // Intrepid2::FunctionSpaceTools::applyFieldSigns<Scalar>(basis_vector,orientations);
1496 
1497  for (int c=0; c<num_cell; c++)
1498  for (int b=0; b<num_basis; b++)
1499  for (int p=0; p<num_ip; p++)
1500  for (int d=0; d<num_dim; d++)
1501  basis_vector(c, b, p, d) *= orientations(c, b);
1502 
1503  if(compute_derivatives) {
1504  // Intrepid2::FunctionSpaceTools::applyFieldSigns<Scalar>(curl_basis_vector,orientations);
1505  for (int c=0; c<num_cell; c++)
1506  for (int b=0; b<num_basis; b++)
1507  for (int p=0; p<num_ip; p++)
1508  for (int d=0; d<num_dim; d++)
1509  curl_basis_vector(c, b, p,d) *= orientations(c, b);
1510  }
1511 
1512  // setup the orientations for the test space
1513  if(build_weighted) {
1514  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_basis_vector.get_view(),orientations.get_view());
1515 
1516  if(compute_derivatives)
1517  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_curl_basis_vector.get_view(),orientations.get_view());
1518  }
1519  }
1520  else if(elmtspace==PureBasis::HDIV) {
1521  // setup the orientations for the trial space
1522  // Intrepid2::FunctionSpaceTools::applyFieldSigns<Scalar>(basis_vector,orientations);
1523 
1524  for (int c=0; c<num_cell; c++)
1525  for (int b=0; b<num_basis; b++)
1526  for (int p=0; p<num_ip; p++)
1527  for (int d=0; d<num_dim; d++)
1528  basis_vector(c, b, p, d) *= orientations(c, b);
1529 
1530  if(compute_derivatives) {
1531  // Intrepid2::FunctionSpaceTools::applyFieldSigns<Scalar>(div_basis,orientations);
1532 
1533  for (int c=0; c<num_cell; c++)
1534  for (int b=0; b<num_basis; b++)
1535  for (int p=0; p<num_ip; p++)
1536  div_basis(c, b, p) *= orientations(c, b);
1537  }
1538 
1539  // setup the orientations for the test space
1540  if(build_weighted) {
1541  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_basis_vector.get_view(),orientations.get_view());
1542 
1543  if(compute_derivatives)
1544  Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::applyFieldSigns(weighted_div_basis.get_view(),orientations.get_view());
1545  }
1546  }
1547 }
1548 
1549 template <typename Scalar>
1551 { return basis_layout->getBasis()->getElementSpace(); }
1552 
1553 template <typename Scalar>
1556  bool computeDerivatives)
1557 {
1558  MDFieldArrayFactory af(prefix,alloc_arrays);
1559 
1560  compute_derivatives = computeDerivatives;
1561  basis_layout = layout;
1562  Teuchos::RCP<const panzer::PureBasis> basisDesc = layout->getBasis();
1563 
1564  // for convience pull out basis and quadrature information
1565  int num_quad = layout->numPoints();
1566  int dim = basisDesc->dimension();
1567  int card = basisDesc->cardinality();
1568  int numcells = basisDesc->numCells();
1569  panzer::PureBasis::EElementSpace elmtspace = basisDesc->getElementSpace();
1571 
1572  intrepid_basis = basisDesc->getIntrepid2Basis<PHX::Device::execution_space,Scalar,Scalar>();
1573 
1574  // allocate field containers
1575  // field sizes defined by http://trilinos.sandia.gov/packages/docs/dev/packages/intrepid/doc/html/basis_page.html#basis_md_array_sec
1576 
1577  // compute basis fields
1578  if(elmtspace==panzer::PureBasis::HGRAD) {
1579  // HGRAD is a nodal field
1580 
1581  // build values
1583  basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("basis_ref",card,num_quad); // F, P
1584  basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("basis",numcells,card,num_quad);
1585 
1586  if(build_weighted)
1587  weighted_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("weighted_basis",numcells,card,num_quad);
1588 
1589  // build gradients
1591 
1592  if(compute_derivatives) {
1593  grad_basis_ref = af.buildStaticArray<Scalar,BASIS,IP,Dim>("grad_basis_ref",card,num_quad,dim); // F, P, D
1594  grad_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("grad_basis",numcells,card,num_quad,dim);
1595 
1596  if(build_weighted)
1597  weighted_grad_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("weighted_grad_basis",numcells,card,num_quad,dim);
1598  }
1599 
1600  // build curl
1602 
1603  // nothing - HGRAD does not support CURL operation
1604  }
1605  else if(elmtspace==panzer::PureBasis::HCURL) {
1606  // HCURL is a vector field
1607 
1608  // build values
1610 
1611  basis_ref_vector = af.buildStaticArray<Scalar,BASIS,IP,Dim>("basis_ref",card,num_quad,dim); // F, P, D
1612  basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("basis",numcells,card,num_quad,dim);
1613 
1614  if(build_weighted)
1615  weighted_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("weighted_basis",numcells,card,num_quad,dim);
1616 
1617  // build gradients
1619 
1620  // nothing - HCURL does not support GRAD operation
1621 
1622  // build curl
1624 
1625  if(compute_derivatives) {
1626  if(dim==2) {
1627  // curl of HCURL basis is not dimension dependent
1628  curl_basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("curl_basis_ref",card,num_quad); // F, P
1629  curl_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("curl_basis",numcells,card,num_quad);
1630 
1631  if(build_weighted)
1632  weighted_curl_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("weighted_curl_basis",numcells,card,num_quad);
1633  }
1634  else if(dim==3){
1635  curl_basis_ref_vector = af.buildStaticArray<Scalar,BASIS,IP,Dim>("curl_basis_ref",card,num_quad,dim); // F, P, D
1636  curl_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("curl_basis",numcells,card,num_quad,dim);
1637 
1638  if(build_weighted)
1639  weighted_curl_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("weighted_curl_basis",numcells,card,num_quad,dim);
1640  }
1641  else { TEUCHOS_ASSERT(false); } // what do I do with 1D?
1642  }
1643  }
1644  else if(elmtspace==panzer::PureBasis::HDIV) {
1645  // HDIV is a vector field
1646 
1647  // build values
1649 
1650  basis_ref_vector = af.buildStaticArray<Scalar,BASIS,IP,Dim>("basis_ref",card,num_quad,dim); // F, P, D
1651  basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("basis",numcells,card,num_quad,dim);
1652 
1653  if(build_weighted)
1654  weighted_basis_vector = af.buildStaticArray<Scalar,Cell,BASIS,IP,Dim>("weighted_basis",numcells,card,num_quad,dim);
1655 
1656  // build gradients
1658 
1659  // nothing - HCURL does not support GRAD operation
1660 
1661  // build curl
1663 
1664  // nothing - HDIV does not support CURL operation
1665 
1666  // build div
1668 
1669  if(compute_derivatives) {
1670  div_basis_ref = af.buildStaticArray<Scalar,BASIS,IP>("div_basis_ref",card,num_quad); // F, P
1671  div_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP>("div_basis",numcells,card,num_quad);
1672 
1673  if(build_weighted)
1674  weighted_div_basis = af.buildStaticArray<Scalar,Cell,BASIS,IP>("weighted_div_basis",numcells,card,num_quad);
1675  }
1676  }
1677  else if(elmtspace==panzer::PureBasis::CONST || elmtspace==panzer::PureBasis::HVOL) {
1678  // CONST is a nodal field
1679 
1680  // build values
1682  basis_ref_scalar = af.buildStaticArray<Scalar,BASIS,IP>("basis_ref",card,num_quad); // F, P
1683  basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("basis",numcells,card,num_quad);
1684 
1685  if(build_weighted)
1686  weighted_basis_scalar = af.buildStaticArray<Scalar,Cell,BASIS,IP>("weighted_basis",numcells,card,num_quad);
1687 
1688  // build gradients
1690 
1691  // nothing - CONST does not support GRAD operation
1692 
1693  // build curl
1695 
1696  // nothing - CONST does not support CURL operation
1697 
1698  // build div
1700 
1701  // nothing - CONST does not support DIV operation
1702  }
1703  else { TEUCHOS_ASSERT(false); }
1704 
1705  basis_coordinates_ref = af.buildStaticArray<Scalar,BASIS,Dim>("basis_coordinates_ref",card,dim);
1706  basis_coordinates = af.buildStaticArray<Scalar,Cell,BASIS,Dim>("basis_coordinates",numcells,card,dim);
1707 }
1708 
1709 } // namespace panzer
void evaluateValues_HCurl(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack...dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
void evaluateValues_HDiv(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
void evaluateValuesCV(const PHX::MDField< Scalar, Cell, IP, Dim > &cell_cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv)
void evaluateReferenceValues(const PHX::MDField< Scalar, IP, Dim > &cub_points, bool compute_derivatives, bool use_vertex_coordinates)
EElementSpace getElementSpace() const
void evaluateValues_HGrad(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
int cardinality() const
Returns the number of basis coefficients.
PHX::MDField< Scalar, T0 > buildStaticArray(const std::string &str, int d0) const
Teuchos::RCP< const PureBasis > getBasis() const
int dimension() const
Returns the dimension of the basis from the topology.
void applyOrientations(const PHX::MDField< const Scalar, Cell, BASIS > &orientations)
Method to apply orientations to a basis values container.
Teuchos::RCP< const shards::CellTopology > getCellTopology() const
#define TEUCHOS_TEST_FOR_EXCEPT_MSG(throw_exception_test, msg)
void evaluateBasisCoordinates(const PHX::MDField< Scalar, Cell, NODE, Dim > &vertex_coordinates, const int in_num_cells=-1)
ArrayTraits< Scalar, PHX::MDField< Scalar > >::size_type size_type
PureBasis::EElementSpace getElementSpace() const
Teuchos::RCP< Intrepid2::Basis< PHX::Device::execution_space, double, double > > getIntrepid2Basis() const
void evaluateValues_Const(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
void setupArrays(const Teuchos::RCP< const panzer::BasisIRLayout > &basis, bool computeDerivatives=true)
Sizes/allocates memory for arrays.
int numCells() const
Returns the number of cells in the data layouts.
void evaluateValues_HVol(const PHX::MDField< Scalar, Cell, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const PHX::MDField< Scalar, Cell, IP > &weighted_measure, const int in_num_cells)
Description and data layouts associated with a particular basis.
#define TEUCHOS_ASSERT(assertion_test)
void evaluateValues(const PHX::MDField< Scalar, IP, Dim > &cub_points, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac, const PHX::MDField< Scalar, Cell, IP > &jac_det, const PHX::MDField< Scalar, Cell, IP, Dim, Dim > &jac_inv, const int in_num_cells=-1)
PHX::MDField< Scalar > buildArray(const std::string &str, int d0) const
PHX::MDField< Scalar > ArrayDynamic