Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_DOF_BasisToBasis_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_BASIS_TO_BASIS_IMPL_HPP
44#define PANZER_DOF_BASIS_TO_BASIS_IMPL_HPP
45
48#include "Intrepid2_FunctionSpaceTools.hpp"
49#include "Phalanx_DataLayout_MDALayout.hpp"
50#include "Intrepid2_Basis.hpp"
51#include "Teuchos_Assert.hpp"
52
53namespace panzer {
54
55//**********************************************************************
56template <typename EvalT, typename TRAITST>
58DOF_BasisToBasis(const std::string & fieldName,
59 const PureBasis & sourceBasis,
60 const PureBasis & targetBasis)
61{
62 TEUCHOS_ASSERT(sourceBasis.numCells() == targetBasis.numCells());
63
64 // **************
65 // Declare fields
66 // **************
67 dof_source_coeff = PHX::MDField<const ScalarT>(fieldName,sourceBasis.functional);
68 dof_target_coeff = PHX::MDField<ScalarT>(fieldName,targetBasis.functional);
69
70 this->addDependentField(dof_source_coeff);
71 this->addEvaluatedField(dof_target_coeff);
72
73 // **************
74 // Get coordinate points for reference cell on target basis
75 // **************
76 Kokkos::DynRankView<double,PHX::Device>intrpCoords =
77 Kokkos::DynRankView<double,PHX::Device>("intrpCoords",targetBasis.cardinality(),targetBasis.dimension());
78
79 targetBasis.getIntrepid2Basis<PHX::exec_space,double,double>()->getDofCoords(intrpCoords);
80
81 // **************
82 // Evaluate source basis values at target basis coordinates
83 // **************
84 Kokkos::DynRankView<double,PHX::Device> basisRef =
85 Kokkos::DynRankView<double,PHX::Device>("basisRef",sourceBasis.cardinality(),targetBasis.cardinality());
86
87 sourceBasis.getIntrepid2Basis()->getValues(basisRef, intrpCoords, Intrepid2::OPERATOR_VALUE);
88
89 // **************
90 // Copy the reference basis values for all cells in workset
91 // **************
92 basis = Kokkos::DynRankView<double,PHX::Device>("basis",sourceBasis.numCells(),sourceBasis.cardinality(),targetBasis.cardinality());
93 Intrepid2::FunctionSpaceTools<PHX::exec_space>::HGRADtransformVALUE(basis,basisRef);
94
95 std::string n = "DOF_BasisToBasis: " + dof_target_coeff.fieldTag().name();
96 this->setName(n);
97}
98
99//**********************************************************************
100template <typename EvalT, typename TRAITST>
101void DOF_BasisToBasis<EvalT,TRAITST>::evaluateFields(typename TRAITST::EvalData workset)
102{
103 // Zero out arrays (intrepid does a sum!)
104 dof_target_coeff.deep_copy(ScalarT(0.0));
105
106 if(workset.num_cells>0) {
107
108 // evaluate function at specified points
109 Intrepid2::FunctionSpaceTools<PHX::exec_space>::evaluate(dof_target_coeff.get_view(),
110 dof_source_coeff.get_view(),
111 basis);
112 }
113}
114
115//**********************************************************************
116
117}
118
119#endif
void evaluateFields(typename TRAITST::EvalData workset)
DOF_BasisToBasis(const std::string &fieldName, const PureBasis &sourceBasis, const PureBasis &targetBasis)
Ctor.
Description and data layouts associated with a particular basis.
int numCells() const
Returns the number of cells in the data layouts.
Teuchos::RCP< Intrepid2::Basis< PHX::Device::execution_space, double, double > > getIntrepid2Basis() const
int cardinality() const
Returns the number of basis coefficients.
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis> or <Cell,Basis>
int dimension() const
Returns the dimension of the basis from the topology.