Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_DOF_PointValues_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_IMPL_HPP
44#define PANZER_DOF_IMPL_HPP
45
46#include <algorithm>
53
54#include "Intrepid2_FunctionSpaceTools.hpp"
55
56namespace panzer {
57
58//**********************************************************************
59//* DOF_PointValues evaluator
60//**********************************************************************
61
62//**********************************************************************
63// MOST EVALUATION TYPES
64//**********************************************************************
65
66//**********************************************************************
67template<typename EvalT, typename TRAITS>
69DOF_PointValues(const Teuchos::ParameterList & p)
70{
71 const std::string fieldName = p.get<std::string>("Name");
72 basis = p.get< Teuchos::RCP<const PureBasis> >("Basis");
73 Teuchos::RCP<const PointRule> pointRule = p.get< Teuchos::RCP<const PointRule> >("Point Rule");
74 is_vector_basis = basis->isVectorBasis();
75
76 std::string evalName = fieldName+"_"+pointRule->getName();
77 if(p.isType<bool>("Use DOF Name")) {
78 if(p.get<bool>("Use DOF Name"))
79 evalName = fieldName;
80 }
81
82 dof_basis = PHX::MDField<const ScalarT,Cell,Point>(fieldName, basis->functional);
83
84 this->addDependentField(dof_basis);
85
86 // setup all basis fields that are required
87 Teuchos::RCP<BasisIRLayout> layout = Teuchos::rcp(new BasisIRLayout(basis,*pointRule));
88 basisValues = Teuchos::rcp(new BasisValues2<double>(basis->name()+"_"+pointRule->getName()+"_"));
89 basisValues->setupArrays(layout,false);
90
91 // the field manager will allocate all of these field
92 // swap between scalar basis value, or vector basis value
93 if(basis->isScalarBasis()) {
94 dof_ip_scalar = PHX::MDField<ScalarT,Cell,Point>(
95 evalName,
96 pointRule->dl_scalar);
97 this->addEvaluatedField(dof_ip_scalar);
98 this->addNonConstDependentField(basisValues->basis_ref_scalar);
99 this->addNonConstDependentField(basisValues->basis_scalar);
100 }
101 else if(basis->isVectorBasis()) {
102 dof_ip_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(
103 evalName,
104 pointRule->dl_vector);
105 this->addEvaluatedField(dof_ip_vector);
106 this->addNonConstDependentField(basisValues->basis_ref_vector);
107 this->addNonConstDependentField(basisValues->basis_vector);
108 }
109 else
110 { TEUCHOS_ASSERT(false); }
111
112 std::string n = "DOF_PointValues: " + dof_basis.fieldTag().name();
113 this->setName(n);
114}
115
116//**********************************************************************
117template<typename EvalT, typename TRAITS>
119postRegistrationSetup(typename TRAITS::SetupData /* sd */,
121{
122 if(!is_vector_basis) {
123 this->utils.setFieldData(basisValues->basis_ref_scalar,fm);
124 this->utils.setFieldData(basisValues->basis_scalar,fm);
125 }
126 else {
127 this->utils.setFieldData(basisValues->basis_ref_vector,fm);
128 this->utils.setFieldData(basisValues->basis_vector,fm);
129 }
130}
131
132//**********************************************************************
133template<typename EvalT, typename TRAITS>
135evaluateFields(typename TRAITS::EvalData workset)
136{
137 const int vector_size = panzer::HP::inst().vectorSize<ScalarT>();
138
139 if(is_vector_basis) {
140 int spaceDim = basisValues->basis_vector.extent(3);
141 if(spaceDim==3) {
142 dof_functors::EvaluateDOFWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(dof_basis.get_static_view(),dof_ip_vector.get_static_view(),basisValues->basis_vector);
143 Kokkos::parallel_for(Kokkos::TeamPolicy<PHX::Device>(workset.num_cells,Kokkos::AUTO(),vector_size),functor);
144 }
145 else {
146 dof_functors::EvaluateDOFWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,2> functor(dof_basis.get_static_view(),dof_ip_vector.get_static_view(),basisValues->basis_vector);
147 Kokkos::parallel_for(Kokkos::TeamPolicy<PHX::Device>(workset.num_cells,Kokkos::AUTO(),vector_size),functor);
148 }
149 }
150 else {
151 dof_functors::EvaluateDOFWithSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(dof_basis,dof_ip_scalar,basisValues->basis_scalar);
152 Kokkos::parallel_for(workset.num_cells,functor);
153 }
154}
155
156//**********************************************************************
157
158//**********************************************************************
159// JACOBIAN EVALUATION TYPES
160//**********************************************************************
161
162//**********************************************************************
163template<typename TRAITS>
165DOF_PointValues(const Teuchos::ParameterList & p)
166{
167 const std::string fieldName = p.get<std::string>("Name");
168 basis = p.get< Teuchos::RCP<const PureBasis> >("Basis");
169 Teuchos::RCP<const PointRule> pointRule = p.get< Teuchos::RCP<const PointRule> >("Point Rule");
170 is_vector_basis = basis->isVectorBasis();
171
172 if(p.isType<Teuchos::RCP<const std::vector<int> > >("Jacobian Offsets Vector")) {
173 const std::vector<int> & offsets = *p.get<Teuchos::RCP<const std::vector<int> > >("Jacobian Offsets Vector");
174
175 // allocate and copy offsets vector to Kokkos array
176 offsets_array = PHX::View<int*>("offsets",offsets.size());
177 for(std::size_t i=0;i<offsets.size();i++)
178 offsets_array(i) = offsets[i];
179
180 accelerate_jacobian = true; // short cut for identity matrix
181 }
182 else
183 accelerate_jacobian = false; // don't short cut for identity matrix
184
185 std::string evalName = fieldName+"_"+pointRule->getName();
186 if(p.isType<bool>("Use DOF Name")) {
187 if(p.get<bool>("Use DOF Name"))
188 evalName = fieldName;
189 }
190
191 dof_basis = PHX::MDField<const ScalarT,Cell,Point>(fieldName, basis->functional);
192
193 this->addDependentField(dof_basis);
194
195 // setup all basis fields that are required
196 Teuchos::RCP<BasisIRLayout> layout = Teuchos::rcp(new BasisIRLayout(basis,*pointRule));
197 basisValues = Teuchos::rcp(new BasisValues2<double>(basis->name()+"_"+pointRule->getName()+"_"));
198 basisValues->setupArrays(layout,false);
199
200 // the field manager will allocate all of these field
201 // swap between scalar basis value, or vector basis value
202 if(basis->isScalarBasis()) {
203 dof_ip_scalar = PHX::MDField<ScalarT,Cell,Point>(
204 evalName,
205 pointRule->dl_scalar);
206 this->addEvaluatedField(dof_ip_scalar);
207 this->addNonConstDependentField(basisValues->basis_ref_scalar);
208 this->addNonConstDependentField(basisValues->basis_scalar);
209 }
210 else if(basis->isVectorBasis()) {
211 dof_ip_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(
212 evalName,
213 pointRule->dl_vector);
214 this->addEvaluatedField(dof_ip_vector);
215 this->addNonConstDependentField(basisValues->basis_ref_vector);
216 this->addNonConstDependentField(basisValues->basis_vector);
217 }
218 else
219 { TEUCHOS_ASSERT(false); }
220
221 std::string n = "DOF_PointValues: " + dof_basis.fieldTag().name() + " Jacobian";
222 this->setName(n);
223}
224
225//**********************************************************************
226template<typename TRAITS>
228postRegistrationSetup(typename TRAITS::SetupData /* sd */,
230{
231 if(!is_vector_basis) {
232 this->utils.setFieldData(basisValues->basis_ref_scalar,fm);
233 this->utils.setFieldData(basisValues->basis_scalar,fm);
234 }
235 else {
236 this->utils.setFieldData(basisValues->basis_ref_vector,fm);
237 this->utils.setFieldData(basisValues->basis_vector,fm);
238 }
239}
240
241//**********************************************************************
242template<typename TRAITS>
244evaluateFields(typename TRAITS::EvalData workset)
245{
246 const int vector_size = panzer::HP::inst().vectorSize<ScalarT>();
247
248 if(is_vector_basis) {
249 if(accelerate_jacobian) {
250 int spaceDim = basisValues->basis_vector.extent(3);
251 if(spaceDim==3) {
252 dof_functors::EvaluateDOFFastSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(dof_basis,dof_ip_vector,offsets_array,basisValues->basis_vector);
253 Kokkos::parallel_for(workset.num_cells,functor);
254 }
255 else {
256 dof_functors::EvaluateDOFFastSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,2> functor(dof_basis,dof_ip_vector,offsets_array,basisValues->basis_vector);
257 Kokkos::parallel_for(workset.num_cells,functor);
258 }
259 }
260 else {
261 int spaceDim = basisValues->basis_vector.extent(3);
262 if(spaceDim==3) {
263 dof_functors::EvaluateDOFWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(dof_basis.get_static_view(),dof_ip_vector.get_static_view(),basisValues->basis_vector);
264 Kokkos::parallel_for(Kokkos::TeamPolicy<PHX::Device>(workset.num_cells,Kokkos::AUTO(),vector_size),functor);
265 }
266 else {
267 dof_functors::EvaluateDOFWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,2> functor(dof_basis.get_static_view(),dof_ip_vector.get_static_view(),basisValues->basis_vector);
268 Kokkos::parallel_for(Kokkos::TeamPolicy<PHX::Device>(workset.num_cells,Kokkos::AUTO(),vector_size),functor);
269 }
270 }
271 }
272 else {
273 if(accelerate_jacobian) {
274 dof_functors::EvaluateDOFFastSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(dof_basis,dof_ip_scalar,offsets_array,basisValues->basis_scalar);
275 Kokkos::parallel_for(workset.num_cells,functor);
276 }
277 else {
278 dof_functors::EvaluateDOFWithSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(dof_basis,dof_ip_scalar,basisValues->basis_scalar);
279 Kokkos::parallel_for(workset.num_cells,functor);
280 }
281 }
282}
283
284}
285
286#endif
PHX::View< const int * > offsets
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
void evaluateFields(typename TRAITS::EvalData d)
DOF_PointValues(const Teuchos::ParameterList &p)
int vectorSize() const
Returns the vector size. Specialized for AD scalar types.
static HP & inst()
Private ctor.