Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_GatherTangents_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_TANGENTS_IMPL_HPP
44#define PANZER_GATHER_TANGENTS_IMPL_HPP
45
46#include "Teuchos_Assert.hpp"
47#include "Phalanx_DataLayout.hpp"
48
49#include "Intrepid2_Kernels.hpp"
50#include "Intrepid2_OrientationTools.hpp"
51
52#include "Panzer_PureBasis.hpp"
53#include "Kokkos_ViewFactory.hpp"
54
55#include "Teuchos_FancyOStream.hpp"
56
57template<typename EvalT,typename Traits>
60 const Teuchos::ParameterList& p)
61{
62 dof_name_ = (p.get< std::string >("DOF Name"));
63
64 if(p.isType< Teuchos::RCP<PureBasis> >("Basis"))
65 basis_ = p.get< Teuchos::RCP<PureBasis> >("Basis");
66 else
67 basis_ = p.get< Teuchos::RCP<const PureBasis> >("Basis");
68
69 pointRule_ = p.get<Teuchos::RCP<const PointRule> >("Point Rule");
70
71 Teuchos::RCP<PHX::DataLayout> basis_layout = basis_->functional;
72 Teuchos::RCP<PHX::DataLayout> vector_layout_vector = basis_->functional_grad;
73
74 // some sanity checks
75 TEUCHOS_ASSERT(basis_->isVectorBasis());
76
77 // setup the orientation field
78 std::string orientationFieldName = basis_->name() + " Orientation";
79 // setup all fields to be evaluated and constructed
80 pointValues_ = panzer::PointValues2<double> (pointRule_->getName()+"_",false);
81 pointValues_.setupArrays(pointRule_);
82
83 // the field manager will allocate all of these field
84 constJac_ = pointValues_.jac;
85 this->addDependentField(constJac_);
86
87 gatherFieldTangents_ = PHX::MDField<ScalarT,Cell,NODE,Dim>(dof_name_+"_Tangents",vector_layout_vector);
88 this->addEvaluatedField(gatherFieldTangents_);
89
90 this->setName("Gather Tangents");
91}
92
93// **********************************************************************
94template<typename EvalT,typename Traits>
98{
99 auto orientations = d.orientations_;
100 orientations_ = Kokkos::View<Intrepid2::Orientation*>("orientations_",orientations->size());
101 auto orientations_host = Kokkos::create_mirror_view(orientations_);
102 for (size_t i=0; i < orientations->size(); ++i)
103 orientations_host(i) = (*orientations)[i];
104 Kokkos::deep_copy(orientations_,orientations_host);
105
106 this->utils.setFieldData(pointValues_.jac,fm);
107 const shards::CellTopology & parentCell = *basis_->getCellTopology();
108 const int edgeDim = 1;
109 edgeParam_ = Intrepid2::RefSubcellParametrization<PHX::Device>::get(edgeDim, parentCell.getKey());
110
111 const int numEdges = gatherFieldTangents_.extent(1);
112 keys_ = Kokkos::View<unsigned int*>("parentCell.keys",numEdges);
113 auto keys_host = Kokkos::create_mirror_view(keys_);
114 for (int i=0; i < numEdges; ++i)
115 keys_host(i) = parentCell.getKey(edgeDim,i);
116 Kokkos::deep_copy(keys_,keys_host);
117}
118
119// **********************************************************************
120template<typename EvalT,typename Traits>
122evaluateFields(typename Traits::EvalData workset)
123{
124
125 if(workset.num_cells<=0)
126 return;
127 else {
128 const shards::CellTopology & parentCell = *basis_->getCellTopology();
129 const int cellDim = parentCell.getDimension();
130 const int numEdges = gatherFieldTangents_.extent(1);
131
132 auto workspace_tmp = Kokkos::createDynRankView(gatherFieldTangents_.get_static_view(),"workspace", workset.num_cells,4, cellDim);
133 const auto worksetJacobians = pointValues_.jac.get_view();
134 const auto cell_local_ids = workset.getLocalCellIDs();
135 auto gatherFieldTangents = gatherFieldTangents_.get_static_view();
136 auto keys = keys_;
137 auto orientations = orientations_;
138 auto edgeParam = edgeParam_;
139
140 // Loop over workset faces and edge points
141 Kokkos::parallel_for("panzer::GatherTangets",workset.num_cells,KOKKOS_LAMBDA(const int c){
142 int edgeOrts[12] = {};
143 orientations(cell_local_ids(c)).getEdgeOrientation(edgeOrts, numEdges);
144
145 auto workspace = Kokkos::subview(workspace_tmp,c,Kokkos::ALL(),Kokkos::ALL());
146 for(int pt = 0; pt < numEdges; pt++) {
147 auto phyEdgeTan = Kokkos::subview(gatherFieldTangents, c, pt, Kokkos::ALL());
148 auto ortEdgeTan = Kokkos::subview(workspace_tmp,c,0,Kokkos::ALL());
149
150 Intrepid2::Impl::OrientationTools::getRefSubcellTangents(
151 workspace,
152 edgeParam,
153 keys(pt),
154 pt,
155 edgeOrts[pt]);
156
157 auto J = Kokkos::subview(worksetJacobians, c, pt, Kokkos::ALL(), Kokkos::ALL());
158 Intrepid2::Kernels::Serial::matvec_product(phyEdgeTan, J, ortEdgeTan);
159
160 }// for pt
161 });// for pCell
162
163 }
164
165}
166
167#endif
void evaluateFields(typename Traits::EvalData d)
GatherTangents(const Teuchos::ParameterList &p)
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)
Kokkos::View< const int *, PHX::Device > getLocalCellIDs() const
Get the local cell IDs for the workset.
int num_cells
DEPRECATED - use: numCells()
Teuchos::RCP< const std::vector< Intrepid2::Orientation > > orientations_