43#ifndef PANZER_RESPONSE_SCATTER_EVALUATOR_FUNCTIONAL_HPP
44#define PANZER_RESPONSE_SCATTER_EVALUATOR_FUNCTIONAL_HPP
49#include "PanzerDiscFE_config.hpp"
56#include "Phalanx_Evaluator_Macros.hpp"
57#include "Phalanx_MDField.hpp"
67 virtual void scatterDerivative(
const PHX::MDField<const panzer::Traits::Jacobian::ScalarT,panzer::Cell> & cellIntegral,
70 const std::vector<Teuchos::ArrayRCP<double> > & dgdx)
const = 0;
72#ifdef Panzer_BUILD_HESSIAN_SUPPORT
73 virtual void scatterHessian(
const PHX::MDField<const panzer::Traits::Hessian::ScalarT,panzer::Cell> & cellIntegral,
76 const std::vector<Teuchos::ArrayRCP<double> > & d2gdx2)
const = 0;
80template <
typename LO,
typename GO>
85 if(globalIndexer!=Teuchos::null)
86 ugis_.push_back(globalIndexer);
92 void scatterDerivative(
const PHX::MDField<const panzer::Traits::Jacobian::ScalarT,panzer::Cell> & cellIntegral,
95 const std::vector<Teuchos::ArrayRCP<double> > & dgdx)
const;
97#ifdef Panzer_BUILD_HESSIAN_SUPPORT
98 void scatterHessian(
const PHX::MDField<const panzer::Traits::Hessian::ScalarT,panzer::Cell> & cellIntegral,
101 const std::vector<Teuchos::ArrayRCP<double> > & d2gdx2)
const;
106 std::vector<Teuchos::RCP<const panzer::GlobalIndexer> >
ugis_;
112template<
typename EvalT,
typename Traits>
114 public PHX::EvaluatorDerived<EvalT, Traits> {
119 const Teuchos::RCP<FunctionalScatterBase> & functionalScatter);
121 const Teuchos::RCP<FunctionalScatterBase> & functionalScatter);
138template <
typename LO,
typename GO>
142 const std::vector<Teuchos::ArrayRCP<double> > & dgdx)
const
145 std::string blockId = wda(workset).block_id;
147 std::vector<int> blockOffsets;
150 TEUCHOS_ASSERT(dgdx.size()==ugis_.size());
152 auto cellIntegral_h = Kokkos::create_mirror_view(cellIntegral.get_view());
153 Kokkos::deep_copy(cellIntegral_h, cellIntegral.get_view());
155 const std::vector<std::size_t> & localCellIds = wda(workset).cell_local_ids;
157 for(std::size_t b=0;b<ugis_.size();b++) {
158 int start = blockOffsets[b];
160 auto LIDs = ugis_[b]->getLIDs();
161 auto LIDs_h = Kokkos::create_mirror_view(LIDs);
162 Kokkos::deep_copy(LIDs_h, LIDs);
164 Teuchos::ArrayRCP<double> dgdx_b = dgdx[b];
167 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
168 std::size_t cellLocalId = localCellIds[worksetCellIndex];
171 for(std::size_t i=0;i<LIDs_h.extent(1);i++) {
172 dgdx_b[LIDs_h(cellLocalId, i)] += cellIntegral_h(worksetCellIndex).dx(start+i);
178#ifdef Panzer_BUILD_HESSIAN_SUPPORT
179template <
typename LO,
typename GO>
183 const std::vector<Teuchos::ArrayRCP<double> > & d2gdx2)
const
185 PHX::View<const LO*> LIDs;
188 std::string blockId = wda(workset).block_id;
190 std::vector<int> blockOffsets;
193 TEUCHOS_ASSERT(d2gdx2.size()==ugis_.size());
196 const std::vector<std::size_t> & localCellIds = wda(workset).cell_local_ids;
197 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
198 std::size_t cellLocalId = localCellIds[worksetCellIndex];
200 for(std::size_t b=0;b<ugis_.size();b++) {
201 int start = blockOffsets[b];
203 LIDs = ugis_[b]->getElementLIDs(cellLocalId);
205 Teuchos::ArrayRCP<double> d2gdx2_b = d2gdx2[b];
208 for(std::size_t i=0;i<LIDs.size();i++) {
209 d2gdx2_b[LIDs[i]] += cellIntegral(worksetCellIndex).dx(start+i).dx(0);
Data for determining cell topology and dimensionality.
Wrapper to PHX::EvaluatorWithBaseImpl that implements Panzer-specific helpers.
virtual void scatterDerivative(const PHX::MDField< const panzer::Traits::Jacobian::ScalarT, panzer::Cell > &cellIntegral, panzer::Traits::EvalData workset, WorksetDetailsAccessor &wda, const std::vector< Teuchos::ArrayRCP< double > > &dgdx) const =0
virtual ~FunctionalScatterBase()
virtual void scatterHessian(const PHX::MDField< const panzer::Traits::Hessian::ScalarT, panzer::Cell > &cellIntegral, panzer::Traits::EvalData workset, WorksetDetailsAccessor &wda, const std::vector< Teuchos::ArrayRCP< double > > &d2gdx2) const =0
std::vector< Teuchos::RCP< const panzer::GlobalIndexer > > ugis_
FunctionalScatter(const std::vector< Teuchos::RCP< const panzer::GlobalIndexer > > &ugis)
FunctionalScatter(const Teuchos::RCP< const panzer::GlobalIndexer > &globalIndexer)
void scatterDerivative(const PHX::MDField< const panzer::Traits::Jacobian::ScalarT, panzer::Cell > &cellIntegral, panzer::Traits::EvalData workset, WorksetDetailsAccessor &wda, const std::vector< Teuchos::ArrayRCP< double > > &dgdx) const
void scatterHessian(const PHX::MDField< const panzer::Traits::Hessian::ScalarT, panzer::Cell > &cellIntegral, panzer::Traits::EvalData workset, WorksetDetailsAccessor &wda, const std::vector< Teuchos::ArrayRCP< double > > &d2gdx2) const
Teuchos::RCP< FunctionalScatterBase > scatterObj_
Teuchos::RCP< Response_Functional< EvalT > > responseObj_
Teuchos::RCP< PHX::FieldTag > scatterHolder_
void preEvaluate(typename Traits::PreEvalData d)
void evaluateFields(typename Traits::EvalData d)
PHX::MDField< const ScalarT, panzer::Cell > cellIntegral_
std::string responseName_
ResponseScatterEvaluator_Functional(const std::string &name, const CellData &cd, const Teuchos::RCP< FunctionalScatterBase > &functionalScatter)
A constructor with concrete arguments instead of a parameter list.
void computeBlockOffsets(const std::string &blockId, const std::vector< Teuchos::RCP< GlobalIndexer > > &ugis, std::vector< int > &blockOffsets)