Panzer  Version of the Day
Panzer_CommonArrayFactories_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_COMMON_ARRAY_FACTORIES_IMPL_HPP
44 #define PANZER_COMMON_ARRAY_FACTORIES_IMPL_HPP
45 
46 #include "Teuchos_RCP.hpp"
47 
48 #include "Phalanx_DataLayout_MDALayout.hpp"
49 
50 #include "Phalanx_KokkosUtilities.hpp"
51 #include "Phalanx_KokkosViewFactory.hpp"
52 
53 namespace panzer {
54 
55 // Implementation for intrepid container factory
56 template <typename Scalar,typename T0>
57 Kokkos::DynRankView<Scalar,PHX::Device> Intrepid2FieldContainerFactory::
58 buildArray(const std::string & str,int d0) const
59 {
60  static_assert(std::is_same<Scalar,double>::value,"ERROR: CommonArryFactory for DynRankView only supports double scalar type!");
61  return Kokkos::DynRankView<Scalar,PHX::Device>(str,d0);
62 }
63 
64 template <typename Scalar,typename T0,typename T1>
65 Kokkos::DynRankView<Scalar,PHX::Device> Intrepid2FieldContainerFactory::
66 buildArray(const std::string & str,int d0,int d1) const
67 {
68  static_assert(std::is_same<Scalar,double>::value,"ERROR: CommonArryFactory for DynRankView only supports double scalar type!");
69  return Kokkos::DynRankView<Scalar,PHX::Device>(str,d0,d1);
70 }
71 
72 template <typename Scalar,typename T0,typename T1,typename T2>
73 Kokkos::DynRankView<Scalar,PHX::Device> Intrepid2FieldContainerFactory::
74 buildArray(const std::string & str,int d0,int d1,int d2) const
75 {
76  static_assert(std::is_same<Scalar,double>::value,"ERROR: CommonArryFactory for DynRankView only supports double scalar type!");
77  return Kokkos::DynRankView<Scalar,PHX::Device>(str,d0,d1,d2);
78 }
79 
80 template <typename Scalar,typename T0,typename T1,typename T2,typename T3>
81 Kokkos::DynRankView<Scalar,PHX::Device> Intrepid2FieldContainerFactory::
82 buildArray(const std::string & str,int d0,int d1,int d2,int d3) const
83 {
84  static_assert(std::is_same<Scalar,double>::value,"ERROR: CommonArryFactory for DynRankView only supports double scalar type!");
85  return Kokkos::DynRankView<Scalar,PHX::Device>(str,d0,d1,d2,d3);
86 }
87 
88 template <typename Scalar,typename T0,typename T1,typename T2,typename T3,typename T4>
89 Kokkos::DynRankView<Scalar,PHX::Device> Intrepid2FieldContainerFactory::
90 buildArray(const std::string & str,int d0,int d1,int d2,int d3,int d4) const
91 {
92  static_assert(std::is_same<Scalar,double>::value,"ERROR: CommonArryFactory for DynRankView only supports double scalar type!");
93  return Kokkos::DynRankView<Scalar,PHX::Device>(str,d0,d1,d2,d3,d4);
94 }
95 
96 // Implementation for MDField array factory
97 template <typename Scalar,typename T0>
98 PHX::MDField<Scalar> MDFieldArrayFactory::
99 buildArray(const std::string & str,int d0) const
100 {
101  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
102 
103  PHX::MDField<Scalar> field = PHX::MDField<Scalar>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0>(d0)));
104 
105  if(allocArray_)
106  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
107 
108  return field;
109 }
110 
111 template <typename Scalar,typename T0,typename T1>
112 PHX::MDField<Scalar> MDFieldArrayFactory::
113 buildArray(const std::string & str,int d0,int d1) const
114 {
115  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
116 
117  PHX::MDField<Scalar> field = PHX::MDField<Scalar>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1 >(d0,d1)));
118 
119  if(allocArray_)
120  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
121 
122  return field;
123 }
124 
125 template <typename Scalar,typename T0,typename T1,typename T2>
126 PHX::MDField<Scalar> MDFieldArrayFactory::
127 buildArray(const std::string & str,int d0,int d1,int d2) const
128 {
129  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
130 
131  PHX::MDField<Scalar> field = PHX::MDField<Scalar>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1,T2>(d0,d1,d2)));
132 
133  if(allocArray_)
134  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
135 
136  return field;
137 }
138 
139 template <typename Scalar,typename T0,typename T1,typename T2,typename T3>
140 PHX::MDField<Scalar> MDFieldArrayFactory::
141 buildArray(const std::string & str,int d0,int d1,int d2,int d3) const
142 {
143  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
144 
145  PHX::MDField<Scalar> field = PHX::MDField<Scalar>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1,T2,T3>(d0,d1,d2,d3)));
146 
147  if(allocArray_)
148  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
149 
150  return field;
151 }
152 
153 template <typename Scalar,typename T0,typename T1,typename T2,typename T3,typename T4>
154 PHX::MDField<Scalar> MDFieldArrayFactory::
155 buildArray(const std::string & str,int d0,int d1,int d2,int d3,int d4) const
156 {
157  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
158 
159  PHX::MDField<Scalar> field = PHX::MDField<Scalar>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1,T2,T3,T4>(d0,d1,d2,d3,d4)));
160 
161  if(allocArray_)
162  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
163 
164  return field;
165 }
166 
167 // Implementation for MDField array factory
168 template <typename Scalar,typename T0>
169 PHX::MDField<Scalar,T0> MDFieldArrayFactory::
170 buildStaticArray(const std::string & str,int d0) const
171 {
172  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
173 
174  PHX::MDField<Scalar,T0> field = PHX::MDField<Scalar,T0>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0>(d0)));
175 
176  if(allocArray_)
177  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
178 
179  return field;
180 }
181 
182 template <typename Scalar,typename T0,typename T1>
183 PHX::MDField<Scalar,T0,T1> MDFieldArrayFactory::
184 buildStaticArray(const std::string & str,int d0,int d1) const
185 {
186  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
187 
188  PHX::MDField<Scalar,T0,T1> field = PHX::MDField<Scalar,T0,T1>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1>(d0,d1)));
189 
190  if(allocArray_)
191  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
192 
193  return field;
194 }
195 
196 template <typename Scalar,typename T0,typename T1,typename T2>
197 PHX::MDField<Scalar,T0,T1,T2> MDFieldArrayFactory::
198 buildStaticArray(const std::string & str,int d0,int d1,int d2) const
199 {
200  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
201 
202  PHX::MDField<Scalar,T0,T1,T2> field = PHX::MDField<Scalar,T0,T1,T2>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1,T2>(d0,d1,d2)));
203 
204  if(allocArray_)
205  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
206 
207  return field;
208 }
209 
210 template <typename Scalar,typename T0,typename T1,typename T2,typename T3>
211 PHX::MDField<Scalar,T0,T1,T2,T3> MDFieldArrayFactory::
212 buildStaticArray(const std::string & str,int d0,int d1,int d2,int d3) const
213 {
214  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
215 
216  PHX::MDField<Scalar,T0,T1,T2,T3> field = PHX::MDField<Scalar,T0,T1,T2,T3>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1,T2,T3>(d0,d1,d2,d3)));
217 
218  if(allocArray_)
219  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
220 
221  return field;
222 }
223 
224 template <typename Scalar,typename T0,typename T1,typename T2,typename T3,typename T4>
225 PHX::MDField<Scalar,T0,T1,T2,T3,T4> MDFieldArrayFactory::
226 buildStaticArray(const std::string & str,int d0,int d1,int d2,int d3,int d4) const
227 {
228  typedef PHX::KokkosViewFactory<Scalar,PHX::Device> ViewFactory;
229 
230  PHX::MDField<Scalar,T0,T1,T2,T3,T4> field = PHX::MDField<Scalar,T0,T1,T2,T3,T4>(prefix_+str,Teuchos::rcp(new PHX::MDALayout<T0,T1,T2,T3,T4>(d0,d1,d2,d3,d4)));
231 
232  if(allocArray_)
233  field.setFieldData(ViewFactory::buildView(field.fieldTag(),ddims_));
234 
235  return field;
236 }
237 
238 } // end namespace panzer
239 
240 #endif
std::vector< PHX::index_size_type > ddims_
PHX::MDField< Scalar, T0 > buildStaticArray(const std::string &str, int d0) const
PHX::MDField< Scalar > buildArray(const std::string &str, int d0) const
PHX::MDField< const ScalarT, Cell, IP > field
Kokkos::DynRankView< Scalar, PHX::Device > buildArray(const std::string &str, int d0) const