Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_DOFCurl_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 #ifndef PANZER_DOF_CURL_IMPL_HPP
44 #define PANZER_DOF_CURL_IMPL_HPP
45 
47 #include "Panzer_BasisIRLayout.hpp"
49 #include "Intrepid2_FunctionSpaceTools.hpp"
50 #include "Phalanx_KokkosDeviceTypes.hpp"
51 
52 namespace panzer {
53 
54 namespace {
55 
56 //**********************************************************************
57 template <typename ScalarT,typename Array,int spaceDim>
58 class EvaluateCurlWithSens_Vector {
59  PHX::MDField<const ScalarT,Cell,Point> dof_value;
60  PHX::MDField<ScalarT,Cell,Point,Dim> dof_curl;
61  Array curl_basis;
62 
63  int numFields;
64  int numPoints;
65 
66 public:
67  typedef typename PHX::Device execution_space;
68 
69  EvaluateCurlWithSens_Vector(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
70  PHX::MDField<ScalarT,Cell,Point,Dim> in_dof_curl,
71  Array in_curl_basis)
72  : dof_value(in_dof_value), dof_curl(in_dof_curl), curl_basis(in_curl_basis)
73  {
74  numFields = curl_basis.extent(1);
75  numPoints = curl_basis.extent(2);
76  }
77  KOKKOS_INLINE_FUNCTION
78  void operator()(const unsigned int cell) const
79  {
80  for (int pt=0; pt<numPoints; pt++) {
81  for (int d=0; d<spaceDim; d++) {
82  // first initialize to the right thing (prevents over writing with 0)
83  // then loop over one less basis function
84  dof_curl(cell,pt,d) = dof_value(cell, 0) * curl_basis(cell, 0, pt, d);
85  for (int bf=1; bf<numFields; bf++)
86  dof_curl(cell,pt,d) += dof_value(cell, bf) * curl_basis(cell, bf, pt, d);
87  }
88  }
89  }
90 };
91 
92 template <typename ScalarT,typename ArrayT>
93 void evaluateCurl_withSens_vector(int numCells,
94  PHX::MDField<ScalarT,Cell,Point,Dim> & dof_curl,
95  PHX::MDField<const ScalarT,Cell,Point> & dof_value,
96  const ArrayT & curl_basis)
97 {
98  if(numCells>0) {
99  // evaluate at quadrature points
100  int numFields = curl_basis.extent(1);
101  int numPoints = curl_basis.extent(2);
102  int spaceDim = curl_basis.extent(3);
103 
104  for (int cell=0; cell<numCells; cell++) {
105  for (int pt=0; pt<numPoints; pt++) {
106  for (int d=0; d<spaceDim; d++) {
107  // first initialize to the right thing (prevents over writing with 0)
108  // then loop over one less basis function
109  dof_curl(cell,pt,d) = dof_value(cell, 0) * curl_basis(cell, 0, pt, d);
110  for (int bf=1; bf<numFields; bf++)
111  dof_curl(cell,pt,d) += dof_value(cell, bf) * curl_basis(cell, bf, pt, d);
112  }
113  }
114  }
115  }
116 }
117 
118 //**********************************************************************
119 template <typename ScalarT,typename Array>
120 class EvaluateCurlWithSens_Scalar {
121  PHX::MDField<const ScalarT,Cell,Point> dof_value;
122  PHX::MDField<ScalarT,Cell,Point> dof_curl;
123  Array curl_basis;
124 
125  int numFields;
126  int numPoints;
127 
128 public:
129  typedef typename PHX::Device execution_space;
130 
131  EvaluateCurlWithSens_Scalar(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
132  PHX::MDField<ScalarT,Cell,Point> in_dof_curl,
133  Array in_curl_basis)
134  : dof_value(in_dof_value), dof_curl(in_dof_curl), curl_basis(in_curl_basis)
135  {
136  numFields = curl_basis.extent(1);
137  numPoints = curl_basis.extent(2);
138  }
139  KOKKOS_INLINE_FUNCTION
140  void operator()(const unsigned int cell) const
141  {
142  for (int pt=0; pt<numPoints; pt++) {
143  // first initialize to the right thing (prevents over writing with 0)
144  // then loop over one less basis function
145  dof_curl(cell,pt) = dof_value(cell, 0) * curl_basis(cell, 0, pt);
146  for (int bf=1; bf<numFields; bf++)
147  dof_curl(cell,pt) += dof_value(cell, bf) * curl_basis(cell, bf, pt);
148  }
149  }
150 };
151 
152 template <typename ScalarT,typename ArrayT>
153 void evaluateCurl_withSens_scalar(int numCells,
154  PHX::MDField<ScalarT,Cell,Point> & dof_curl,
155  PHX::MDField<const ScalarT,Cell,Point> & dof_value,
156  const ArrayT & curl_basis)
157 {
158  if(numCells>0) {
159  // evaluate at quadrature points
160  int numFields = curl_basis.extent(1);
161  int numPoints = curl_basis.extent(2);
162 
163  for (int cell=0; cell<numCells; cell++) {
164  for (int pt=0; pt<numPoints; pt++) {
165  // first initialize to the right thing (prevents over writing with 0)
166  // then loop over one less basis function
167  dof_curl(cell,pt) = dof_value(cell, 0) * curl_basis(cell, 0, pt);
168  for (int bf=1; bf<numFields; bf++)
169  dof_curl(cell,pt) += dof_value(cell, bf) * curl_basis(cell, bf, pt);
170  }
171  }
172  }
173 }
174 
175 //**********************************************************************
176 template <typename ScalarT,typename Array,int spaceDim>
177 class EvaluateCurlFastSens_Vector {
178  PHX::MDField<const ScalarT,Cell,Point> dof_value;
179  PHX::MDField<ScalarT,Cell,Point,Dim> dof_curl;
180  Kokkos::View<const int*,PHX::Device> offsets;
181  Array curl_basis;
182 
183  int numFields;
184  int numPoints;
185 
186 public:
187  typedef typename PHX::Device execution_space;
188 
189  EvaluateCurlFastSens_Vector(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
190  PHX::MDField<ScalarT,Cell,Point,Dim> in_dof_curl,
191  Kokkos::View<const int*,PHX::Device> in_offsets,
192  Array in_curl_basis)
193  : dof_value(in_dof_value), dof_curl(in_dof_curl), offsets(in_offsets), curl_basis(in_curl_basis)
194  {
195  numFields = curl_basis.extent(1);
196  numPoints = curl_basis.extent(2);
197  }
198  KOKKOS_INLINE_FUNCTION
199  void operator()(const unsigned int cell) const
200  {
201  for (int pt=0; pt<numPoints; pt++) {
202  for (int d=0; d<spaceDim; d++) {
203  // first initialize to the right thing (prevents over writing with 0)
204  // then loop over one less basis function
205  dof_curl(cell,pt,d) = dof_value(cell, 0).val() * curl_basis(cell, 0, pt, d);
206  dof_curl(cell,pt,d).fastAccessDx(offsets(0)) = dof_value(cell, 0).fastAccessDx(offsets(0)) * curl_basis(cell, 0, pt, d);
207  for (int bf=1; bf<numFields; bf++) {
208  dof_curl(cell,pt,d).val() += dof_value(cell, bf).val() * curl_basis(cell, bf, pt, d);
209  dof_curl(cell,pt,d).fastAccessDx(offsets(bf)) += dof_value(cell, bf).fastAccessDx(offsets(bf)) * curl_basis(cell, bf, pt, d);
210  }
211  }
212  }
213  }
214 };
215 template <typename ScalarT,typename ArrayT>
216 void evaluateCurl_fastSens_vector(int numCells,
217  PHX::MDField<ScalarT,Cell,Point,Dim> & dof_curl,
218  PHX::MDField<const ScalarT,Cell,Point> & dof_value,
219  const std::vector<int> & offsets,
220  const ArrayT & curl_basis)
221 {
222  if(numCells>0) {
223  int numFields = curl_basis.extent(1);
224  int numPoints = curl_basis.extent(2);
225  int spaceDim = curl_basis.extent(3);
226 
227  for (int cell=0; cell<numCells; cell++) {
228  for (int pt=0; pt<numPoints; pt++) {
229  for (int d=0; d<spaceDim; d++) {
230  // first initialize to the right thing (prevents over writing with 0)
231  // then loop over one less basis function
232  dof_curl(cell,pt,d) = ScalarT(numFields, dof_value(cell, 0).val() * curl_basis(cell, 0, pt, d));
233  dof_curl(cell,pt,d).fastAccessDx(offsets[0]) = dof_value(cell, 0).fastAccessDx(offsets[0]) * curl_basis(cell, 0, pt, d);
234  for (int bf=1; bf<numFields; bf++) {
235  dof_curl(cell,pt,d).val() += dof_value(cell, bf).val() * curl_basis(cell, bf, pt, d);
236  dof_curl(cell,pt,d).fastAccessDx(offsets[bf]) += dof_value(cell, bf).fastAccessDx(offsets[bf]) * curl_basis(cell, bf, pt, d);
237  }
238  }
239  }
240  }
241  }
242 }
243 
244 //**********************************************************************
245 template <typename ScalarT,typename Array>
246 class EvaluateCurlFastSens_Scalar {
247  PHX::MDField<const ScalarT,Cell,Point> dof_value;
248  PHX::MDField<ScalarT,Cell,Point> dof_curl;
249  Kokkos::View<const int*,PHX::Device> offsets;
250  Array curl_basis;
251 
252  int numFields;
253  int numPoints;
254 
255 public:
256  typedef typename PHX::Device execution_space;
257 
258  EvaluateCurlFastSens_Scalar(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
259  PHX::MDField<ScalarT,Cell,Point> in_dof_curl,
260  Kokkos::View<const int*,PHX::Device> in_offsets,
261  Array in_curl_basis)
262  : dof_value(in_dof_value), dof_curl(in_dof_curl), offsets(in_offsets), curl_basis(in_curl_basis)
263  {
264  numFields = curl_basis.extent(1);
265  numPoints = curl_basis.extent(2);
266  }
267  KOKKOS_INLINE_FUNCTION
268  void operator()(const unsigned int cell) const
269  {
270  for (int pt=0; pt<numPoints; pt++) {
271  // first initialize to the right thing (prevents over writing with 0)
272  // then loop over one less basis function
273  dof_curl(cell,pt) = dof_value(cell, 0).val() * curl_basis(cell, 0, pt);
274  dof_curl(cell,pt).fastAccessDx(offsets(0)) = dof_value(cell, 0).fastAccessDx(offsets(0)) * curl_basis(cell, 0, pt);
275  for (int bf=1; bf<numFields; bf++) {
276  dof_curl(cell,pt).val() += dof_value(cell, bf).val() * curl_basis(cell, bf, pt);
277  dof_curl(cell,pt).fastAccessDx(offsets(bf)) += dof_value(cell, bf).fastAccessDx(offsets(bf)) * curl_basis(cell, bf, pt);
278  }
279  }
280  }
281 };
282 template <typename ScalarT,typename ArrayT>
283 void evaluateCurl_fastSens_scalar(int numCells,
284  PHX::MDField<ScalarT,Cell,Point> & dof_curl,
285  PHX::MDField<const ScalarT,Cell,Point> & dof_value,
286  const std::vector<int> & offsets,
287  const ArrayT & curl_basis)
288 {
289  if(numCells>0) {
290  int numFields = curl_basis.extent(1);
291  int numPoints = curl_basis.extent(2);
292 
293  for (int cell=0; cell<numCells; cell++) {
294  for (int pt=0; pt<numPoints; pt++) {
295  // first initialize to the right thing (prevents over writing with 0)
296  // then loop over one less basis function
297  dof_curl(cell,pt) = ScalarT(numFields, dof_value(cell, 0).val() * curl_basis(cell, 0, pt));
298  dof_curl(cell,pt).fastAccessDx(offsets[0]) = dof_value(cell, 0).fastAccessDx(offsets[0]) * curl_basis(cell, 0, pt);
299  for (int bf=1; bf<numFields; bf++) {
300  dof_curl(cell,pt).val() += dof_value(cell, bf).val() * curl_basis(cell, bf, pt);
301  dof_curl(cell,pt).fastAccessDx(offsets[bf]) += dof_value(cell, bf).fastAccessDx(offsets[bf]) * curl_basis(cell, bf, pt);
302  }
303  }
304  }
305  }
306 }
307 
308 //**********************************************************************
309 
310 }
311 
312 //**********************************************************************
313 // MOST EVALUATION TYPES
314 //**********************************************************************
315 
316 //**********************************************************************
317 template<typename EvalT, typename TRAITS>
320  use_descriptors_(false),
321  dof_value( p.get<std::string>("Name"),
322  p.get< Teuchos::RCP<panzer::BasisIRLayout> >("Basis")->functional),
323  basis_name(p.get< Teuchos::RCP<panzer::BasisIRLayout> >("Basis")->name())
324 {
326  = p.get< Teuchos::RCP<BasisIRLayout> >("Basis")->getBasis();
327 
328  // Verify that this basis supports the curl operation
329  TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsCurl(),std::logic_error,
330  "DOFCurl: Basis of type \"" << basis->name() << "\" does not support CURL");
331  TEUCHOS_TEST_FOR_EXCEPTION(!basis->requiresOrientations(),std::logic_error,
332  "DOFCurl: Basis of type \"" << basis->name() << "\" in DOF Curl should require orientations. So we are throwing.");
333 
334  // build dof_curl
335  basis_dimension = basis->dimension();
336  if(basis_dimension==2) {
337  dof_curl_scalar = PHX::MDField<ScalarT,Cell,Point>(p.get<std::string>("Curl Name"),
338  p.get< Teuchos::RCP<panzer::IntegrationRule> >("IR")->dl_scalar );
339  this->addEvaluatedField(dof_curl_scalar);
340  }
341  else if(basis_dimension==3) {
342  dof_curl_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(p.get<std::string>("Curl Name"),
343  p.get< Teuchos::RCP<panzer::IntegrationRule> >("IR")->dl_vector );
344  this->addEvaluatedField(dof_curl_vector);
345  }
346  else
347  { TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,"DOFCurl only works for 2D and 3D basis functions"); }
348 
349  // add to evaluation graph
350  this->addDependentField(dof_value);
351 
352  std::string n = "DOFCurl: " + (basis_dimension==2 ? dof_curl_scalar.fieldTag().name() : dof_curl_vector.fieldTag().name())+ " ()";
353  this->setName(n);
354 }
355 
356 //**********************************************************************
357 template<typename EvalT, typename TRAITS>
359 DOFCurl(const PHX::FieldTag & input,
360  const PHX::FieldTag & output,
361  const panzer::BasisDescriptor & bd,
363  int basis_dim)
364  : use_descriptors_(true)
365  , bd_(bd)
366  , id_(id)
367  , dof_value(input)
368 {
369  TEUCHOS_ASSERT(bd_.getType()=="HCurl");
370 
371  basis_dimension = basis_dim; // user specified
372 
373  // build dof_curl
374  if(basis_dimension==2) {
375  dof_curl_scalar = output;
376  this->addEvaluatedField(dof_curl_scalar);
377  }
378  else if(basis_dimension==3) {
379  dof_curl_vector = output;
380  this->addEvaluatedField(dof_curl_vector);
381  }
382  else
383  { TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,"DOFCurl only works for 2D and 3D basis functions"); }
384 
385  // add to evaluation graph
386  this->addDependentField(dof_value);
387 
388  std::string n = "DOFCurl: " + (basis_dimension==2 ? dof_curl_scalar.fieldTag().name() : dof_curl_vector.fieldTag().name())+ " ()";
389  this->setName(n);
390 }
391 
392 //**********************************************************************
393 template<typename EvalT, typename TRAITS>
395 postRegistrationSetup(typename TRAITS::SetupData sd,
397 {
398  this->utils.setFieldData(dof_value,fm);
399  if(basis_dimension==3)
400  this->utils.setFieldData(dof_curl_vector,fm);
401  else
402  this->utils.setFieldData(dof_curl_scalar,fm);
403 
404  if(not use_descriptors_)
405  basis_index = panzer::getBasisIndex(basis_name, (*sd.worksets_)[0], this->wda);
406 }
407 
408 //**********************************************************************
409 template<typename EvalT, typename TRAITS>
411 evaluateFields(typename TRAITS::EvalData workset)
412 {
413  const panzer::BasisValues2<double> & basisValues = use_descriptors_ ? this->wda(workset).getBasisValues(bd_,id_)
414  : *this->wda(workset).bases[basis_index];
415 
416  if(basis_dimension==3) {
417  EvaluateCurlWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(dof_value,dof_curl_vector,basisValues.curl_basis_vector);
418  Kokkos::parallel_for(workset.num_cells,functor);
419  }
420  else {
421  EvaluateCurlWithSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(dof_value,dof_curl_scalar,basisValues.curl_basis_scalar);
422  Kokkos::parallel_for(workset.num_cells,functor);
423  }
424 }
425 
426 //**********************************************************************
427 
428 //**********************************************************************
429 // JACOBIAN EVALUATION TYPES
430 //**********************************************************************
431 
432 //**********************************************************************
433 template<typename TRAITS>
436  use_descriptors_(false),
437  dof_value( p.get<std::string>("Name"),
438  p.get< Teuchos::RCP<panzer::BasisIRLayout> >("Basis")->functional),
439  basis_name(p.get< Teuchos::RCP<panzer::BasisIRLayout> >("Basis")->name())
440 {
442  = p.get< Teuchos::RCP<BasisIRLayout> >("Basis")->getBasis();
443 
444  // do you specialize because you know where the basis functions are and can
445  // skip a large number of AD calculations?
446  if(p.isType<Teuchos::RCP<const std::vector<int> > >("Jacobian Offsets Vector")) {
447  offsets = *p.get<Teuchos::RCP<const std::vector<int> > >("Jacobian Offsets Vector");
448 
449  // allocate and copy offsets vector to Kokkos array
450  Kokkos::View<int*,PHX::Device> offsets_array_nc
451  = Kokkos::View<int*,PHX::Device>("offsets",offsets.size());
452  for(std::size_t i=0;i<offsets.size();i++)
453  offsets_array_nc(i) = offsets[i];
454  offsets_array = offsets_array_nc;
455 
456  accelerate_jacobian = true; // short cut for identity matrix
457  }
458  else
459  accelerate_jacobian = false; // don't short cut for identity matrix
460 
461  // Verify that this basis supports the curl operation
462  TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsCurl(),std::logic_error,
463  "DOFCurl: Basis of type \"" << basis->name() << "\" does not support CURL");
464  TEUCHOS_TEST_FOR_EXCEPTION(!basis->requiresOrientations(),std::logic_error,
465  "DOFCurl: Basis of type \"" << basis->name() << "\" in DOF Curl should require orientations. So we are throwing.");
466 
467  // build dof_curl
468  basis_dimension = basis->dimension();
469  if(basis_dimension==2) {
470  dof_curl_scalar = PHX::MDField<ScalarT,Cell,Point>(p.get<std::string>("Curl Name"),
471  p.get< Teuchos::RCP<panzer::IntegrationRule> >("IR")->dl_scalar );
472  this->addEvaluatedField(dof_curl_scalar);
473  }
474  else if(basis_dimension==3) {
475  dof_curl_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(p.get<std::string>("Curl Name"),
476  p.get< Teuchos::RCP<panzer::IntegrationRule> >("IR")->dl_vector );
477  this->addEvaluatedField(dof_curl_vector);
478  }
479  else
480  { TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,"DOFCurl only works for 2D and 3D basis functions"); }
481 
482  // add to evaluation graph
483  this->addDependentField(dof_value);
484 
485  std::string n = "DOFCurl: " + (basis_dimension==2 ? dof_curl_scalar.fieldTag().name() : dof_curl_vector.fieldTag().name())+ " (Jacobian)";
486  this->setName(n);
487 }
488 
489 //**********************************************************************
490 template<typename TRAITS>
492 DOFCurl(const PHX::FieldTag & input,
493  const PHX::FieldTag & output,
494  const panzer::BasisDescriptor & bd,
496  int basis_dim)
497  : use_descriptors_(true)
498  , bd_(bd)
499  , id_(id)
500  , dof_value(input)
501 {
502  TEUCHOS_ASSERT(bd_.getType()=="HCurl");
503 
504  basis_dimension = basis_dim; // user specified
505 
506  accelerate_jacobian = false; // don't short cut for identity matrix
507 
508  // build dof_curl
509  if(basis_dimension==2) {
510  dof_curl_scalar = output;
511  this->addEvaluatedField(dof_curl_scalar);
512  }
513  else if(basis_dimension==3) {
514  dof_curl_vector = output;
515  this->addEvaluatedField(dof_curl_vector);
516  }
517  else
518  { TEUCHOS_TEST_FOR_EXCEPTION(true,std::logic_error,"DOFCurl only works for 2D and 3D basis functions"); }
519 
520  // add to evaluation graph
521  this->addDependentField(dof_value);
522 
523  std::string n = "DOFCurl: " + (basis_dimension==2 ? dof_curl_scalar.fieldTag().name() : dof_curl_vector.fieldTag().name())+ " (Jacobian)";
524  this->setName(n);
525 }
526 
527 //**********************************************************************
528 template<typename TRAITS>
530 postRegistrationSetup(typename TRAITS::SetupData sd,
532 {
533  this->utils.setFieldData(dof_value,fm);
534  if(basis_dimension==3)
535  this->utils.setFieldData(dof_curl_vector,fm);
536  else
537  this->utils.setFieldData(dof_curl_scalar,fm);
538 
539  if(not use_descriptors_)
540  basis_index = panzer::getBasisIndex(basis_name, (*sd.worksets_)[0], this->wda);
541 }
542 
543 template<typename TRAITS>
545 evaluateFields(typename TRAITS::EvalData workset)
546 {
547  const panzer::BasisValues2<double> & basisValues = use_descriptors_ ? this->wda(workset).getBasisValues(bd_,id_)
548  : *this->wda(workset).bases[basis_index];
549 
550  if(!accelerate_jacobian) {
551  if(basis_dimension==3) {
552  EvaluateCurlWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(dof_value,dof_curl_vector,basisValues.curl_basis_vector);
553  Kokkos::parallel_for(workset.num_cells,functor);
554  }
555  else {
556  EvaluateCurlWithSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(dof_value,dof_curl_scalar,basisValues.curl_basis_scalar);
557  Kokkos::parallel_for(workset.num_cells,functor);
558  }
559 
560  return;
561  }
562  else {
563 
564  if(basis_dimension==3) {
565  EvaluateCurlFastSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(dof_value,dof_curl_vector,offsets_array,basisValues.curl_basis_vector);
566  Kokkos::parallel_for(workset.num_cells,functor);
567  }
568  else {
569  EvaluateCurlFastSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(dof_value,dof_curl_scalar,offsets_array,basisValues.curl_basis_scalar);
570  Kokkos::parallel_for(workset.num_cells,functor);
571  }
572  }
573 }
574 
575 }
576 
577 #endif
std::string basis_name
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
panzer::BasisDescriptor bd_
T & get(const std::string &name, T def_value)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Array_CellBasisIPDim curl_basis_vector
int numPoints
PHX::MDField< ScalarT, Cell, Point > dof_curl_scalar
int numFields
PHX::MDField< ScalarT, Cell, Point, Dim > dof_curl_vector
const std::string & getType() const
Get type of basis.
std::vector< std::string >::size_type getBasisIndex(std::string basis_name, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
Returns the index in the workset bases for a particular BasisIRLayout name.
std::size_t basis_index
Array curl_basis
DOFCurl(const Teuchos::ParameterList &p)
panzer::IntegrationDescriptor id_
bool isType(const std::string &name) const
#define TEUCHOS_ASSERT(assertion_test)
void evaluateFields(typename TRAITS::EvalData d)
PHX::MDField< const ScalarT, Cell, Point > dof_value
PHX::MDField< ScalarT, Cell, Point, Dim > dof_curl
Kokkos::View< const int *, PHX::Device > offsets
Array_CellBasisIP curl_basis_scalar