43#ifndef PANZER_GLOBAL_STATISTICS_IMPL_HPP
44#define PANZER_GLOBAL_STATISTICS_IMPL_HPP
46#include "Intrepid2_FunctionSpaceTools.hpp"
52#include "Phalanx_DataLayout_MDALayout.hpp"
53#include "Teuchos_ScalarTraits.hpp"
54#include "Teuchos_CommHelpers.hpp"
60template<
typename EvalT,
typename Traits>
63 const Teuchos::ParameterList& p)
65 comm = p.get< Teuchos::RCP<const Teuchos::Comm<int> > >(
"Comm");
67 global_data = p.get<Teuchos::RCP<panzer::GlobalData> >(
"Global Data");
72 std::string names_string = p.get<std::string>(
"Names");
73 std::vector<std::string> names;
76 Teuchos::RCP<panzer::IntegrationRule> ir = p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR");
79 for (
typename std::vector<std::string>::const_iterator name = names.begin(); name != names.end(); ++name)
80 field_values.push_back(PHX::MDField<const ScalarT,Cell,IP>(*name, ir->dl_scalar));
82 Teuchos::RCP<PHX::MDALayout<Cell> > cell_dl = Teuchos::rcp(
new PHX::MDALayout<Cell>(ir->dl_scalar->extent(0)));
83 volumes = PHX::MDField<ScalarT,Cell>(
"Cell Volumes",cell_dl);
85 tmp = PHX::MDField<ScalarT,Cell>(
"GlobalStatistics:tmp:"+names_string,cell_dl);
86 ones = PHX::MDField<ScalarT,Cell,IP>(
"GlobalStatistics:ones:"+names_string,ir->dl_scalar);
88 this->addEvaluatedField(volumes);
89 this->addEvaluatedField(tmp);
90 this->addEvaluatedField(ones);
91 for (
typename std::vector<PHX::MDField<const ScalarT,Cell,IP> >::const_iterator
field = field_values.begin();
93 this->addDependentField(*
field);
96 averages.resize(field_values.size());
97 maxs.resize(field_values.size());
98 mins.resize(field_values.size());
99 global_maxs.resize(field_values.size());
100 global_mins.resize(field_values.size());
101 global_averages.resize(field_values.size());
103 ir_order = ir->cubature_degree;
105 std::string n =
"GlobalStatistics: " + names_string;
110template<
typename EvalT,
typename Traits>
118 auto l_ones = ones.get_static_view();
119 Kokkos::parallel_for(
"GlobalStatistics", l_ones.extent(0), KOKKOS_LAMBDA(
int cell) {
120 for (std::size_t ip = 0; ip < l_ones.extent(1); ++ip)
121 l_ones(cell,ip) = 1.0;
126template<
typename EvalT,
typename Traits>
135 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::integrate(volumes.get_view(),
137 (this->wda(workset).int_rules[ir_index])->weighted_measure.get_view());
138 auto volumes_h = Kokkos::create_mirror_view(as_view(volumes));
139 Kokkos::deep_copy(volumes_h, as_view(volumes));
141 for (index_t cell = 0; cell < workset.
num_cells; ++cell)
142 total_volume += volumes_h(cell);
144 typename std::vector<PHX::MDField<ScalarT,Cell,IP> >::size_type field_index = 0;
145 for (
typename std::vector<PHX::MDField<const ScalarT,Cell,IP> >::iterator
field = field_values.begin();
146 field != field_values.end(); ++
field,++field_index) {
148 Intrepid2::FunctionSpaceTools<PHX::Device::execution_space>::integrate(tmp.get_view(),
150 (this->wda(workset).int_rules[ir_index])->weighted_measure.get_view());
151 auto tmp_h = Kokkos::create_mirror_view(tmp.get_static_view());
152 auto field_h = Kokkos::create_mirror_view(
field->get_static_view());
153 Kokkos::deep_copy(tmp_h, tmp.get_static_view());
154 Kokkos::deep_copy(field_h,
field->get_static_view());
157 for (index_t cell = 0; cell < workset.
num_cells; ++cell) {
158 averages[field_index] += tmp_h(cell);
160 for (
typename PHX::MDField<ScalarT,Cell,IP>::size_type ip = 0; ip < (
field->extent(1)); ++ip) {
161 maxs[field_index] = std::max( field_h(cell,ip), maxs[field_index]);
162 mins[field_index] = std::min( field_h(cell,ip), mins[field_index]);
170template<
typename EvalT,
typename Traits>
176 total_volume = Teuchos::ScalarTraits<ScalarT>::zero();
178 for (
typename std::vector<ScalarT>::iterator
field = averages.begin();
field != averages.end(); ++
field)
179 *
field = Teuchos::ScalarTraits<ScalarT>::zero();
181 for (
typename std::vector<ScalarT>::iterator
field = maxs.begin();
field != maxs.end(); ++
field)
182 *
field = Teuchos::ScalarTraits<ScalarT>::rmin();
184 for (
typename std::vector<ScalarT>::iterator
field = mins.begin();
field != mins.end(); ++
field)
185 *
field = Teuchos::ScalarTraits<ScalarT>::rmax();
189template<
typename EvalT,
typename Traits>
195 this->postprocess(*(global_data->os));
199template<
typename EvalT,
typename TRAITS>
203 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"SHOULD NEVER BE CALLED!");
210 Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM,
static_cast<int>(1), &total_volume, &global_total_volume);
211 Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM,
static_cast<int>(averages.size()), &averages[0], &global_averages[0]);
212 Teuchos::reduceAll(*comm, Teuchos::REDUCE_MAX,
static_cast<int>(maxs.size()), &maxs[0], &global_maxs[0]);
213 Teuchos::reduceAll(*comm, Teuchos::REDUCE_MIN,
static_cast<int>(mins.size()), &mins[0], &global_mins[0]);
215 for (std::vector<ScalarT>::size_type i = 0; i < field_values.size(); ++i)
216 global_averages[i] /= global_total_volume;
218 if (comm->getRank() == 0) {
222 std::size_t precision = 8;
223 os << std::scientific << std::showpoint << std::setprecision(precision) << std::left;
225 std::size_t name_width = 0;
226 for (std::vector<ScalarT>::size_type i = 0; i < field_values.size(); ++i)
227 name_width = std::max(name_width,field_values[i].fieldTag().name().size());
229 std::size_t value_width = precision + 7;
231 os << std::setw(name_width) <<
"Field"
232 <<
" " << std::setw(value_width) <<
"Average"
233 <<
" " << std::setw(value_width) <<
"Maximum (@IP)"
234 <<
" " << std::setw(value_width) <<
"Minimum (@IP)"
237 for (std::vector<ScalarT>::size_type i = 0; i < field_values.size(); ++i) {
238 os << std::setw(name_width) << field_values[i].fieldTag().name()
239 <<
" " << std::setw(value_width) << global_averages[i]
240 <<
" " << std::setw(value_width) << global_maxs[i]
241 <<
" " << std::setw(value_width) << global_mins[i] << std::endl;
249template<
typename EvalT,
typename TRAITS>
252 return tmp.fieldTag();
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.
GlobalStatistics(const Teuchos::ParameterList &p)
void postprocess(std::ostream &os)
void postEvaluate(typename Traits::PostEvalData d)
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
const PHX::FieldTag & getRequiredFieldTag()
void preEvaluate(typename Traits::PreEvalData d)
void evaluateFields(typename Traits::EvalData d)
int num_cells
DEPRECATED - use: numCells()
void StringTokenizer(std::vector< std::string > &tokens, const std::string &str, const std::string delimiters, bool trim)
Tokenize a string, put tokens in a vector.
std::vector< int >::size_type getIntegrationRuleIndex(int ir_degree, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_