43#ifndef PANZER_BCSTRATEGY_WEAKDIRICHLET_DEFAULT_IMPL_IMPL_HPP
44#define PANZER_BCSTRATEGY_WEAKDIRICHLET_DEFAULT_IMPL_IMPL_HPP
46#include "Teuchos_ParameterList.hpp"
47#include "Teuchos_RCP.hpp"
48#include "Teuchos_Assert.hpp"
49#include "Phalanx_DataLayout_MDALayout.hpp"
50#include "Phalanx_FieldManager.hpp"
55#include "Phalanx_MDField.hpp"
56#include "Phalanx_DataLayout.hpp"
57#include "Phalanx_DataLayout_MDALayout.hpp"
62#include "Panzer_WeakDirichlet_Residual.hpp"
64#ifdef PANZER_HAVE_EPETRA_STACK
65#include "Panzer_GatherSolution_Epetra.hpp"
66#include "Panzer_ScatterResidual_Epetra.hpp"
69#include "Panzer_Normals.hpp"
72template <
typename EvalT>
75 const Teuchos::RCP<panzer::GlobalData>& global_data) :
83template <
typename EvalT>
91template <
typename EvalT>
96 const Teuchos::ParameterList& user_data)
const
98 buildAndRegisterGatherAndOrientationEvaluators(fm,pb,lof,user_data);
99 buildAndRegisterScatterEvaluators(fm,pb,lof,user_data);
103template <
typename EvalT>
108 const Teuchos::ParameterList& user_data)
const
110 using Teuchos::ParameterList;
122 for (vector<std::tuple<std::string,std::string,std::string,
int,Teuchos::RCP<panzer::PureBasis>,Teuchos::RCP<panzer::IntegrationRule> > >::const_iterator eq =
123 m_residual_contributions.begin(); eq != m_residual_contributions.end(); ++eq) {
125 const string& residual_name = std::get<0>(*eq);
126 const string& dof_name = std::get<1>(*eq);
127 const string& flux_name = std::get<2>(*eq);
129 const RCP<const panzer::PureBasis> basis = std::get<4>(*eq);
130 const RCP<const panzer::IntegrationRule> ir = std::get<5>(*eq);
136 ParameterList p(s.str());
137 p.set<std::string>(
"Name",
"Side Normal");
139 p.set< Teuchos::RCP<panzer::IntegrationRule> >(
"IR", Teuchos::rcp_const_cast<panzer::IntegrationRule>(ir));
140 p.set<
bool>(
"Normalize",
true);
144 this->
template registerEvaluator<EvalT>(fm, op);
149 ParameterList p(
"Neumann Residual: " + residual_name +
" to DOF: " + dof_name);
150 p.set(
"Residual Name", residual_name);
151 p.set(
"DOF Name",dof_name);
152 p.set(
"Flux Name", flux_name);
153 p.set(
"Normal Name",
"Side Normal");
154 p.set(
"Basis", basis);
157 RCP< PHX::Evaluator<panzer::Traits> > op =
160 this->
template registerEvaluator<EvalT>(fm, op);
167template <
typename EvalT>
172 const Teuchos::ParameterList& )
const
174 using Teuchos::ParameterList;
183 for (vector<std::tuple<std::string,std::string,std::string,
int,Teuchos::RCP<panzer::PureBasis>,Teuchos::RCP<panzer::IntegrationRule> > >::const_iterator eq =
184 m_residual_contributions.begin(); eq != m_residual_contributions.end(); ++eq) {
186 const string& residual_name = std::get<0>(*eq);
187 const string& dof_name = std::get<1>(*eq);
188 const RCP<const panzer::PureBasis> basis = std::get<4>(*eq);
189 const RCP<const panzer::IntegrationRule> ir = std::get<5>(*eq);
193 ParameterList p(
"Scatter: "+ residual_name +
" to " + dof_name);
196 string scatter_field_name =
"Dummy Scatter: " + this->m_bc.identifier() + residual_name;
197 p.set(
"Scatter Name", scatter_field_name);
198 p.set(
"Basis", basis);
200 RCP<vector<string> > residual_names = rcp(
new vector<string>);
201 residual_names->push_back(residual_name);
202 p.set(
"Dependent Names", residual_names);
204 RCP<map<string,string> > names_map = rcp(
new map<string,string>);
205 names_map->insert(std::pair<string,string>(residual_name,dof_name));
206 p.set(
"Dependent Map", names_map);
208 RCP< PHX::Evaluator<panzer::Traits> > op = lof.
buildScatter<EvalT>(p);
210 this->
template registerEvaluator<EvalT>(fm, op);
215 PHX::Tag<typename EvalT::ScalarT> tag(scatter_field_name,
216 rcp(
new PHX::MDALayout<Dummy>(0)));
217 fm.template requireField<EvalT>(tag);
226template <
typename EvalT>
230 m_required_dof_names.push_back(required_dof_name);
234template <
typename EvalT>
237 const std::string dof_name,
238 const std::string flux_name,
239 const int integration_order,
242 Teuchos::RCP<panzer::PureBasis> basis = this->getBasis(dof_name,side_pb);
244 Teuchos::RCP<panzer::IntegrationRule> ir = buildIntegrationRule(integration_order,side_pb);
246 m_residual_contributions.push_back(std::make_tuple(residual_name,
255template <
typename EvalT>
256const std::vector<std::tuple<std::string,std::string,std::string,int,Teuchos::RCP<panzer::PureBasis>,Teuchos::RCP<panzer::IntegrationRule> > >
259 return m_residual_contributions;
263template <
typename EvalT>
264Teuchos::RCP<panzer::PureBasis>
268 const std::vector<std::pair<std::string,Teuchos::RCP<panzer::PureBasis> > >& dofBasisPair = side_pb.
getProvidedDOFs();
269 Teuchos::RCP<panzer::PureBasis> basis;
270 for (std::vector<std::pair<std::string,Teuchos::RCP<panzer::PureBasis> > >::const_iterator it =
271 dofBasisPair.begin(); it != dofBasisPair.end(); ++it) {
272 if (it->first == dof_name)
276 TEUCHOS_TEST_FOR_EXCEPTION(is_null(basis), std::runtime_error,
277 "Error the name \"" << dof_name
278 <<
"\" is not a valid DOF for the boundary condition:\n"
279 << this->m_bc <<
"\n");
285template <
typename EvalT>
286Teuchos::RCP<panzer::IntegrationRule>
296template <
typename EvalT>
const panzer::BC bc() const
Returns the boundary condition data for this object.
virtual void buildAndRegisterGatherScatterEvaluators(PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &side_pb, const panzer::LinearObjFactory< panzer::Traits > &lof, const Teuchos::ParameterList &user_data) const
virtual void buildAndRegisterScatterEvaluators(PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &side_pb, const LinearObjFactory< panzer::Traits > &lof, const Teuchos::ParameterList &user_data) const
const std::vector< std::tuple< std::string, std::string, std::string, int, Teuchos::RCP< panzer::PureBasis >, Teuchos::RCP< panzer::IntegrationRule > > > getResidualContributionData() const
Returns information for the residual contribution integrations associated with this Neumann BC.
virtual ~BCStrategy_WeakDirichlet_DefaultImpl()
virtual void buildAndRegisterGatherAndOrientationEvaluators(PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &side_pb, const LinearObjFactory< panzer::Traits > &lof, const Teuchos::ParameterList &user_data) const
Teuchos::RCP< panzer::PureBasis > getBasis(const std::string dof_name, const panzer::PhysicsBlock &side_pb) const
Finds the basis for the corresponding dof_name in the physics block.
Teuchos::RCP< panzer::IntegrationRule > buildIntegrationRule(const int integration_order, const panzer::PhysicsBlock &side_pb) const
Allocates and returns the integration rule associated with an integration order and side physics bloc...
virtual void addResidualContribution(const std::string residual_name, const std::string dof_name, const std::string flux_name, const int integration_order, const panzer::PhysicsBlock &side_pb)
Adds a residual contribution for a neumann condition to a particular equation.
BCStrategy_WeakDirichlet_DefaultImpl(const panzer::BC &bc, const Teuchos::RCP< panzer::GlobalData > &global_data)
virtual void requireDOFGather(const std::string required_dof_name)
Requires that a gather evaluator for the DOF be constructed.
Stores input information for a boundary condition.
Default implementation for accessing the GlobalData object.
Teuchos::RCP< PHX::Evaluator< Traits > > buildScatter(const Teuchos::ParameterList &pl) const
Use preconstructed scatter evaluators.
Object that contains information on the physics and discretization of a block of elements with the SA...
void buildAndRegisterGatherAndOrientationEvaluators(PHX::FieldManager< panzer::Traits > &fm, const panzer::LinearObjFactory< panzer::Traits > &lof, const Teuchos::ParameterList &user_data) const
const panzer::CellData & cellData() const
const std::vector< StrPureBasisPair > & getProvidedDOFs() const
Evaluates a Weak Dirichlet BC residual contribution.