Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_GatherIntegrationCoordinates_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_GATHER_INTEGRATION_COORDINATES_IMPL_HPP
44 #define PANZER_GATHER_INTEGRATION_COORDINATES_IMPL_HPP
45 
46 #include "Teuchos_Assert.hpp"
47 #include "Phalanx_DataLayout.hpp"
48 
51 
52 #include "Teuchos_FancyOStream.hpp"
53 
54 template<typename EvalT,typename TRAITS>
55 std::string
57 fieldName(int degree)
58 {
59  std::stringstream ss;
60  ss << "IR_" << degree << " IntegrationCoordinates";
61  return ss.str();
62 }
63 
64 template<typename EvalT,typename TRAITS>
67 {
68  quadDegree_ = quad.cubature_degree;
69 
70  quadCoordinates_ = PHX::MDField<ScalarT,Cell,Point,Dim>(fieldName(quadDegree_),quad.dl_vector);
71 
72  this->addEvaluatedField(quadCoordinates_);
73  this->addUnsharedField(Teuchos::rcp_const_cast<PHX::FieldTag>(quadCoordinates_.fieldTagPtr()));
74 
75  this->setName("GatherIntegrationCoordinates: "+fieldName(quadDegree_));
76 }
77 
78 // **********************************************************************
79 template<typename EvalT,typename TRAITS>
81 postRegistrationSetup(typename TRAITS::SetupData sd,
82  PHX::FieldManager<TRAITS>& /* fm */)
83 {
84  quadIndex_ = panzer::getIntegrationRuleIndex(quadDegree_, (*sd.worksets_)[0], this->wda);
85  // make sure to zero out derivative array as we only set the value for AD types in the loop below
86  Kokkos::deep_copy(quadCoordinates_.get_static_view(),0.0);
87 }
88 
89 // **********************************************************************
90 template<typename EvalT,typename TRAITS>
92 evaluateFields(typename TRAITS::EvalData workset)
93 {
94  // const Kokkos::DynRankView<double,PHX::Device> & quadCoords = this->wda(workset).int_rules[quadIndex_]->ip_coordinates;
95  const IntegrationValues2<double> & iv = *this->wda(workset).int_rules[quadIndex_];
96  auto s_ip_coordinates = iv.ip_coordinates.get_static_view();
97  auto d_quadCoordinates = quadCoordinates_.get_static_view();
98 
99  // just copy the array
100  Kokkos::MDRangePolicy<PHX::Device,Kokkos::Rank<3>> policy({0,0,0},{int(workset.num_cells),s_ip_coordinates.extent_int(1),s_ip_coordinates.extent_int(2)});
101  Kokkos::parallel_for("GatherIntegrationCoords", policy, KOKKOS_LAMBDA (const int i, const int j, const int k) {
102  auto s_ip_coordinates_tmp = s_ip_coordinates;
103  auto d_quadCoordinates_tmp = d_quadCoordinates;
104  if constexpr(Sacado::IsADType<typename EvalT::ScalarT>::value) {
105  d_quadCoordinates_tmp(i,j,k).val() = s_ip_coordinates_tmp(i,j,k);
106  }
107  else {
108  d_quadCoordinates_tmp(i,j,k) = s_ip_coordinates_tmp(i,j,k);
109  }
110  });
111 }
112 
113 // **********************************************************************
114 #endif
std::vector< int >::size_type getIntegrationRuleIndex(int ir_degree, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
Teuchos::RCP< PHX::DataLayout > dl_vector
Data layout for vector fields.
KOKKOS_FORCEINLINE_FUNCTION array_type get_static_view()
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)