Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_Integrator_CurlBasisDotVector_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_Integrator_CurlBasisDotVector_impl_hpp
44 #define Panzer_Integrator_CurlBasisDotVector_impl_hpp
45 
47 //
48 // Include Files
49 //
51 
52 // Intrepid2
53 #include "Intrepid2_FunctionSpaceTools.hpp"
54 
55 // Panzer
56 #include "Panzer_BasisIRLayout.hpp"
60 
61 // Phalanx
62 #include "Phalanx_KokkosDeviceTypes.hpp"
63 
64 namespace panzer
65 {
67  //
68  // Main Constructor
69  //
71  template<typename EvalT, typename Traits>
75  const std::string& resName,
76  const std::string& valName,
77  const panzer::BasisIRLayout& basis,
78  const panzer::IntegrationRule& ir,
79  const double& multiplier, /* = 1 */
80  const std::vector<std::string>& fmNames /* =
81  std::vector<std::string>() */)
82  :
83  evalStyle_(evalStyle),
84  useDescriptors_(false),
85  multiplier_(multiplier),
86  basisName_(basis.name())
87  {
88  using Kokkos::View;
89  using panzer::BASIS;
90  using panzer::Cell;
91  using panzer::Dim;
93  using panzer::IP;
94  using panzer::PureBasis;
95  using PHX::MDField;
96  using std::invalid_argument;
97  using std::logic_error;
98  using std::string;
99  using Teuchos::RCP;
100 
101  // Ensure the input makes sense.
102  TEUCHOS_TEST_FOR_EXCEPTION(resName == "", invalid_argument, "Error: " \
103  "Integrator_CurlBasisDotVector called with an empty residual name.")
104  TEUCHOS_TEST_FOR_EXCEPTION(valName == "", invalid_argument, "Error: " \
105  "Integrator_CurlBasisDotVector called with an empty value name.")
106  RCP<const PureBasis> tmpBasis = basis.getBasis();
107  TEUCHOS_TEST_FOR_EXCEPTION(not tmpBasis->isVectorBasis(), logic_error,
108  "Error: Integrator_CurlBasisDotVector: Basis of type \""
109  << tmpBasis->name() << "\" is not a vector basis.")
110  TEUCHOS_TEST_FOR_EXCEPTION(not tmpBasis->requiresOrientations(),
111  logic_error, "Error: Integrator_CurlBasisDotVector: Basis of type \""
112  << tmpBasis->name() << "\" does not require orientations, though it " \
113  "should for its use in this Evaluator.")
114  TEUCHOS_TEST_FOR_EXCEPTION(not tmpBasis->supportsCurl(), logic_error,
115  "Error: Integrator_CurlBasisDotVector: Basis of type \""
116  << tmpBasis->name() << "\" does not support curl.")
117  TEUCHOS_TEST_FOR_EXCEPTION(not ((tmpBasis->dimension() == 2) or
118  (tmpBasis->dimension() == 3)), logic_error,
119  "Error: Integrator_CurlBasisDotVector requires either a two- or " \
120  "three-dimensional basis. The basis \"" << tmpBasis->name()
121  << "\" is neither.")
122 
123  // Use a scalar field only if we're dealing with a two-dimensional case.
124  spaceDim_ = tmpBasis->dimension();
125 
126  // Create the field for the vector-values quantity we're integrating.
127  if (spaceDim_ == 2)
128  {
129  vector2D_ = MDField<const ScalarT, Cell, IP>(valName, ir.dl_scalar);
130  this->addDependentField(vector2D_);
131  }
132  else // if (spaceDim_ == 3)
133  {
134  vector3D_ = MDField<const ScalarT, Cell, IP, Dim>(valName, ir.dl_vector);
135  this->addDependentField(vector3D_);
136  } // end if spaceDim_ is 2 or 3
137 
138  // Create the field that we're either contributing to or evaluating
139  // (storing).
140  field_ = MDField<ScalarT, Cell, BASIS>(resName, basis.functional);
141  if (evalStyle == EvaluatorStyle::CONTRIBUTES)
142  this->addContributedField(field_);
143  else // if (evalStyle == EvaluatorStyle::EVALUATES)
144  this->addEvaluatedField(field_);
145 
146  // Add the dependent field multipliers, if there are any.
147  int i(0);
148  fieldMults_.resize(fmNames.size());
149  kokkosFieldMults_ = View<View<const ScalarT**,typename PHX::DevLayout<ScalarT>::type,PHX::Device>*>(
150  "CurlBasisDotVector::KokkosFieldMultipliers", fmNames.size());
151  for (const auto& name : fmNames)
152  {
153  fieldMults_[i++] = MDField<const ScalarT, Cell, IP>(name, ir.dl_scalar);
154  this->addDependentField(fieldMults_[i - 1]);
155  } // end loop over the field multipliers
156 
157  // Set the name of this object.
158  string n("Integrator_CurlBasisDotVector (");
159  if (evalStyle == EvaluatorStyle::CONTRIBUTES)
160  n += "CONTRIBUTES";
161  else // if (evalStyle == EvaluatorStyle::EVALUATES)
162  n += "EVALUATES";
163  n += "): " + field_.fieldTag().name();
164  this->setName(n);
165  } // end of Main Constructor
166 
168  //
169  // ParameterList Constructor
170  //
172  template<typename EvalT, typename Traits>
175  const Teuchos::ParameterList& p)
176  :
178  panzer::EvaluatorStyle::EVALUATES,
179  p.get<std::string>("Residual Name"),
180  p.get<std::string>("Value Name"),
181  (*p.get<Teuchos::RCP<panzer::BasisIRLayout>>("Basis")),
182  (*p.get<Teuchos::RCP<panzer::IntegrationRule>>("IR")),
183  p.get<double>("Multiplier"),
184  p.isType<Teuchos::RCP<const std::vector<std::string>>>
185  ("Field Multipliers") ?
186  (*p.get<Teuchos::RCP<const std::vector<std::string>>>
187  ("Field Multipliers")) : std::vector<std::string>())
188  {
190  using Teuchos::RCP;
191 
192  // Ensure that the input ParameterList didn't contain any bogus entries.
193  RCP<ParameterList> validParams = this->getValidParameters();
194  p.validateParameters(*validParams);
195  } // end of ParameterList Constructor
196 
198  //
199  // FieldTag Constructor
200  //
202  template<typename EvalT, typename TRAITS>
206  const PHX::FieldTag& resTag,
207  const PHX::FieldTag& valTag,
208  const panzer::BasisDescriptor& bd,
210  const int& spaceDim, /* = 3 */
211  const double& multiplier, /* = 1 */
212  const std::vector<PHX::FieldTag>& multipliers /* =
213  std::vector<PHX::FieldTag>() */)
214  :
215  evalStyle_(evalStyle),
216  useDescriptors_(true),
217  bd_(bd),
218  id_(id),
219  multiplier_(multiplier),
220  spaceDim_(spaceDim)
221  {
222  using Kokkos::View;
224  using std::logic_error;
225  using std::string;
226 
227  // Ensure the input makes sense.
228  TEUCHOS_TEST_FOR_EXCEPTION(bd_.getType() != "HCurl", logic_error,
229  "Error: Integrator_CurlBasisDotVector: Basis of type \""
230  << bd_.getType() << "\" does not support curl.")
231  TEUCHOS_TEST_FOR_EXCEPTION(not ((spaceDim == 2) or (spaceDim == 3)),
232  logic_error, "Error: Integrator_CurlBasisDotVector works on either " \
233  "two- or three-dimensional problems. You've provided spaceDim = "
234  << spaceDim << ".")
235 
236  // Create the field for the vector-valued quantity we're integrating.
237  if (spaceDim_ == 2)
238  {
239  vector2D_ = valTag;
240  this->addDependentField(vector2D_);
241  }
242  else // if (spaceDim_ == 3)
243  {
244  vector3D_ = valTag;
245  this->addDependentField(vector3D_);
246  } // end if spaceDim_ is 2 or 3
247 
248  // Create the field that we're either contributing to or evaluating
249  // (storing).
250  field_ = resTag;
251  if (evalStyle == EvaluatorStyle::CONTRIBUTES)
252  this->addContributedField(field_);
253  else // if (evalStyle == EvaluatorStyle::EVALUATES)
254  this->addEvaluatedField(field_);
255 
256  // Add the dependent field multipliers, if there are any.
257  int i(0);
258  fieldMults_.resize(multipliers.size());
259  kokkosFieldMults_ = View<View<const ScalarT**,typename PHX::DevLayout<ScalarT>::type,PHX::Device>*>(
260  "CurlBasisDotVector::KokkosFieldMultipliers", multipliers.size());
261  for (const auto& fm : multipliers)
262  {
263  fieldMults_[i++] = fm;
264  this->addDependentField(fieldMults_[i - 1]);
265  } // end loop over the field multipliers
266 
267  // Set the name of this object.
268  string n("Integrator_CurlBasisDotVector (");
269  if (evalStyle == EvaluatorStyle::CONTRIBUTES)
270  n += "CONTRIBUTES";
271  else // if (evalStyle == EvaluatorStyle::EVALUATES)
272  n += "EVALUATES";
273  n += "): " + field_.fieldTag().name();
274  this->setName(n);
275  } // end of Other Constructor
276 
278  //
279  // postRegistrationSetup()
280  //
282  template<typename EvalT, typename Traits>
283  void
286  typename Traits::SetupData sd,
288  {
289  using panzer::getBasisIndex;
290  using PHX::MDField;
291  using std::vector;
292 
293  // Get the Kokkos::Views of the field multipliers.
294  for (size_t i(0); i < fieldMults_.size(); ++i)
295  kokkosFieldMults_(i) = fieldMults_[i].get_static_view();
296 
297  // Determine the index in the Workset bases for our particular basis
298  // name.
299  if (not useDescriptors_)
300  basisIndex_ = getBasisIndex(basisName_, (*sd.worksets_)[0], this->wda);
301 
302  // Set up the field that will be used to build of the result of this
303  // integration.
304  MDFieldArrayFactory af("",
305  fm.template getKokkosExtendedDataTypeDimensions<EvalT>(), true);
306  if (spaceDim_ == 2)
307  result2D_ = af.buildStaticArray<ScalarT, Cell, IP>(
308  "Integrator_CurlBasisDotVector: 2-D Result", vector2D_.extent(0),
309  vector2D_.extent(1));
310  else // if (spaceDim_ == 3)
311  result3D_ = af.buildStaticArray<ScalarT, Cell, IP, Dim>(
312  "Integrator_CurlBasisDotVector: 3-D Result", vector3D_.extent(0),
313  vector3D_.extent(1), 3);
314  } // end of postRegistrationSetup()
315 
317  //
318  // Anonymous namespace containing classes for performing the integration.
319  //
321  namespace
322  {
334  template<typename ScalarT>
335  class PreMultiply2D
336  {
337  public:
338 
343  struct ScalarMultiplierTag
344  {
345  }; // end of struct ScalarMultiplierTag
346 
351  struct FieldMultiplierTag
352  {
353  }; // end of struct FieldMultiplierTag
354 
366  KOKKOS_INLINE_FUNCTION
367  void operator()(const ScalarMultiplierTag, const unsigned cell) const
368  {
369  int numQP(result.extent(1));
370  for (int qp(0); qp < numQP; ++qp)
371  result(cell, qp) = multiplier * vectorField(cell, qp);
372  } // end of ScalarMultiplierTag operator()()
373 
386  KOKKOS_INLINE_FUNCTION
387  void operator()(const FieldMultiplierTag, const unsigned cell) const
388  {
389  int numQP(result.extent(1));
390  for (int qp(0); qp < numQP; ++qp)
391  result(cell, qp) *= fieldMult(cell, qp);
392  } // end of FieldMultiplierTag operator()()
393 
398  using execution_space = PHX::Device;
399 
404  PHX::MDField<ScalarT, panzer::Cell, panzer::IP> result;
405 
410  double multiplier;
411 
416  PHX::MDField<const ScalarT, panzer::Cell, panzer::IP> vectorField;
417 
422  PHX::MDField<const ScalarT, panzer::Cell, panzer::IP> fieldMult;
423  }; // end of class PreMultiply2D
424 
436  template<typename ScalarT>
437  class PreMultiply3D
438  {
439  public:
440 
445  struct ScalarMultiplierTag
446  {
447  }; // end of struct ScalarMultiplierTag
448 
453  struct FieldMultiplierTag
454  {
455  }; // end of struct FieldMultiplierTag
456 
468  KOKKOS_INLINE_FUNCTION
469  void operator()(const ScalarMultiplierTag, const unsigned cell) const
470  {
471  int numQP(result.extent(1)), numDim(result.extent(2));
472  for (int qp(0); qp < numQP; ++qp)
473  for (int dim(0); dim < numDim; ++dim)
474  result(cell, qp, dim) = multiplier * vectorField(cell, qp, dim);
475  } // end of ScalarMultiplierTag operator()()
476 
489  KOKKOS_INLINE_FUNCTION
490  void operator()(const FieldMultiplierTag, const unsigned cell) const
491  {
492  int numQP(result.extent(1)), numDim(result.extent(2));
493  for (int qp(0); qp < numQP; ++qp)
494  for (int dim(0); dim < numDim; ++dim)
495  result(cell, qp, dim) *= fieldMult(cell, qp);
496  } // end of FieldMultiplierTag operator()()
497 
502  using execution_space = PHX::Device;
503 
508  PHX::MDField<ScalarT, panzer::Cell, panzer::IP, panzer::Dim> result;
509 
514  double multiplier;
515 
520  PHX::MDField<const ScalarT, panzer::Cell, panzer::IP, panzer::Dim>
521  vectorField;
522 
527  PHX::MDField<const ScalarT, panzer::Cell, panzer::IP> fieldMult;
528  }; // end of class PreMultiply3D
529 
541  template<typename ScalarT, int spaceDim>
542  class Integrate3D
543  {
544  public:
545 
558  KOKKOS_INLINE_FUNCTION
559  void operator()(const unsigned cell) const
560  {
562  int numBases(weightedCurlBasis.extent(1)),
563  numQP(weightedCurlBasis.extent(2));
564  for (int basis(0); basis < numBases; ++basis)
565  {
567  field(cell, basis) = 0.0;
568  for (int qp(0); qp < numQP; ++qp)
569  for (int dim(0); dim < spaceDim; ++dim)
570  field(cell, basis) += result(cell, qp, dim) *
571  weightedCurlBasis(cell, basis, qp, dim);
572  } // end loop over the basis functions
573  } // end of operator()
574 
579  using execution_space = PHX::Device;
580 
584  PHX::MDField<ScalarT, panzer::Cell, panzer::IP, panzer::Dim> result;
585 
590  PHX::MDField<ScalarT, panzer::Cell, panzer::BASIS> field;
591 
595  PHX::MDField<double, panzer::Cell, panzer::BASIS, panzer::IP,
596  panzer::Dim> weightedCurlBasis;
597 
603  }; // end of class Integrate3D
604 
616  template<typename ScalarT>
617  class Integrate2D
618  {
619  public:
620 
632  KOKKOS_INLINE_FUNCTION
633  void operator()(const unsigned cell) const
634  {
636  int numBases(weightedCurlBasis.extent(1)),
637  numQP(weightedCurlBasis.extent(2));
638  for (int basis(0); basis < numBases; ++basis)
639  {
641  field(cell, basis) = 0.0;
642  for (int qp(0); qp < numQP; ++qp)
643  field(cell, basis) += result(cell, qp) *
644  weightedCurlBasis(cell, basis, qp);
645  } // end loop over the basis functions
646  } // end of operator()
647 
652  using execution_space = PHX::Device;
653 
657  PHX::MDField<ScalarT, panzer::Cell, panzer::IP> result;
658 
663  PHX::MDField<ScalarT, panzer::Cell, panzer::BASIS> field;
664 
668  PHX::MDField<double, panzer::Cell, panzer::BASIS, panzer::IP>
670 
676  }; // end of class Integrate2D
677 
678  } // end of anonymous namespace
679 
681  //
682  // evaluateFields()
683  //
685  template<typename EvalT, typename Traits>
686  void
689  typename Traits::EvalData workset)
690  {
691  using Kokkos::parallel_for;
692  using Kokkos::RangePolicy;
693  using panzer::BasisValues2;
694  using PHX::Device;
695  using PHX::MDField;
696  using std::vector;
697 
698  // Grab the basis information.
699  const BasisValues2<double>& bv = useDescriptors_ ?
700  this->wda(workset).getBasisValues(bd_, id_) :
701  *this->wda(workset).bases[basisIndex_];
702 
703  // If we're dealing with a two- or three-dimensional problem...
704  if (spaceDim_ == 2)
705  {
706  // Create an object to multiply the integrand by the scalar multiplier
707  // and any field multipliers out in front of the integral.
708  using PreMultiply = PreMultiply2D<ScalarT>;
709  using ScalarMultiplierTag = typename PreMultiply::ScalarMultiplierTag;
710  using FieldMultiplierTag = typename PreMultiply::FieldMultiplierTag;
711  PreMultiply preMultiply;
712  preMultiply.result = result2D_;
713  preMultiply.multiplier = multiplier_;
714  preMultiply.vectorField = vector2D_;
715 
716  // Multiply the integrand by the scalar multiplier out in front of the
717  // integral.
718  parallel_for(RangePolicy<Device, ScalarMultiplierTag>(0,
719  workset.num_cells), preMultiply);
720 
721  // Multiply the integrand by any field multipliers out in front of the
722  // integral.
723  for (const auto& field : fieldMults_)
724  {
725  preMultiply.fieldMult = field;
726  parallel_for(RangePolicy<Device, FieldMultiplierTag>(0,
727  workset.num_cells), preMultiply);
728  } // end loop over the field multipliers
729 
730  // Create an object to do the actual integration and then do it.
731  Integrate2D<ScalarT> integrate;
732  integrate.result = result2D_;
733  integrate.field = field_;
734  integrate.weightedCurlBasis = bv.weighted_curl_basis_scalar;
735  integrate.evalStyle = evalStyle_;
736  parallel_for(workset.num_cells, integrate);
737  }
738  else // if (spaceDim_ == 3)
739  {
740  // Create an object to multiply the integrand by the scalar multiplier
741  // and any field multipliers out in front of the integral.
742  using PreMultiply = PreMultiply3D<ScalarT>;
743  using ScalarMultiplierTag = typename PreMultiply::ScalarMultiplierTag;
744  using FieldMultiplierTag = typename PreMultiply::FieldMultiplierTag;
745  PreMultiply preMultiply;
746  preMultiply.result = result3D_;
747  preMultiply.multiplier = multiplier_;
748  preMultiply.vectorField = vector3D_;
749 
750  // Multiply the integrand by the scalar multiplier out in front of the
751  // integral.
752  parallel_for(RangePolicy<Device, ScalarMultiplierTag>(0,
753  workset.num_cells), preMultiply);
754 
755  // Multiply the integrand by any field multipliers out in front of the
756  // integral.
757  for (const auto& field : fieldMults_)
758  {
759  preMultiply.fieldMult = field;
760  parallel_for(RangePolicy<Device, FieldMultiplierTag>(0,
761  workset.num_cells), preMultiply);
762  } // end loop over the field multipliers
763 
764  // Create an object to do the actual integration and then do it.
765  Integrate3D<ScalarT, 3> integrate;
766  integrate.result = result3D_;
767  integrate.field = field_;
768  integrate.weightedCurlBasis = bv.weighted_curl_basis_vector;
769  integrate.evalStyle = evalStyle_;
770  parallel_for(workset.num_cells, integrate);
771  } // end if spaceDim_ is 2 or 3
772  } // end of evaluateFields()
773 
775  //
776  // getValidParameters()
777  //
779  template<typename EvalT, typename TRAITS>
783  {
784  using panzer::BasisIRLayout;
786  using std::string;
787  using std::vector;
789  using Teuchos::RCP;
790  using Teuchos::rcp;
791 
792  RCP<ParameterList> p = rcp(new ParameterList);
793  p->set<string>("Residual Name", "?");
794  p->set<string>("Value Name", "?");
795  RCP<BasisIRLayout> basis;
796  p->set("Basis", basis);
798  p->set("IR", ir);
799  p->set<double>("Multiplier", 1.0);
801  p->set("Field Multipliers", fms);
802  return p;
803  } // end of getValidParameters()
804 
805 } // end of namespace panzer
806 
807 #endif // Panzer_Integrator_CurlBasisDotVector_impl_hpp
Array_CellBasisIPDim weighted_curl_basis_vector
Kokkos::View< Kokkos::View< const ScalarT **, typename PHX::DevLayout< ScalarT >::type, PHX::Device > * > kokkosFieldMults_
The Kokkos::View representation of the (possibly empty) list of fields that are multipliers out in fr...
void evaluateFields(typename Traits::EvalData d)
Evaluate Fields.
std::vector< PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > > fieldMults_
The (possibly empty) list of fields that are multipliers out in front of the integral ( ...
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > vectorField
A field representing the vector-valued function we&#39;re integrating ( ).
Teuchos::RCP< Teuchos::ParameterList > getValidParameters() const
Get Valid Parameters.
EvaluatorStyle
An indication of how an Evaluator will behave.
PHX::MDField< Scalar, T0 > buildStaticArray(const std::string &str, int d0) const
Teuchos::RCP< const PureBasis > getBasis() const
panzer::BasisDescriptor bd_
The BasisDescriptor for the basis to use.
PHX::MDField< ScalarT, panzer::Cell, panzer::IP > result
A field that will be used to build up the result of the integral we&#39;re performing.
double multiplier
The scalar multiplier out in front of the integral ( ).
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
const std::string & getType() const
Get type of basis.
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
panzer::EvaluatorStyle evalStyle
The EvaluatorStyle of the parent Integrator_CurlBasisDotVector object.
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP, panzer::Dim > vector3D_
A field representing the vector-valued function we&#39;re integrating ( ) in a three-dimensional problem...
int spaceDim_
The spatial dimension of the vector-valued function we&#39;re integrating, either 2 or 3...
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.
PHX::MDField< double, panzer::Cell, panzer::BASIS, panzer::IP, panzer::Dim > weightedCurlBasis
The vector basis information necessary for integration.
Array_CellBasisIP weighted_curl_basis_scalar
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
void postRegistrationSetup(typename Traits::SetupData sd, PHX::FieldManager< Traits > &fm)
Post-Registration Setup.
const std::vector< std::pair< int, LocalOrdinal > > &pid_and_lid const
Teuchos::RCP< PHX::DataLayout > dl_vector
Data layout for vector fields.
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field_
A field to which we&#39;ll contribute, or in which we&#39;ll store, the result of computing this integral...
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we&#39;ll contribute, or in which we&#39;ll store, the result of computing this integral...
Integrator_CurlBasisDotVector(const panzer::EvaluatorStyle &evalStyle, const std::string &resName, const std::string &valName, const panzer::BasisIRLayout &basis, const panzer::IntegrationRule &ir, const double &multiplier=1, const std::vector< std::string > &fmNames=std::vector< std::string >())
Main Constructor.
Description and data layouts associated with a particular basis.
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > fieldMult
One of the field multipliers ( , , etc.) out in front of the integral.
Teuchos::RCP< PHX::DataLayout > functional
&lt;Cell,Basis&gt;
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > vector2D_
A field representing the vector-valued function we&#39;re integrating ( ) in a two-dimensional problem...
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_
typename EvalT::ScalarT ScalarT
The scalar type.