Panzer  Version of the Day
Panzer_Sum_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_SUM_IMPL_HPP
44 #define PANZER_SUM_IMPL_HPP
45 
46 #include <cstddef>
47 #include <string>
48 #include <vector>
49 
50 namespace panzer {
51 
52 //**********************************************************************
53 template<typename EvalT, typename Traits>
56  const Teuchos::ParameterList& p)
57 {
58  std::string sum_name = p.get<std::string>("Sum Name");
59  Teuchos::RCP<std::vector<std::string> > value_names =
60  p.get<Teuchos::RCP<std::vector<std::string> > >("Values Names");
61  Teuchos::RCP<PHX::DataLayout> data_layout =
62  p.get< Teuchos::RCP<PHX::DataLayout> >("Data Layout");
63 
64  TEUCHOS_ASSERT(static_cast<int>(value_names->size()) < MAX_VALUES);
65 
66  // check if the user wants to scale each term independently
67  auto local_scalars = PHX::View<double *>("scalars",value_names->size());
68  auto local_scalars_host = Kokkos::create_mirror_view(local_scalars);
69  if(p.isType<Teuchos::RCP<const std::vector<double> > >("Scalars")) {
70  auto scalars_v = *p.get<Teuchos::RCP<const std::vector<double> > >("Scalars");
71 
72  // safety/sanity check
73  TEUCHOS_ASSERT(scalars_v.size()==value_names->size());
74 
75  for (std::size_t i=0; i < value_names->size(); ++i)
76  local_scalars_host(i) = scalars_v[i];
77  }
78  else {
79  for (std::size_t i=0; i < value_names->size(); ++i)
80  local_scalars_host(i) = 1.0;
81  }
82  Kokkos::deep_copy(local_scalars,local_scalars_host);
83 
84  scalars = local_scalars;
85 
86  sum = PHX::MDField<ScalarT>(sum_name, data_layout);
87 
88  this->addEvaluatedField(sum);
89 
90  for (std::size_t i=0; i < value_names->size(); ++i) {
91  values[i] = PHX::MDField<const ScalarT>( (*value_names)[i], data_layout);
92  this->addDependentField(values[i]);
93  }
94  /*
95  values.resize(value_names->size());
96  for (std::size_t i=0; i < value_names->size(); ++i) {
97  values[i] = PHX::MDField<const ScalarT>( (*value_names)[i], data_layout);
98  this->addDependentField(values[i]);
99  }
100  */
101 
102  std::string n = "Sum Evaluator";
103  this->setName(n);
104 }
105 
106 //**********************************************************************
107 template<typename EvalT, typename Traits>
108 void
111  typename Traits::SetupData /* worksets */,
112  PHX::FieldManager<Traits>& /* fm */)
113 {
114  cell_data_size = sum.size() / sum.fieldTag().dataLayout().extent(0);
115 }
116 
117 
118 //**********************************************************************
119 template<typename EvalT, typename TRAITS>
120 template<unsigned int RANK>
121 KOKKOS_INLINE_FUNCTION
123  auto num_vals = scalars.extent(0);
124 
125 
126  if (RANK == 1 )
127  {
128  for (std::size_t iv = 0; iv < num_vals; ++iv)
129  sum(i) += scalars(iv)*(values[iv](i));
130  }
131  else if (RANK == 2)
132  {
133  const size_t dim_1 = sum.extent(1);
134  for (std::size_t j = 0; j < dim_1; ++j)
135  for (std::size_t iv = 0; iv < num_vals; ++iv)
136  sum(i,j) += scalars(iv)*(values[iv](i,j));
137  }
138  else if (RANK == 3)
139  {
140  const size_t dim_1 = sum.extent(1),dim_2 = sum.extent(2);
141  for (std::size_t j = 0; j < dim_1; ++j)
142  for (std::size_t k = 0; k < dim_2; ++k)
143  for (std::size_t iv = 0; iv < num_vals; ++iv)
144  sum(i,j,k) += scalars(iv)*(values[iv](i,j,k));
145  }
146  else if (RANK == 4)
147  {
148  const size_t dim_1 = sum.extent(1),dim_2 = sum.extent(2),dim_3 = sum.extent(3);
149  for (std::size_t j = 0; j < dim_1; ++j)
150  for (std::size_t k = 0; k < dim_2; ++k)
151  for (std::size_t l = 0; l < dim_3; ++l)
152  for (std::size_t iv = 0; iv < num_vals; ++iv)
153  sum(i,j,k,l) += scalars(iv)*(values[iv](i,j,k,l));
154  }
155  else if (RANK == 5)
156  {
157  const size_t dim_1 = sum.extent(1),dim_2 = sum.extent(2),dim_3 = sum.extent(3),dim_4 = sum.extent(4);
158  for (std::size_t j = 0; j < dim_1; ++j)
159  for (std::size_t k = 0; k < dim_2; ++k)
160  for (std::size_t l = 0; l < dim_3; ++l)
161  for (std::size_t m = 0; m < dim_4; ++m)
162  for (std::size_t iv = 0; iv < num_vals; ++iv)
163  sum(i,j,k,l,m) += scalars(iv)*(values[iv](i,j,k,l,m));
164  }
165  else if (RANK == 6)
166  {
167  const size_t dim_1 = sum.extent(1),dim_2 = sum.extent(2),dim_3 = sum.extent(3),dim_4 = sum.extent(4),dim_5 = sum.extent(5);
168  for (std::size_t j = 0; j < dim_1; ++j)
169  for (std::size_t k = 0; k < dim_2; ++k)
170  for (std::size_t l = 0; l < dim_3; ++l)
171  for (std::size_t m = 0; m < dim_4; ++m)
172  for (std::size_t n = 0; n < dim_5; ++n)
173  for (std::size_t iv = 0; iv < num_vals; ++iv)
174  sum(i,j,k,l,m,n) += scalars(iv)*(values[iv](i,j,k,l,m,n));
175  }
176 }
177 
178 //**********************************************************************
179 template<typename EvalT, typename Traits>
180 void
183  typename Traits::EvalData /* workset */)
184 {
185 
186  sum.deep_copy(ScalarT(0.0));
187 
188  size_t rank = sum.rank();
189  const size_t length = sum.extent(0);
190  if (rank == 1 )
191  {
192  Kokkos::parallel_for(Kokkos::RangePolicy<PanzerSumTag<1> >(0, length), *this);
193  }
194  else if (rank == 2)
195  {
196  Kokkos::parallel_for(Kokkos::RangePolicy<PanzerSumTag<2> >(0, length), *this);
197  }
198  else if (rank == 3)
199  {
200  Kokkos::parallel_for(Kokkos::RangePolicy<PanzerSumTag<3> >(0, length), *this);
201  }
202  else if (rank == 4)
203  {
204  Kokkos::parallel_for(Kokkos::RangePolicy<PanzerSumTag<4> >(0, length), *this);
205  }
206  else if (rank == 5)
207  {
208  Kokkos::parallel_for(Kokkos::RangePolicy<PanzerSumTag<5> >(0, length), *this);
209  }
210  else if (rank == 6)
211  {
212  Kokkos::parallel_for(Kokkos::RangePolicy<PanzerSumTag<6> >(0, length), *this);
213  }
214  else
215  {
216  TEUCHOS_TEST_FOR_EXCEPTION(true, std::runtime_error, "ERROR: rank of sum is higher than supported");
217  }
218 }
219 
220 
221 //**********************************************************************
222 
223 template<typename EvalT, typename TRAITS,typename Tag0>
225 SumStatic(const Teuchos::ParameterList& p)
226 {
227  std::string sum_name = p.get<std::string>("Sum Name");
228  Teuchos::RCP<std::vector<std::string> > value_names =
229  p.get<Teuchos::RCP<std::vector<std::string> > >("Values Names");
230  Teuchos::RCP<PHX::DataLayout> data_layout =
231  p.get< Teuchos::RCP<PHX::DataLayout> >("Data Layout");
232 
233  // sanity check
234  TEUCHOS_ASSERT(data_layout->rank()==1);
235 
236  sum = PHX::MDField<ScalarT,Tag0>(sum_name, data_layout);
237 
238  this->addEvaluatedField(sum);
239 
240  values.resize(value_names->size());
241  for (std::size_t i=0; i < value_names->size(); ++i) {
242  values[i] = PHX::MDField<const ScalarT,Tag0>( (*value_names)[i], data_layout);
243  this->addDependentField(values[i]);
244  }
245 
246  std::string n = "SumStatic Rank 1 Evaluator";
247  this->setName(n);
248 }
249 
250 //**********************************************************************
251 
252 template<typename EvalT, typename TRAITS,typename Tag0>
254 evaluateFields(typename TRAITS::EvalData /* d */)
255 {
256  sum.deep_copy(ScalarT(0.0));
257  for (std::size_t i = 0; i < sum.extent(0); ++i)
258  for (std::size_t d = 0; d < values.size(); ++d)
259  sum(i) += (values[d])(i);
260 }
261 
262 //**********************************************************************
263 //**********************************************************************
264 
265 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1>
267 SumStatic(const Teuchos::ParameterList& p)
268 {
269  std::string sum_name = p.get<std::string>("Sum Name");
270  Teuchos::RCP<std::vector<std::string> > value_names =
271  p.get<Teuchos::RCP<std::vector<std::string> > >("Values Names");
272  Teuchos::RCP<PHX::DataLayout> data_layout =
273  p.get< Teuchos::RCP<PHX::DataLayout> >("Data Layout");
274 
275  // check if the user wants to scale each term independently
276  if(p.isType<Teuchos::RCP<const std::vector<double> > >("Scalars")) {
277  Teuchos::RCP<const std::vector<double> > scalar_values
278  = p.get<Teuchos::RCP<const std::vector<double> > >("Scalars");
279 
280  // safety/sanity check
281  TEUCHOS_ASSERT(scalar_values->size()==value_names->size());
282  useScalars = true;
283 
284  PHX::View<double*> scalars_nc = PHX::View<double*>("scalars",scalar_values->size());
285 
286  for(std::size_t i=0;i<scalar_values->size();i++)
287  scalars_nc(i) = (*scalar_values)[i];
288 
289  scalars = scalars_nc;
290  }
291  else {
292  useScalars = false;
293  }
294 
295  // sanity check
296  TEUCHOS_ASSERT(data_layout->rank()==2);
297 
298  sum = PHX::MDField<ScalarT,Tag0,Tag1>(sum_name, data_layout);
299 
300  this->addEvaluatedField(sum);
301 
302  values.resize(value_names->size());
303  for (std::size_t i=0; i < value_names->size(); ++i) {
304  values[i] = PHX::MDField<const ScalarT,Tag0,Tag1>( (*value_names)[i], data_layout);
305  this->addDependentField(values[i]);
306  }
307  numValues = value_names->size();
308  TEUCHOS_ASSERT(numValues<=MAX_VALUES);
309 
310  std::string n = "SumStatic Rank 2 Evaluator";
311  this->setName(n);
312 }
313 
314 //**********************************************************************
315 
316 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1>
318 SumStatic(const std::vector<PHX::Tag<typename EvalT::ScalarT>> & inputs,
319  const std::vector<double> & scalar_values,
320  const PHX::Tag<typename EvalT::ScalarT> & output)
321 {
322  TEUCHOS_ASSERT(scalar_values.size()==inputs.size());
323 
324  // check if the user wants to scale each term independently
325  if(scalars.size()==0) {
326  useScalars = false;
327  }
328  else {
329  useScalars = true;
330 
331  PHX::View<double*> scalars_nc = PHX::View<double*>("scalars",scalar_values.size());
332 
333  for(std::size_t i=0;i<scalar_values.size();i++)
334  scalars_nc(i) = scalar_values[i];
335 
336  scalars = scalars_nc;
337  }
338 
339  // sanity check
340  TEUCHOS_ASSERT(inputs.size()<=MAX_VALUES);
341 
342  sum = output;
343  this->addEvaluatedField(sum);
344 
345  values.resize(inputs.size());
346  for (std::size_t i=0; i < inputs.size(); ++i) {
347  values[i] = inputs[i];
348  this->addDependentField(values[i]);
349  }
350 
351  numValues = inputs.size();
352 
353  std::string n = "SumStatic Rank 2 Evaluator";
354  this->setName(n);
355 }
356 
357 //**********************************************************************
358 
359 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1>
361 postRegistrationSetup(typename TRAITS::SetupData /* d */,
362  PHX::FieldManager<TRAITS>& /* fm */)
363 {
364  for (std::size_t i=0; i < values.size(); ++i)
365  value_views[i] = values[i].get_static_view();
366 }
367 
368 //**********************************************************************
369 
370 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1>
372 evaluateFields(typename TRAITS::EvalData /* d */)
373 {
374  sum.deep_copy(ScalarT(0.0));
375 
376  // Kokkos::parallel_for(sum.extent(0), *this);
377  if(useScalars)
378  Kokkos::parallel_for(Kokkos::RangePolicy<PHX::Device,ScalarsTag>(0,sum.extent(0)), *this);
379  else
380  Kokkos::parallel_for(Kokkos::RangePolicy<PHX::Device,NoScalarsTag>(0,sum.extent(0)), *this);
381 }
382 
383 //**********************************************************************
384 
385 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1>
386 KOKKOS_INLINE_FUNCTION
388 operator()(const ScalarsTag, const unsigned c ) const
389 {
390  for (int i=0;i<numValues;i++) {
391  for (int j = 0; j < sum.extent_int(1); ++j)
392  sum(c,j) += scalars(i)*value_views[i](c,j);
393  }
394 }
395 
396 //**********************************************************************
397 
398 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1>
399 KOKKOS_INLINE_FUNCTION
401 operator()(const NoScalarsTag, const unsigned c ) const
402 {
403  for (int i=0;i<numValues;i++) {
404  for (int j = 0; j < sum.extent_int(1); ++j)
405  sum(c,j) += value_views[i](c,j);
406  }
407 }
408 
409 
410 //**********************************************************************
411 //**********************************************************************
412 
413 /*
414 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1,typename Tag2>
415 SumStatic<EvalT,TRAITS,Tag0,Tag1,Tag2>::
416 SumStatic(const Teuchos::ParameterList& p)
417 {
418  std::string sum_name = p.get<std::string>("Sum Name");
419  Teuchos::RCP<std::vector<std::string> > value_names =
420  p.get<Teuchos::RCP<std::vector<std::string> > >("Values Names");
421  Teuchos::RCP<PHX::DataLayout> data_layout =
422  p.get< Teuchos::RCP<PHX::DataLayout> >("Data Layout");
423 
424  // sanity check
425  TEUCHOS_ASSERT(data_layout->rank()==3);
426 
427  sum = PHX::MDField<ScalarT,Tag0,Tag1,Tag2>(sum_name, data_layout);
428 
429  this->addEvaluatedField(sum);
430 
431  values.resize(value_names->size());
432  for (std::size_t i=0; i < value_names->size(); ++i) {
433  values[i] = PHX::MDField<ScalarT,Tag0,Tag1,Tag2>( (*value_names)[i], data_layout);
434  this->addDependentField(values[i]);
435  }
436 
437  std::string n = "Sum Evaluator";
438  this->setName(n);
439 }
440 */
441 
442 //**********************************************************************
443 
444 /*
445 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1,typename Tag2>
446 void SumStatic<EvalT,TRAITS,Tag0,Tag1,Tag2>::
447 evaluateFields(typename TRAITS::EvalData d)
448 {
449  sum.deep_copy(ScalarT(0.0));
450 
451  for (std::size_t d = 0; d < values.size(); ++d)
452  for (std::size_t i = 0; i < sum.extent(0); ++i)
453  for (std::size_t j = 0; j < sum.extent(1); ++j)
454  for (std::size_t k = 0; k < sum.extent(2); ++k)
455  sum(i,j,k) += (values[d])(i);
456 }
457 */
458 
459 //**********************************************************************
460 //**********************************************************************
461 
462 template<typename EvalT, typename TRAITS,typename Tag0,typename Tag1,typename Tag2>
463 Teuchos::RCP<PHX::Evaluator<TRAITS> >
464 buildStaticSumEvaluator(const std::string & sum_name,
465  const std::vector<std::string> & value_names,
466  const Teuchos::RCP<PHX::DataLayout> & data_layout)
467 {
468  Teuchos::ParameterList p;
469  p.set<std::string>("Sum Name",sum_name);
470  p.set<Teuchos::RCP<std::vector<std::string> > >("Values Names",Teuchos::rcp(new std::vector<std::string>(value_names)));
471  p.set< Teuchos::RCP<PHX::DataLayout> >("Data Layout",data_layout);
472 
473  return Teuchos::rcp(new SumStatic<EvalT,TRAITS,Tag0,Tag1,Tag2>(p));
474 }
475 
476 }
477 
478 #endif
Teuchos::RCP< PHX::Evaluator< TRAITS > > buildStaticSumEvaluator(const std::string &sum_name, const std::vector< std::string > &value_names, const Teuchos::RCP< PHX::DataLayout > &data_layout)
Sum(const Teuchos::ParameterList &p)
EvalT::ScalarT ScalarT
Definition: Panzer_Sum.hpp:120
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
SumStatic(const Teuchos::ParameterList &p)
KOKKOS_INLINE_FUNCTION void operator()(PanzerSumTag< RANK >, const int &i) const
void evaluateFields(typename Traits::EvalData d)
typename EvalT::ScalarT ScalarT
Definition: Panzer_Sum.hpp:88
void evaluateFields(typename TRAITS::EvalData d)