Intrepid2
Intrepid2_HVOL_TET_Cn_FEM.hpp
Go to the documentation of this file.
1// @HEADER
2// ************************************************************************
3//
4// Intrepid2 Package
5// Copyright (2007) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38// Mauro Perego (mperego@sandia.gov)
39//
40// ************************************************************************
41// @HEADER
42
48#ifndef __INTREPID2_HVOL_TET_CN_FEM_HPP__
49#define __INTREPID2_HVOL_TET_CN_FEM_HPP__
50
51#include "Intrepid2_Basis.hpp"
52
54#include "Teuchos_LAPACK.hpp"
55
56
57namespace Intrepid2 {
58
74 namespace Impl {
75
80 public:
81 typedef struct Tetrahedron<4> cell_topology_type;
85 template<EOperator opType>
86 struct Serial {
87 template<typename outputValueViewType,
88 typename inputPointViewType,
89 typename workViewType,
90 typename vinvViewType>
92 static void
96 const vinvViewType vinv );
97
98
100 static ordinal_type
101 getWorkSizePerPoint(ordinal_type order) {
103 switch (opType) {
104 case OPERATOR_GRAD:
105 case OPERATOR_CURL:
106 case OPERATOR_D1:
107 return 7*cardinality;
108 default:
110 }
111 }
112 };
113
114 template<typename DeviceType, ordinal_type numPtsPerEval,
115 typename outputValueValueType, class ...outputValueProperties,
116 typename inputPointValueType, class ...inputPointProperties,
117 typename vinvValueType, class ...vinvProperties>
118 static void
119 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
120 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
121 const Kokkos::DynRankView<vinvValueType, vinvProperties...> vinv,
122 const EOperator operatorType);
123
127 template<typename outputValueViewType,
128 typename inputPointViewType,
129 typename vinvViewType,
130 typename workViewType,
131 EOperator opType,
132 ordinal_type numPtsEval>
133 struct Functor {
134 outputValueViewType _outputValues;
135 const inputPointViewType _inputPoints;
136 const vinvViewType _vinv;
137 workViewType _work;
138
141 inputPointViewType inputPoints_,
142 vinvViewType vinv_,
144 : _outputValues(outputValues_), _inputPoints(inputPoints_),
145 _vinv(vinv_), _work(work_) {}
146
148 void operator()(const size_type iter) const {
149 const auto ptBegin = Util<ordinal_type>::min(iter*numPtsEval, _inputPoints.extent(0));
150 const auto ptEnd = Util<ordinal_type>::min(ptBegin+numPtsEval, _inputPoints.extent(0));
151
152 const auto ptRange = Kokkos::pair<ordinal_type,ordinal_type>(ptBegin, ptEnd);
153 const auto input = Kokkos::subview( _inputPoints, ptRange, Kokkos::ALL() );
154
155 typename workViewType::pointer_type ptr = _work.data() + _work.extent(0)*ptBegin*get_dimension_scalar(_work);
156
157 auto vcprop = Kokkos::common_view_alloc_prop(_work);
158 workViewType work(Kokkos::view_wrap(ptr,vcprop), (ptEnd-ptBegin)*_work.extent(0));
159
160 switch (opType) {
161 case OPERATOR_VALUE : {
162 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange );
163 Serial<opType>::getValues( output, input, work, _vinv );
164 break;
165 }
166 case OPERATOR_GRAD :
167 case OPERATOR_D1 :
168 case OPERATOR_D2 :
169 //case OPERATOR_D3 :
170 {
171 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange, Kokkos::ALL() );
172 Serial<opType>::getValues( output, input, work, _vinv );
173 break;
174 }
175 default: {
176 INTREPID2_TEST_FOR_ABORT( true,
177 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::Functor) operator is not supported");
178
179 }
180 }
181 }
182 };
183 };
184 }
185
186 template<typename DeviceType = void,
187 typename outputValueType = double,
188 typename pointValueType = double>
190 : public Basis<DeviceType,outputValueType,pointValueType> {
191 public:
195
198 Basis_HVOL_TET_Cn_FEM(const ordinal_type order,
199 const EPointType pointType = POINTTYPE_EQUISPACED);
200
201
205
207
208 using Basis<DeviceType,outputValueType,pointValueType>::getValues;
209
210 virtual
211 void
212 getValues( OutputViewType outputValues,
213 const PointViewType inputPoints,
214 const EOperator operatorType = OPERATOR_VALUE) const override {
215#ifdef HAVE_INTREPID2_DEBUG
217 inputPoints,
218 operatorType,
219 this->getBaseCellTopology(),
220 this->getCardinality() );
221#endif
222 constexpr ordinal_type numPtsPerEval = Parameters::MaxNumPtsPerBasisEval;
223 Impl::Basis_HVOL_TET_Cn_FEM::
224 getValues<DeviceType,numPtsPerEval>( outputValues,
225 inputPoints,
226 this->vinv_,
227 operatorType);
228 }
229
230 virtual
231 void
232 getDofCoords( ScalarViewType dofCoords ) const override {
233#ifdef HAVE_INTREPID2_DEBUG
234 // Verify rank of output array.
235 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.rank() != 2, std::invalid_argument,
236 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getDofCoords) rank = 2 required for dofCoords array");
237 // Verify 0th dimension of output array.
238 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->getCardinality(), std::invalid_argument,
239 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
240 // Verify 1st dimension of output array.
241 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->getBaseCellTopology().getDimension(), std::invalid_argument,
242 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
243#endif
244 Kokkos::deep_copy(dofCoords, this->dofCoords_);
245 }
246
247 virtual
248 void
249 getDofCoeffs( ScalarViewType dofCoeffs ) const override {
250#ifdef HAVE_INTREPID2_DEBUG
251 // Verify rank of output array.
252 INTREPID2_TEST_FOR_EXCEPTION( dofCoeffs.rank() != 1, std::invalid_argument,
253 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getdofCoeffs) rank = 1 required for dofCoeffs array");
254 // Verify 0th dimension of output array.
255 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->getCardinality(), std::invalid_argument,
256 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getdofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
257#endif
258 Kokkos::deep_copy(dofCoeffs, 1.0);
259 }
260
261 void
262 getVandermondeInverse( ScalarViewType vinv ) const {
263 // has to be same rank and dimensions
264 Kokkos::deep_copy(vinv, this->vinv_);
265 }
266
267 virtual
268 const char*
269 getName() const override {
270 return "Intrepid2_HVOL_TET_Cn_FEM";
271 }
272
273 virtual
274 bool
275 requireOrientation() const override {
276 return false;
277 }
278
283
284 private:
285
288 Kokkos::DynRankView<scalarType,DeviceType> vinv_;
289 EPointType pointType_;
290
291 };
292
293}// namespace Intrepid2
294
296
297#endif
Header file for the abstract base class Intrepid2::Basis.
BasisPtr< typename Kokkos::HostSpace::device_type, OutputType, PointType > HostBasisPtr
Pointer to a Basis whose device type is on the host (Kokkos::HostSpace::device_type),...
void getValues_HVOL_Args(const outputValueViewType outputValues, const inputPointViewType inputPoints, const EOperator operatorType, const shards::CellTopology cellTopo, const ordinal_type basisCard)
Runtime check of the arguments for the getValues method in an HVOL-conforming FEM basis....
Definition file for FEM basis functions of degree n for H(vol) functions on TET.
Header file for Intrepid2::PointTools class to provide utilities for barycentric coordinates,...
Implementation of the default HVOL-compatible Lagrange basis of arbitrary degree on Tetrahedron cell.
virtual const char * getName() const override
Returns basis name.
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Coefficients for computing degrees of freedom for Lagrangian basis If P is an element of the space sp...
virtual bool requireOrientation() const override
True if orientation is required.
Basis_HVOL_TET_Cn_FEM(const ordinal_type order, const EPointType pointType=POINTTYPE_EQUISPACED)
Constructor.
virtual HostBasisPtr< outputValueType, pointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
Kokkos::DynRankView< scalarType, DeviceType > vinv_
inverse of Generalized Vandermonde matrix, whose columns store the expansion coefficients of the noda...
virtual void getDofCoords(ScalarViewType dofCoords) const override
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
virtual void getValues(OutputViewType outputValues, const PointViewType inputPoints, const EOperator operatorType=OPERATOR_VALUE) const override
Evaluation of a FEM basis on a reference cell.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
Kokkos::View< ordinal_type ***, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray3DHost
View type for 3d host array.
ordinal_type basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
ordinal_type getCardinality() const
Returns cardinality of the basis.
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
Device DeviceType
(Kokkos) Device type on which Basis is templated. Does not necessarily return true for Kokkos::is_dev...
Kokkos::DynRankView< scalarType, DeviceType > dofCoords_
Coordinates of degrees-of-freedom for basis functions defined in physical space.
Kokkos::View< ordinal_type **, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray2DHost
View type for 2d host array.
ScalarTraits< pointValueType >::scalar_type scalarType
Scalar type for point values.
shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation https://trili...
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
See Intrepid2::Basis_HVOL_TET_Cn_FEM.
static constexpr ordinal_type MaxNumPtsPerBasisEval
The maximum number of points to eval in serial mode.
small utility functions