Panzer  Version of the Day
Panzer_BCStrategy_Dirichlet_DefaultImpl_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_BCSTRATEGY_DIRICHLET_DEFAULT_IMPL_IMPL_HPP
44 #define PANZER_BCSTRATEGY_DIRICHLET_DEFAULT_IMPL_IMPL_HPP
45 
46 #include "Teuchos_ParameterList.hpp"
47 #include "Teuchos_RCP.hpp"
48 #include "Teuchos_Assert.hpp"
49 
50 #include "Phalanx_DataLayout_MDALayout.hpp"
51 #include "Phalanx_FieldManager.hpp"
52 
53 #include "Phalanx_MDField.hpp"
54 #include "Phalanx_DataLayout.hpp"
55 #include "Phalanx_DataLayout_MDALayout.hpp"
56 
57 // Evaluators
58 #include "Panzer_PhysicsBlock.hpp"
59 #include "Panzer_PureBasis.hpp"
60 #include "Panzer_Dirichlet_Residual.hpp"
63 #include "Panzer_GatherSolution_Epetra.hpp"
64 #include "Panzer_GatherBasisCoordinates.hpp"
65 #include "Panzer_ScatterDirichletResidual_Epetra.hpp"
66 #include "Panzer_BasisValues_Evaluator.hpp"
67 #include "Panzer_PointValues_Evaluator.hpp"
68 #include "Panzer_DOF.hpp"
70 
71 // ***********************************************************************
72 template <typename EvalT>
75  const Teuchos::RCP<panzer::GlobalData>& global_data,
76  const bool in_check_apply_bc) :
77  panzer::BCStrategy<EvalT>(bc),
79  check_apply_bc(in_check_apply_bc),
80  descriptor_map_built(false)
81 {
82 
83 }
84 
85 // ***********************************************************************
86 template <typename EvalT>
89 {
90 
91 }
92 
93 // ***********************************************************************
94 template <typename EvalT>
97  const panzer::PhysicsBlock& pb,
99  const Teuchos::ParameterList& user_data) const
100 {
101  // for deprecated interface support
102  buildDescriptorMapFromVectors();
103 
104  buildAndRegisterGatherAndOrientationEvaluators(fm,pb,lof,user_data);
105  buildAndRegisterScatterEvaluators(fm,pb,lof,user_data);
106 }
107 
108 // ***********************************************************************
109 
110 template <typename EvalT>
113  const panzer::PhysicsBlock& pb,
115  const Teuchos::ParameterList& /* user_data */) const
116 {
117  using Teuchos::ParameterList;
118  using Teuchos::RCP;
119  using Teuchos::rcp;
120  using std::vector;
121  using std::map;
122  using std::string;
123  using std::pair;
124 
125  // for deprecated interface support
126  buildDescriptorMapFromVectors();
127 
128  // Scatter
129  // for (map<string,string>::const_iterator res_to_dof = residual_to_dof_names_map.begin();
130  // res_to_dof != residual_to_dof_names_map.end(); ++res_to_dof) {
131  for(DescriptorIterator itr=m_provided_dofs_desc.begin();
132  itr!=m_provided_dofs_desc.end(); ++itr) {
133 
134  const DOFDescriptor & desc = itr->second;
135 
136  // there is no residual to scatter
137  if(!desc.residualName.first)
138  continue;
139 
140  // std::string residualName = res_to_dof->first;
141  // std::string dofName = res_to_dof->second;
142  std::string dofName = desc.dofName;
143  std::string residualName = desc.residualName.second;
144 
145  ParameterList p("Scatter: "+residualName + " to " + dofName);
146 
147  // Set name
148  string scatter_field_name = "Dummy Scatter: " + this->m_bc.identifier() + residualName;
149  p.set("Scatter Name", scatter_field_name);
150 
151  // Set basis
152  const vector<pair<string,RCP<panzer::PureBasis> > >& dofBasisPair = pb.getProvidedDOFs();
153  RCP<panzer::PureBasis> basis;
154  for (vector<pair<string,RCP<panzer::PureBasis> > >::const_iterator it =
155  dofBasisPair.begin(); it != dofBasisPair.end(); ++it) {
156  if (it->first == dofName)
157  basis = it->second;
158  }
159 
160  TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::is_null(basis), std::runtime_error,
161  "Error the name \"" << dofName
162  << "\" is not a valid DOF for the boundary condition:\n"
163  << this->m_bc << "\n");
164 
165  p.set("Basis", basis);
166 
167  RCP<vector<string> > residual_names = rcp(new vector<string>);
168  residual_names->push_back(residualName);
169  p.set("Dependent Names", residual_names);
170 
171  RCP<map<string,string> > names_map = rcp(new map<string,string>);
172  names_map->insert(std::make_pair(residualName,dofName));
173  p.set("Dependent Map", names_map);
174 
175  TEUCHOS_TEST_FOR_EXCEPTION(!pb.cellData().isSide(), std::logic_error,
176  "Error - physics block is not a side set!");
177 
178  p.set<int>("Side Subcell Dimension",
179  pb.getBaseCellTopology().getDimension() - 1);
180  p.set<int>("Local Side ID", pb.cellData().side());
181 
182  p.set("Check Apply BC",check_apply_bc);
183 
184  RCP< PHX::Evaluator<panzer::Traits> > op = lof.buildScatterDirichlet<EvalT>(p);
185 
186  this->template registerEvaluator<EvalT>(fm, op);
187 
188  // Require variables
189  {
190  using panzer::Dummy;
191  PHX::Tag<typename EvalT::ScalarT> tag(scatter_field_name,
192  rcp(new PHX::MDALayout<Dummy>(0)));
193  fm.template requireField<EvalT>(tag);
194  }
195 
196  }
197 }
198 
199 // ***********************************************************************
200 
201 template <typename EvalT>
204  const panzer::PhysicsBlock& pb,
206  const Teuchos::ParameterList& /* user_data */) const
207 {
208  using Teuchos::ParameterList;
209  using Teuchos::RCP;
210  using Teuchos::rcp;
211  using std::vector;
212  using std::map;
213  using std::string;
214  using std::pair;
215 
216  // for deprecated interface support
217  buildDescriptorMapFromVectors();
218 
219  // volume cell data object (used for handling vector valued fields)
221 
222  // **************************
223  // Coordinates for basis functions (no integration points needed)
224  // **************************
225  {
226  const std::map<std::string,Teuchos::RCP<panzer::PureBasis> > & bases = pb.getBases();
227  for (std::map<std::string,Teuchos::RCP<panzer::PureBasis> >::const_iterator it=bases.begin();
228  it!=bases.end();it++) {
229 
230  Teuchos::RCP<panzer::PureBasis> basis = it->second;
231 
232  // add basis coordinates no matter what
233  {
234  RCP< PHX::Evaluator<panzer::Traits> > basis_op
236  this->template registerEvaluator<EvalT>(fm, basis_op);
237  }
238 
239  // add point values and basis values
240  if(basis->isVectorBasis()) {
241  RCP<const panzer::PointRule> pointRule = rcp(new panzer::PointRule(basis->name()+":BasisPoints",basis->cardinality(),cellData));
242 
243  {
244  RCP< PHX::Evaluator<panzer::Traits> > eval
245  = rcp(new panzer::PointValues_Evaluator<EvalT,panzer::Traits>(pointRule,basis));
246 
247  this->template registerEvaluator<EvalT>(fm, eval);
248  }
249  {
250  // note basis values are not constructed!
251  RCP< PHX::Evaluator<panzer::Traits> > eval
252  = rcp(new panzer::BasisValues_Evaluator<EvalT,panzer::Traits>(pointRule,basis,false));
253 
254  this->template registerEvaluator<EvalT>(fm, eval);
255  }
256  }
257  }
258  }
259 
260  // Gather
261  for(DescriptorIterator itr=m_provided_dofs_desc.begin();
262  itr!=m_provided_dofs_desc.end(); ++itr) {
263 
264  // get the dofName from the descriptor
265  std::string dofName = itr->second.dofName;
266  std::string fieldDof = !itr->second.timeDerivative.first
267  ? itr->second.dofName : itr->second.timeDerivative.second;
268 
269  const vector<pair<string,RCP<panzer::PureBasis> > >& dofBasisPair = pb.getProvidedDOFs();
270  RCP<panzer::PureBasis> basis;
271  for (vector<pair<string,RCP<panzer::PureBasis> > >::const_iterator it =
272  dofBasisPair.begin(); it != dofBasisPair.end(); ++it) {
273  if (it->first == dofName)
274  basis = it->second;
275  }
276 
277  TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::is_null(basis), std::runtime_error,
278  "Error the name \"" << dofName
279  << "\" is not a valid DOF for the boundary condition:\n"
280  << this->m_bc << "\n");
281 
282  {
283  ParameterList p("BC Gather");
284  RCP<vector<string> > gather_field_names_vec = rcp(new vector<string>);
285  RCP<vector<string> > gather_names_vec = rcp(new vector<string>);
286  gather_field_names_vec->push_back(fieldDof);
287  gather_names_vec->push_back(dofName);
288 
289  p.set("DOF Names", gather_field_names_vec);
290  p.set("Indexer Names", gather_names_vec);
291  p.set("Basis", basis);
292  p.set("Use Time Derivative Solution Vector",itr->second.timeDerivative.first);
293 
294  RCP< PHX::Evaluator<panzer::Traits> > op = lof.buildGather<EvalT>(p);
295  this->template registerEvaluator<EvalT>(fm, op);
296  }
297 
298  if(basis->requiresOrientations()) {
299  ParameterList p("Gather Orientation");
300  RCP<vector<string> > gather_field_names_vec = rcp(new vector<string>);
301  RCP<vector<string> > gather_names_vec = rcp(new vector<string>);
302  gather_field_names_vec->push_back(fieldDof);
303  gather_names_vec->push_back(dofName);
304 
305  p.set("DOF Names", gather_field_names_vec);
306  p.set("Indexer Names", gather_names_vec);
307  p.set("Basis", basis);
308 
309  RCP< PHX::Evaluator<panzer::Traits> > op = lof.buildGatherOrientation<EvalT>(p);
310 
311  this->template registerEvaluator<EvalT>(fm, op);
312  }
313 
314  // evaluator a vector basis at the basis points
315  if(basis->isVectorBasis()) {
316  RCP<const panzer::PointRule> pointRule = rcp(new panzer::PointRule(basis->name()+":BasisPoints",basis->cardinality(),cellData));
317 
318  ParameterList p;
319  p.set("Name",fieldDof);
320  p.set("Basis",basis.getConst());
321  p.set("Point Rule",pointRule);
322 
323  RCP< PHX::Evaluator<panzer::Traits> > eval
325  this->template registerEvaluator<EvalT>(fm, eval);
326  }
327 
328  }
329 
330  // Dirichlet Residual: residual = dof_value - target_value
331  for(DescriptorIterator itr=m_provided_dofs_desc.begin();
332  itr!=m_provided_dofs_desc.end(); ++itr) {
333 
334  const DOFDescriptor & desc = itr->second;
335 
336  // there is no residual to scatter
337  if(!desc.residualName.first)
338  continue;
339 
340  // std::string dofName = desc.dofName;
341  std::string dofName = !itr->second.timeDerivative.first
342  ? itr->second.dofName : itr->second.timeDerivative.second;
343  std::string residualName = desc.residualName.second;
344  std::string targetName = desc.targetName.second;
345 
346  const vector<pair<string,RCP<panzer::PureBasis> > >& dofBasisPair = pb.getProvidedDOFs();
347  RCP<panzer::PureBasis> basis;
348  for (vector<pair<string,RCP<panzer::PureBasis> > >::const_iterator it =
349  dofBasisPair.begin(); it != dofBasisPair.end(); ++it) {
350  if (it->first == itr->second.dofName)
351  basis = it->second;
352  }
353 
354  TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::is_null(basis), std::runtime_error,
355  "Error the name \"" << itr->second.dofName
356  << "\" is not a valid DOF for the boundary condition:\n"
357  << this->m_bc << "\n");
358 
359  if(basis->isScalarBasis() || desc.coefficientResidual) {
360  ParameterList p("Dirichlet Residual: "+residualName + " to " + dofName);
361  p.set("Residual Name", residualName);
362  p.set("DOF Name", dofName);
363  p.set("Value Name", targetName);
364  p.set("Data Layout", basis->functional);
365 
366  RCP< PHX::Evaluator<panzer::Traits> > op =
368 
369  this->template registerEvaluator<EvalT>(fm, op);
370  }
371  // This assumes that dofs on faces are all named "<dof>_face"
372  else if(basis->isVectorBasis()&&basis->supportsDiv()) {
373  RCP<const panzer::PointRule> pointRule = rcp(new panzer::PointRule(basis->name()+":BasisPoints",basis->cardinality(),cellData));
374 
375  ParameterList p;
376  p.set("Residual Name", residualName);
377  p.set("DOF Name", dofName);
378  p.set("Value Name", targetName);
379  p.set("Basis", basis.getConst());
380  p.set("Point Rule", pointRule);
381 
382  RCP< PHX::Evaluator<panzer::Traits> > op =
384 
385  this->template registerEvaluator<EvalT>(fm, op);
386  }
387  else if(basis->isVectorBasis()) {
388  RCP<const panzer::PointRule> pointRule = rcp(new panzer::PointRule(basis->name()+":BasisPoints",basis->cardinality(),cellData));
389 
390  ParameterList p;
391  p.set("Residual Name", residualName);
392  p.set("DOF Name", dofName);
393  p.set("Value Name", targetName);
394  p.set("Basis", basis.getConst());
395  p.set("Point Rule", pointRule);
396 
397  RCP< PHX::Evaluator<panzer::Traits> > op =
399 
400  this->template registerEvaluator<EvalT>(fm, op);
401  }
402 
403  }
404 }
405 
406 // ***********************************************************************
407 
408 template <typename EvalT>
411 {
412  if(descriptor_map_built)
413  return;
414 
415  // build the reverse lookup (this assumes that the map is invertible, though it should be!)
416  std::map<std::string,std::string> dof_names_to_residual_map;
417  for(std::map<std::string,std::string>::const_iterator itr=residual_to_dof_names_map.begin();
418  itr!=residual_to_dof_names_map.end();++itr) {
419  dof_names_to_residual_map[itr->second] = itr->first;
420  }
421 
422  for(std::size_t i=0;i<required_dof_names.size();i++) {
423  std::string dof_name = required_dof_names[i];
424 
425  // add the DOF right away
426  const_cast<BCStrategy_Dirichlet_DefaultImpl<EvalT> *>(this)->addDOF(dof_name);
427 
428  // add target information if its required
429  if(dof_names_to_residual_map.find(dof_name)!=dof_names_to_residual_map.end()) {
430  std::string residual_name = dof_names_to_residual_map[dof_name];
431  std::string target_name = residual_to_target_field_map.find(residual_name)->second;
432 
433  // add the descriptor from the data structures: Note the const_cast is neccessary
434  // only because this is protecting deprecated code
435  const_cast<BCStrategy_Dirichlet_DefaultImpl<EvalT> *>(this)->
436  addTarget(target_name,
437  dof_name,
438  residual_name);
439  }
440  }
441 
442  descriptor_map_built = true;
443 }
444 
445 // ***********************************************************************
446 
447 template <typename EvalT>
449 addDOF(const std::string & dofName)
450 {
451  DescriptorIterator itr = m_provided_dofs_desc.find(dofName);
452  TEUCHOS_ASSERT(itr==m_provided_dofs_desc.end());
453 
454  DOFDescriptor & desc = m_provided_dofs_desc[dofName];
455  desc.dofName = dofName;
456 }
457 
458 // ***********************************************************************
459 
460 template <typename EvalT>
462 addTarget(const std::string & targetName,
463  const std::string & dofName,
464  const std::string & residualName)
465 {
466  typename std::map<std::string,DOFDescriptor>::iterator itr = m_provided_dofs_desc.find(dofName);
467  TEUCHOS_ASSERT(itr!=m_provided_dofs_desc.end());
468 
469  DOFDescriptor & desc = itr->second;
470  desc.dofName = dofName;
471  desc.coefficientResidual = false;
472  desc.targetName = std::make_pair(true,targetName);
473  desc.residualName = (residualName=="") ? std::make_pair(true,"RESIDUAL_"+dofName)
474  : std::make_pair(true,residualName);
475 }
476 
477 template <typename EvalT>
479 addCoefficientTarget(const std::string & targetName,
480  const std::string & dofName,
481  const std::string & residualName)
482 {
483  typename std::map<std::string,DOFDescriptor>::iterator itr = m_provided_dofs_desc.find(dofName);
484  TEUCHOS_ASSERT(itr!=m_provided_dofs_desc.end());
485 
486  DOFDescriptor & desc = itr->second;
487  desc.dofName = dofName;
488  desc.coefficientResidual = true;
489  desc.targetName = std::make_pair(true,targetName);
490  desc.residualName = (residualName=="") ? std::make_pair(true,"RESIDUAL_"+dofName)
491  : std::make_pair(true,residualName);
492 }
493 
494 // ***********************************************************************
495 
496 template <typename EvalT>
498 addDotTarget(const std::string & targetName,
499  const std::string & dofName,
500  const std::string & dotName,
501  const std::string & residualName)
502 {
503  typename std::map<std::string,DOFDescriptor>::iterator itr = m_provided_dofs_desc.find(dofName);
504  TEUCHOS_ASSERT(itr!=m_provided_dofs_desc.end());
505 
506  DOFDescriptor & desc = itr->second;
507  desc.dofName = dofName;
508  desc.targetName = std::make_pair(true,targetName);
509  desc.timeDerivative = (dotName=="") ? std::make_pair(true,"DXDT_"+dofName)
510  : std::make_pair(true,dotName);
511  desc.residualName = (residualName=="") ? std::make_pair(true,"RESIDUAL_"+dofName)
512  : std::make_pair(true,residualName);
513 }
514 
515 // ***********************************************************************
516 
517 #endif
Teuchos::RCP< PHX::Evaluator< Traits > > buildScatterDirichlet(const Teuchos::ParameterList &pl) const
Use preconstructed dirichlet scatter evaluators.
Object that contains information on the physics and discretization of a block of elements with the SA...
Default implementation for accessing the GlobalData object.
Interpolates basis DOF values to IP DOF values.
Teuchos::RCP< PHX::Evaluator< Traits > > buildGather(const Teuchos::ParameterList &pl) const
Use preconstructed gather evaluators.
BCStrategy_Dirichlet_DefaultImpl(const panzer::BC &bc, const Teuchos::RCP< panzer::GlobalData > &global_data, const bool check_apply_bc=false)
Interpolates basis DOF values to IP DOF Curl values.
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< const shards::CellTopology > getCellTopology() const
Get CellTopology for the base cell.
std::map< std::string, DOFDescriptor >::const_iterator DescriptorIterator
For convenience, declare the DOFDescriptor iterator.
Data for determining cell topology and dimensionality.
void addTarget(const std::string &targetName, const std::string &dofName, const std::string &residualName="")
bool isSide() const
void addDotTarget(const std::string &targetName, const std::string &dofName, const std::string &dotName="", const std::string &residualName="")
const panzer::CellData & cellData() const
void addCoefficientTarget(const std::string &targetName, const std::string &dofName, const std::string &residualName="")
virtual void buildAndRegisterGatherScatterEvaluators(PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &pb, const panzer::LinearObjFactory< panzer::Traits > &lof, const Teuchos::ParameterList &user_data) const
const std::vector< StrPureBasisPair > & getProvidedDOFs() const
const std::map< std::string, Teuchos::RCP< panzer::PureBasis > > & getBases() const
Returns the unique set of bases, key is the unique panzer::PureBasis::name() of the basis...
Stores input information for a boundary condition.
Definition: Panzer_BC.hpp:81
Evaluates a Dirichlet BC residual corresponding to a field value.
Gathers coordinates for the basis function from the workset and stores them in the field manager...
const shards::CellTopology getBaseCellTopology() const
Interpolates basis DOF values to IP DOF values.
std::size_t numCells() 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
Teuchos::RCP< PHX::Evaluator< Traits > > buildGatherOrientation(const Teuchos::ParameterList &pl) const
Use preconstructed gather evaluators.