Intrepid2
Intrepid2_ProjectionToolsDefHVOL.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 
49 #ifndef __INTREPID2_PROJECTIONTOOLSDEFHVOL_HPP__
50 #define __INTREPID2_PROJECTIONTOOLSDEFHVOL_HPP__
51 
53 #include "Intrepid2_ArrayTools.hpp"
55 
56 
57 namespace Intrepid2 {
58 namespace Experimental {
59 
60 template<typename DeviceType>
61 template<typename BasisType,
62 typename ortValueType, class ...ortProperties>
63 void
64 ProjectionTools<DeviceType>::getHVolEvaluationPoints(typename BasisType::ScalarViewType ePoints,
65  const Kokkos::DynRankView<ortValueType, ortProperties...> /*orts*/,
66  const BasisType* cellBasis,
68  const EvalPointsType ePointType) {
69  ordinal_type dim = cellBasis->getBaseCellTopology().getDimension();
70  auto refEPoints = Kokkos::create_mirror_view_and_copy(MemSpaceType(),projStruct->getEvalPoints(dim,0,ePointType));
71  auto ePointsRange = projStruct->getPointsRange(ePointType);
72  RealSpaceTools<DeviceType>::clone(Kokkos::subview(ePoints, Kokkos::ALL(), ePointsRange(dim, 0), Kokkos::ALL()), refEPoints);
73 }
74 
75 
76 template<typename DeviceType>
77 template<typename basisCoeffsValueType, class ...basisCoeffsProperties,
78 typename funValsValueType, class ...funValsProperties,
79 typename BasisType,
80 typename ortValueType,class ...ortProperties>
81 void
82 ProjectionTools<DeviceType>::getHVolBasisCoeffs(Kokkos::DynRankView<basisCoeffsValueType,basisCoeffsProperties...> basisCoeffs,
83  const Kokkos::DynRankView<funValsValueType,funValsProperties...> targetAtTargetEPoints,
84  const typename BasisType::ScalarViewType targetEPoints,
85  const Kokkos::DynRankView<ortValueType, ortProperties...> orts,
86  const BasisType* cellBasis,
88 
89  typedef typename BasisType::scalarType scalarType;
90  typedef Kokkos::DynRankView<scalarType,DeviceType> ScalarViewType;
91  ordinal_type dim = cellBasis->getBaseCellTopology().getDimension();
92 
93  ordinal_type basisCardinality = cellBasis->getCardinality();
94 
95  ordinal_type numCells = targetAtTargetEPoints.extent(0);
96 
97  auto refTargetEWeights = projStruct->getTargetEvalWeights(dim,0);
98  auto targetEPointsRange = projStruct->getTargetPointsRange();
99 
100  auto refBasisEWeights = projStruct->getBasisEvalWeights(dim,0);
101  auto basisEPointsRange = projStruct->getBasisPointsRange();
102 
103  ordinal_type numTargetEPoints = range_size(targetEPointsRange(dim,0));
104  ordinal_type numBasisEPoints = range_size(basisEPointsRange(dim,0));
105 
106  ScalarViewType basisAtBasisEPoints("basisAtBasisEPoints", 1, basisCardinality, numBasisEPoints);
107  ScalarViewType basisAtTargetEPoints("basisAtTargetEPoints", basisCardinality, numTargetEPoints);
108 
109  ScalarViewType basisEPoints("basisEPoints",numCells,projStruct->getNumBasisEvalPoints(), dim);
110  getHVolEvaluationPoints(basisEPoints, orts, cellBasis, projStruct, EvalPointsType::BASIS);
111 
112  cellBasis->getValues(Kokkos::subview(basisAtBasisEPoints, 0, Kokkos::ALL(), Kokkos::ALL()), Kokkos::subview(basisEPoints,0, Kokkos::ALL(), Kokkos::ALL()));
113  if(targetEPoints.rank()==3)
114  cellBasis->getValues(basisAtTargetEPoints, Kokkos::subview(targetEPoints, 0, Kokkos::ALL(), Kokkos::ALL()));
115  else
116  cellBasis->getValues(basisAtTargetEPoints, targetEPoints);
117 
118  ScalarViewType weightedBasisAtTargetEPoints("weightedBasisAtTargetEPoints_",numCells, basisCardinality, numTargetEPoints);
119  ScalarViewType weightedBasisAtBasisEPoints("weightedBasisAtBasisEPoints", 1, basisCardinality, numBasisEPoints);
120 
121  auto tagToOrdinal = Kokkos::create_mirror_view_and_copy(MemSpaceType(), cellBasis->getAllDofOrdinal());
122  auto cellDofs = Kokkos::subview(tagToOrdinal, dim, 0, Kokkos::ALL());
123 
124  ScalarViewType
125  massMat0("massMat0", 1, basisCardinality, basisCardinality),
126  massMat("massMat", numCells, basisCardinality, basisCardinality),
127  rhsMat("rhsMat", numCells, basisCardinality );
128 
129  ordinal_type offsetBasis = basisEPointsRange(dim,0).first;
130  ordinal_type offsetTarget = targetEPointsRange(dim,0).first;
131 
132  using HostSpaceType = typename Kokkos::Impl::is_space<DeviceType>::host_mirror_space::execution_space;
133  auto hWeightedBasisAtBasisEPoints = Kokkos::create_mirror_view(weightedBasisAtBasisEPoints);
134  auto hWeightedBasisAtTargetEPoints = Kokkos::create_mirror_view(weightedBasisAtTargetEPoints);
135  auto hBasisAtBasisEPoints = Kokkos::create_mirror_view_and_copy(HostSpaceType(), basisAtBasisEPoints);
136  auto hBasisAtTargetEPoints = Kokkos::create_mirror_view_and_copy(HostSpaceType(), basisAtTargetEPoints);
137 
138  for(ordinal_type j=0; j <basisCardinality; ++j) {
139  ordinal_type idof = cellBasis->getDofOrdinal(dim, 0, j);
140  for(ordinal_type iq=0; iq <ordinal_type(refBasisEWeights.extent(0)); ++iq)
141  hWeightedBasisAtBasisEPoints(0,j,iq) = hBasisAtBasisEPoints(0,idof,offsetBasis+iq) * refBasisEWeights(iq);
142  for(ordinal_type iq=0; iq <ordinal_type(refTargetEWeights.extent(0)); ++iq)
143  hWeightedBasisAtTargetEPoints(0,j,iq) = hBasisAtTargetEPoints(idof,offsetTarget+iq)* refTargetEWeights(iq);
144  }
145  Kokkos::deep_copy(weightedBasisAtBasisEPoints,hWeightedBasisAtBasisEPoints);
146  Kokkos::deep_copy(weightedBasisAtTargetEPoints,hWeightedBasisAtTargetEPoints);
147  FunctionSpaceTools<DeviceType >::integrate(massMat0, basisAtBasisEPoints, weightedBasisAtBasisEPoints);
148  RealSpaceTools<DeviceType>::clone(massMat, Kokkos::subview(massMat0,0,Kokkos::ALL(), Kokkos::ALL()));
149  RealSpaceTools<DeviceType>::clone(weightedBasisAtTargetEPoints, Kokkos::subview(weightedBasisAtTargetEPoints,0,Kokkos::ALL(), Kokkos::ALL()));
150  FunctionSpaceTools<DeviceType >::integrate(rhsMat, targetAtTargetEPoints, weightedBasisAtTargetEPoints);
151 
152  typedef Kokkos::DynRankView<scalarType, Kokkos::LayoutRight, DeviceType> WorkArrayViewType;
153  ScalarViewType t_("t",numCells, basisCardinality);
154  WorkArrayViewType w_("w",numCells,basisCardinality);
155 
156  ElemSystem cellSystem("cellSystem", true);
157  cellSystem.solve(basisCoeffs, massMat, rhsMat, t_, w_, cellDofs, basisCardinality);
158 }
159 
160 }
161 }
162 
163 #endif
164 
Header file for the abstract base class Intrepid2::DefaultCubatureFactory.
const range_tag getBasisPointsRange() const
Returns the range tag of the basis evaluation points subcells.
view_type getBasisEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the basis evaluation weights on a subcell.
Header file for the Intrepid2::FunctionSpaceTools class.
const range_tag getPointsRange(const EvalPointsType type) const
Returns the range tag of the basis/target evaluation points in subcells.
Class to solve a square system A x = b on each cell A is expected to be saddle a point (KKT) matrix o...
ordinal_type getNumBasisEvalPoints()
Returns number of basis evaluation points.
Header file for Intrepid2::ArrayTools class providing utilities for array operations.
static void getHVolBasisCoeffs(Kokkos::DynRankView< basisCoeffsValueType, basisCoeffsProperties... > basisCoeffs, const Kokkos::DynRankView< funValsValueType, funValsProperties... > targetAtEvalPoints, const typename BasisType::ScalarViewType evaluationPoints, const Kokkos::DynRankView< ortValueType, ortProperties... > cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct)
Computes the basis coefficients of the HVol projection of the target function.
static void integrate(Kokkos::DynRankView< outputValueValueType, outputValueProperties... > outputValues, const Kokkos::DynRankView< leftValueValueType, leftValueProperties... > leftValues, const Kokkos::DynRankView< rightValueValueType, rightValueProperties... > rightValues, const bool sumInto=false)
Contracts leftValues and rightValues arrays on the point and possibly space dimensions and stores the...
const range_tag getTargetPointsRange() const
Returns the range of the target function evaluation points on subcells.
void solve(ViewType1 basisCoeffs, ViewType2 elemMat, ViewType2 elemRhs, ViewType2 tau, ViewType3 w, const ViewType4 elemDof, ordinal_type n, ordinal_type m=0)
Solve the system and returns the basis coefficients solve the system either using Kokkos Kernel QR or...
An helper class to compute the evaluation points and weights needed for performing projections...
static void getHVolEvaluationPoints(typename BasisType::ScalarViewType evaluationPoints, const Kokkos::DynRankView< ortValueType, ortProperties... > cellOrientations, const BasisType *cellBasis, ProjectionStruct< DeviceType, typename BasisType::scalarType > *projStruct, const EvalPointsType evalPointType=EvalPointsType::TARGET)
Computes evaluation points for HVol projection.
static void clone(Kokkos::DynRankView< outputValueType, outputProperties... > output, const Kokkos::DynRankView< inputValueType, inputProperties... > input)
Clone input array.
view_type getEvalPoints(const ordinal_type subCellDim, const ordinal_type subCellId, EvalPointsType type) const
Returns the basis/target evaluation points on a subcell.
view_type getTargetEvalWeights(const ordinal_type subCellDim, const ordinal_type subCellId)
Returns the function evaluation weights on a subcell.