Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_STK_ScatterCellAvgQuantity_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_STK_SCATTER_CELL_AVG_QUANTITY_IMPL_HPP
44#define PANZER_STK_SCATTER_CELL_AVG_QUANTITY_IMPL_HPP
45
46#include "Teuchos_Assert.hpp"
47
48#include "Phalanx_config.hpp"
49#include "Phalanx_Evaluator_Macros.hpp"
50#include "Phalanx_MDField.hpp"
51#include "Phalanx_DataLayout.hpp"
52#include "Phalanx_DataLayout_MDALayout.hpp"
53
56
57#include "Teuchos_FancyOStream.hpp"
58#include "Teuchos_ArrayRCP.hpp"
59
60namespace panzer_stk {
61
62template<typename EvalT, typename Traits>
65 const Teuchos::ParameterList& p) :
66 mesh_(p.get<Teuchos::RCP<STK_Interface> >("Mesh"))
67{
68 using panzer::Cell;
69 using panzer::Point;
70
71 std::string scatterName = p.get<std::string>("Scatter Name");
72
73 // if it's there pull the ouptut scaling
74 if (p.isParameter("Variable Scale Factors Map"))
75 varScaleFactors_ = p.get<Teuchos::RCP<std::map<std::string,double>>>("Variable Scale Factors Map");
76
77 const std::vector<std::string> & names =
78 *(p.get< Teuchos::RCP< std::vector<std::string> > >("Field Names"));
79
80 Teuchos::RCP<panzer::IntegrationRule> intRule =
81 p.get< Teuchos::RCP<panzer::IntegrationRule> >("IR");
82
83 // build dependent fields
84 scatterFields_.resize(names.size());
85 stkFields_.resize(names.size());
86 for (std::size_t fd = 0; fd < names.size(); ++fd) {
87 scatterFields_[fd] =
88 PHX::MDField<const ScalarT,Cell,Point>(names[fd],intRule->dl_scalar);
89 this->addDependentField(scatterFields_[fd]);
90 }
91
92 // setup a dummy field to evaluate
93 PHX::Tag<ScalarT> scatterHolder(scatterName,Teuchos::rcp(new PHX::MDALayout<panzer::Dummy>(0)));
94 this->addEvaluatedField(scatterHolder);
95
96 this->setName(scatterName+": STK-Scatter Cell Fields");
97}
98
99template<typename EvalT, typename Traits>
100void
103 typename Traits::SetupData /* d */,
105{
106 for (std::size_t fd = 0; fd < scatterFields_.size(); ++fd) {
107 std::string fieldName = scatterFields_[fd].fieldTag().name();
108
109 stkFields_[fd] = mesh_->getMetaData()->get_field<VariableField>(stk::topology::ELEMENT_RANK, fieldName);
110 }
111}
112
113template<typename EvalT, typename Traits>
114void
117 typename Traits::EvalData workset)
118{
119 panzer::MDFieldArrayFactory af("",true);
120
121 // for convenience pull out some objects from workset
122 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
123 std::string blockId = this->wda(workset).block_id;
124
125 for(std::size_t fieldIndex=0; fieldIndex<scatterFields_.size();fieldIndex++) {
126 auto field = scatterFields_[fieldIndex].get_static_view();
127 auto average = af.buildStaticArray<double,panzer::Cell,panzer::NODE>("",field.extent(0),1).get_static_view();
128 // write to double field
129 Kokkos::parallel_for("ScatterCellAvgQuantity",field.extent(0), KOKKOS_LAMBDA(int i) {
130 for(unsigned j=0; j<field.extent(1);j++)
131 average(i,0) += Sacado::scalarValue(field(i,j));
132 average(i,0) /= field.extent(1);
133 });
134 Kokkos::fence();
135
136 std::string varname = scatterFields_[fieldIndex].fieldTag().name();
137 double scalef = 1.0;
138
139 if (!varScaleFactors_.is_null())
140 {
141 std::map<std::string,double> *tmp_sfs = varScaleFactors_.get();
142 if(tmp_sfs->find(varname) != tmp_sfs->end())
143 scalef = (*tmp_sfs)[varname];
144 }
145
146 mesh_->setCellFieldData(varname,blockId,localCellIds,average,scalef);
147 }
148}
149
150} // end panzer_stk
151
152#endif
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we'll contribute, or in which we'll store, the result of computing this integral.
PHX::MDField< Scalar, T0 > buildStaticArray(const std::string &str, int d0) const
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
panzer_stk::STK_Interface::SolutionFieldType VariableField
Teuchos::RCP< std::map< std::string, double > > varScaleFactors_
std::vector< PHX::MDField< const ScalarT, panzer::Cell, panzer::Point > > scatterFields_