Panzer  Version of the Day
Panzer_WorksetContainer.cpp
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 
44 
46 
50 #include "Panzer_Dimension.hpp"
51 
52 namespace panzer {
53 
56  : worksetSize_(1)
57 {}
58 WorksetContainer::WorksetContainer(const Teuchos::RCP<const WorksetFactoryBase> & factory,
59  const std::map<std::string,WorksetNeeds> & needs)
60  : wkstFactory_(factory), worksetSize_(0)
61 {
62  // thats all!
63  ebToNeeds_ = needs;
64 }
65 
70  : wkstFactory_(wc.wkstFactory_)
71  , worksetSize_(wc.worksetSize_)
72 {
73 }
74 
79 {
80  this->clearVolumeWorksets();
81  this->clearSideWorksets();
82 }
83 
85 {
86  worksets_.clear();
87 }
88 
90 {
91  sideWorksets_.clear();
92 }
93 
95 setNeeds(const std::string & eBlock,const WorksetNeeds & needs)
96 {
97  clear(); // clear out old worksets
98  ebToNeeds_[eBlock] = needs;
99 }
100 
102 const WorksetNeeds & WorksetContainer::lookupNeeds(const std::string & eBlock) const
103 {
104  std::map<std::string,WorksetNeeds>::const_iterator itr = ebToNeeds_.find(eBlock);
105 
106  TEUCHOS_TEST_FOR_EXCEPTION(itr==ebToNeeds_.end(),std::logic_error,
107  "WorksetContainer::lookupNeeds no WorksetNeeds object is associated "
108  "with the element block \""+eBlock+"\".");
109 
110  return itr->second;
111 }
112 
113 Teuchos::RCP<std::vector<Workset> >
115 {
116  Teuchos::RCP<std::vector<Workset> > worksetVector;
117  WorksetMap::iterator itr = worksets_.find(wd);
118  if(itr==worksets_.end()) {
119  // couldn't find workset, build it!
120  WorksetNeeds needs;
121  if(hasNeeds())
122  needs = lookupNeeds(wd.getElementBlock());
123  worksetVector = wkstFactory_->getWorksets(wd,needs);
124 
125  // apply orientations to the just constructed worksets
126  if(worksetVector!=Teuchos::null && wd.applyOrientations()) {
127  applyOrientations(wd.getElementBlock(),*worksetVector);
128  }
129 
130  if(worksetVector!=Teuchos::null)
131  setIdentifiers(wd,*worksetVector);
132 
133  // store vector for reuse in the future
134  worksets_[wd] = worksetVector;
135  }
136  else
137  worksetVector = itr->second;
138 
139  return worksetVector;
140 }
141 
143 Teuchos::RCP<std::map<unsigned,Workset> >
145 {
146  Teuchos::RCP<std::map<unsigned,Workset> > worksetMap;
147 
148  // this is the key for the workset map
149  SideMap::iterator itr = sideWorksets_.find(desc);
150 
151  if(itr==sideWorksets_.end()) {
152  // couldn't find workset, build it!
153  if (desc.connectsElementBlocks()) {
154  worksetMap = wkstFactory_->getSideWorksets(desc, lookupNeeds(desc.getElementBlock(0)),
155  lookupNeeds(desc.getElementBlock(1)));
156  }
157  else {
158  worksetMap = wkstFactory_->getSideWorksets(desc,lookupNeeds(desc.getElementBlock(0)));
159  }
160 
161  // apply orientations to the worksets for this side
162  if(worksetMap!=Teuchos::null)
163  applyOrientations(desc,*worksetMap);
164 
165  if(worksetMap!=Teuchos::null)
166  setIdentifiers(desc,*worksetMap);
167 
168  // store map for reuse in the future
169  sideWorksets_[desc] = worksetMap;
170  }
171  else {
172  worksetMap = itr->second;
173  }
174 
175  return worksetMap;
176 }
177 
178 
180 setGlobalIndexer(const Teuchos::RCP<const panzer::GlobalIndexer> & ugi)
181 {
182  // apply the orientations for stored worksets
183  applyOrientations(ugi);
184 }
185 
187 addBasis(const std::string & type,int order,const std::string & rep_field)
188 {
189  using Teuchos::RCP;
190  using Teuchos::rcp;
191 
192  for(auto itr=ebToNeeds_.begin();itr!=ebToNeeds_.end();++itr) {
193  WorksetNeeds & needs = itr->second;
194  RCP<PureBasis> basis = rcp(new PureBasis(type,order,needs.cellData));
195 
196  // add in the new basis
197  needs.bases.push_back(basis);
198  needs.rep_field_name.push_back(rep_field);
199  }
200 
201  // clear all arrays, lazy evaluation means it will be rebuilt
202  clear();
203 }
204 
206 applyOrientations(const Teuchos::RCP<const panzer::GlobalIndexer> & ugi)
207 {
208  // this gurantees orientations won't accidently be applied twice.
209  TEUCHOS_ASSERT(globalIndexer_==Teuchos::null);
210 
211  globalIndexer_ = ugi;
212 
213  // this should be created once and stored in an appropriate place
214  TEUCHOS_TEST_FOR_EXCEPTION(globalIndexer_ == Teuchos::null, std::logic_error,
215  "global indexer is not set yet");
217 
218  // loop over volume worksets, apply orientations to each
219  for(WorksetMap::iterator itr=worksets_.begin();
220  itr!=worksets_.end();++itr) {
221  std::string eBlock = itr->first.getElementBlock();
222 
223  applyOrientations(eBlock,*itr->second);
224  }
225 
226  // loop over side worksets, apply orientations to each
227  for(SideMap::iterator itr=sideWorksets_.begin();
228  itr!=sideWorksets_.end();itr++) {
229 
230  applyOrientations(itr->first,*itr->second);
231  }
232 }
233 
235 setIdentifiers(const WorksetDescriptor & wd,std::vector<Workset> & worksets)
236 {
237  std::size_t hash = std::hash<WorksetDescriptor>()(wd); // this is really ugly, is this really a C++ standard?
238  for(std::size_t i=0;i<worksets.size();i++)
239  worksets[i].setIdentifier(hash+i);
240 }
241 
243 setIdentifiers(const WorksetDescriptor & wd,std::map<unsigned,Workset> & workset_map)
244 {
245  std::size_t hash = std::hash<WorksetDescriptor>()(wd); // this is really ugly, is this really a C++ standard?
246  std::size_t offset = 0;
247  for(auto itr : workset_map) {
248  // itr.second.setIdentifier(hash+offset);
249  workset_map[itr.first].setIdentifier(hash+offset);
250 
251  offset++;
252  }
253 }
254 
256 applyOrientations(const std::string & eBlock, std::vector<Workset> & worksets) const
257 {
258  using Teuchos::RCP;
259 
261  // this is for volume worksets //
263 
264  // short circuit if no global indexer exists
265  if((globalIndexer_==Teuchos::null) and (wkstFactory_->getOrientationsInterface() == Teuchos::null)) {
266  Teuchos::FancyOStream fout(Teuchos::rcpFromRef(std::cout));
267  fout.setOutputToRootOnly(0);
268 
269  fout << "Panzer Warning: No global indexer assigned to a workset container or factory. "
270  << "Orientation of the basis for edge basis functions cannot be applied, "
271  << "if those basis functions are used, there will be problems!" << std::endl;
272  return;
273  }
274 
275  // this should be matched to global indexer size (not sure how to retrive it)
276  if(globalIndexer_!=Teuchos::null){
277  TEUCHOS_TEST_FOR_EXCEPTION(orientations_ == Teuchos::null, std::logic_error,
278  "intrepid2 orientation is not constructed");
279  }
280 
281  // loop over each basis requiring orientations, then apply them
283 
284  // Note: It may be faster to loop over the basis pairs on the inside (not really sure)
285 
286  WorksetNeeds needs;
287  if(hasNeeds())
288  needs = lookupNeeds(eBlock);
289 
290  if(needs.bases.size()>0) {
291  // sanity check that we aren't missing something (the old and new "needs" should not be used together)
292  TEUCHOS_ASSERT(needs.getBases().size()==0);
293 
294  for(std::size_t w=0;w<needs.bases.size();w++) {
295  const PureBasis & basis = *needs.bases[w];
296 
297  // no need for this if orientations are not required!
298  if(!basis.requiresOrientations())
299  continue;
300 
301  // build accessors for orientation fields
302  std::vector<Intrepid2::Orientation> ortsPerBlock;
303 
304  // loop over worksets compute and apply orientations
305  for(std::size_t i=0;i<worksets.size();i++) {
306  // break out of the workset loop
307  if(worksets[i].num_cells<=0) continue;
308 
309  for(std::size_t j=0;j<worksets[i].numDetails();j++) {
310  WorksetDetails & details = worksets[i](j);
311 
312  ortsPerBlock.clear();
313  for (int k=0;k<worksets[i].num_cells;++k) {
314  ortsPerBlock.push_back((*orientations_)[details.cell_local_ids[k]]);
315  }
316 
317  for(std::size_t basis_index=0;basis_index<details.bases.size();basis_index++) {
318  Teuchos::RCP<const BasisIRLayout> layout = details.bases[basis_index]->basis_layout;
319 
320  // only apply orientations if its relevant to the current needs
321  if(layout->getBasis()->name()!=basis.name())
322  continue;
323 
324  TEUCHOS_ASSERT(layout!=Teuchos::null);
325  TEUCHOS_ASSERT(layout->getBasis()!=Teuchos::null);
326  if(layout->getBasis()->requiresOrientations()) {
327  // apply orientations for this basis
328  auto & bv = *details.bases[basis_index];
329  if(not bv.orientationsApplied())
330  bv.applyOrientations(ortsPerBlock,(int) worksets[i].num_cells);
331  }
332  }
333  }
334  }
335  } // end for w
336  }
337  else if(needs.getBases().size()>0) {
338  // sanity check that we aren't missing something (the old and new "needs" should not be used together)
339  TEUCHOS_ASSERT(needs.bases.size()==0);
340 
341  // This is for forwards compatibility, the needs now use "getBasis" calls as opposed
342  // to director accessors.
343  for(const auto & bd : needs.getBases()) {
344 
345  // build accessors for orientation fields
346  std::vector<Intrepid2::Orientation> ortsPerBlock;
347 
348  // loop over worksets compute and apply orientations
349  for(std::size_t i=0;i<worksets.size();i++) {
350  // break out of the workset loop
351  if(worksets[i].num_cells<=0) continue;
352 
353  for(std::size_t j=0;j<worksets[i].numDetails();j++) {
354  WorksetDetails & details = worksets[i](j);
355 
356  ortsPerBlock.clear();
357  // for (int k=0;k<worksets[i].num_cells;++k) {
358  for (int k=0;k<details.numOwnedCells();++k) {
359  ortsPerBlock.push_back((*orientations_)[details.cell_local_ids[k]]);
360  }
361 
362  for(const auto & id : needs.getIntegrators()) {
363  // apply orientations for this basis
364  auto & bv = details.getBasisValues(bd,id);
365  if(not bv.orientationsApplied())
366  bv.applyOrientations(ortsPerBlock,(int) worksets[i].num_cells);
367  }
368  }
369  }
370  } // end for w
371  }
372 }
373 
374 
376 applyOrientations(const WorksetDescriptor & desc,std::map<unsigned,Workset> & worksets) const
377 {
378  using Teuchos::RCP;
379 
381  // this is for side worksets //
383 
384  // short circuit if no global indexer exists
385  if((globalIndexer_==Teuchos::null) and (wkstFactory_->getOrientationsInterface() == Teuchos::null)) {
386  Teuchos::FancyOStream fout(Teuchos::rcpFromRef(std::cout));
387  fout.setOutputToRootOnly(0);
388 
389  fout << "Panzer Warning: No global indexer assigned to a workset container or factory. "
390  << "Orientation of the basis for edge basis functions cannot be applied, "
391  << "if those basis functions are used, there will be problems!";
392  return;
393  }
394 
395  // loop over each basis requiring orientations, then apply them
397 
398  // Note: It may be faster to loop over the basis pairs on the inside (not really sure)
399  WorksetNeeds needs;
400  if(hasNeeds())
401  needs = lookupNeeds(desc.getElementBlock());
402 
403  if(needs.bases.size()>0) {
404  // sanity check that we aren't missing something (the old and new "needs" should not be used together)
405  TEUCHOS_ASSERT(needs.getBases().size()==0);
406  for(std::size_t i=0;i<needs.bases.size();i++) {
407  const PureBasis & basis = *needs.bases[i];
408 
409  // no need for this if orientations are not required!
410  if(!basis.requiresOrientations()) continue;
411 
412  // build accessors for orientation fields
413  std::vector<Intrepid2::Orientation> ortsPerBlock;
414 
415  // loop over worksets compute and apply orientations
416  for(std::map<unsigned,Workset>::iterator itr=worksets.begin();
417  itr!=worksets.end();++itr) {
418 
419  // break out of the workset loop
420  if(itr->second.num_cells<=0) continue;
421 
422  for(std::size_t j=0;j<itr->second.numDetails();j++) {
423  WorksetDetails & details = itr->second(j);
424 
425  ortsPerBlock.clear();
426  for (int k=0;k<itr->second.num_cells;++k) {
427  ortsPerBlock.push_back((*orientations_)[details.cell_local_ids[k]]);
428  }
429 
430  for(std::size_t basis_index=0;basis_index<details.bases.size();basis_index++) {
431  Teuchos::RCP<const BasisIRLayout> layout = details.bases[basis_index]->basis_layout;
432 
433  // only apply orientations if its relevant to the current needs
434  if(layout->getBasis()->name()!=basis.name())
435  continue;
436 
437  TEUCHOS_ASSERT(layout!=Teuchos::null);
438  TEUCHOS_ASSERT(layout->getBasis()!=Teuchos::null);
439  if(layout->getBasis()->requiresOrientations()) {
440  // apply orientations for this basis
441  auto & bv = *details.bases[basis_index];
442  if(not bv.orientationsApplied())
443  bv.applyOrientations(ortsPerBlock,(int) itr->second.num_cells);
444  }
445  }
446  }
447  }
448  } // end for i
449  }
450  else if(needs.getBases().size()>0) {
451  // sanity check that we aren't missing something (the old and new "needs" should not be used together)
452  TEUCHOS_ASSERT(needs.bases.size()==0);
453 
454  // This is for forwards compatibility, the needs now use "getBasis" calls as opposed
455  // to director accessors.
456  for(const auto & bd : needs.getBases()) {
457 
458  // build accessors for orientation fields
459  std::vector<Intrepid2::Orientation> ortsPerBlock;
460 
461  // loop over worksets compute and apply orientations
462  for(std::map<unsigned,Workset>::iterator itr=worksets.begin();
463  itr!=worksets.end();++itr) {
464 
465  // break out of the workset loop
466  if(itr->second.num_cells<=0) continue;
467 
468  for(std::size_t j=0;j<itr->second.numDetails();j++) {
469  WorksetDetails & details = itr->second(j);
470 
471  ortsPerBlock.clear();
472  for (int k=0;k<itr->second.num_cells;++k) {
473  ortsPerBlock.push_back((*orientations_)[details.cell_local_ids[k]]);
474  }
475 
476  for(const auto & id : needs.getIntegrators()) {
477  // apply orientations for this basis
478  auto & bv = details.getBasisValues(bd,id);
479  if(not bv.orientationsApplied())
480  bv.applyOrientations(ortsPerBlock,(int) itr->second.num_cells);
481  }
482  }
483  }
484  } // end for w
485  }
486 }
487 
489  const std::vector<std::string> & elementBlockNames,
490  std::map<std::string,Teuchos::RCP<std::vector<Workset> > > & volumeWksts)
491 {
492  for(std::size_t i=0;i<elementBlockNames.size();i++) {
493  WorksetDescriptor wd = blockDescriptor(elementBlockNames[i]);
494  volumeWksts[elementBlockNames[i]] = wc.getWorksets(wd);
495  }
496 }
497 
499  const std::vector<BC> & bcs,
500  std::map<BC,Teuchos::RCP<std::map<unsigned,Workset> >,LessBC> & sideWksts)
501 {
502  for(std::size_t i=0;i<bcs.size();i++) {
503  WorksetDescriptor wd(bcs[i].elementBlockID(),bcs[i].sidesetID());
504  Teuchos::RCP<std::map<unsigned,Workset> > wksts = wc.getSideWorksets(wd);
505  if(wksts!=Teuchos::null)
506  sideWksts[bcs[i]] = wksts;
507  }
508 }
509 
510 }
std::vector< Teuchos::RCP< const PureBasis > > bases
std::string name() const
A unique key that is the combination of the basis type and basis order.
void setGlobalIndexer(const Teuchos::RCP< const panzer::GlobalIndexer > &ugi)
Teuchos::RCP< std::vector< Intrepid2::Orientation > > orientations_
WorksetContainer()
Default contructor, starts with no workset factory objects.
const std::vector< panzer::BasisDescriptor > & getBases() const
Get a list of bases being requested.
Teuchos::RCP< const panzer::GlobalIndexer > globalIndexer_
Teuchos::RCP< const WorksetFactoryBase > wkstFactory_
std::vector< std::string > rep_field_name
void buildIntrepidOrientation(std::vector< Intrepid2::Orientation > &orientation, panzer::ConnManager &connMgr)
Class that provides access to worksets on each element block and side set.
bool requiresOrientations() const
WorksetMap worksets_
Maps element blocks to input physics block objects.
const std::vector< panzer::IntegrationDescriptor > & getIntegrators() const
Get a list of integrators being requested.
bool connectsElementBlocks() const
Identifies this workset as an interface between two element blocks.
Teuchos::RCP< std::map< unsigned, Workset > > getSideWorksets(const WorksetDescriptor &desc)
Access, and construction of side worksets.
Teuchos::RCP< std::vector< Workset > > getWorksets(const WorksetDescriptor &wd)
Access to volume worksets.
void setIdentifiers(const WorksetDescriptor &wd, std::vector< Workset > &worksets)
const WorksetNeeds & lookupNeeds(const std::string &eBlock) const
Look up an input physics block, throws an exception if it can not be found.
void applyOrientations(const Teuchos::RCP< const panzer::GlobalIndexer > &ugi)
void addBasis(const std::string &type, int order, const std::string &rep_field)
WorksetDescriptor blockDescriptor(const std::string &eBlock)
void getSideWorksetsFromContainer(WorksetContainer &wc, const std::vector< BC > &bcs, std::map< BC, Teuchos::RCP< std::map< unsigned, Workset > >, LessBC > &sideWksts)
Stores input information for a boundary condition.
Definition: Panzer_BC.hpp:81
Description and data layouts associated with a particular basis.
void setNeeds(const std::string &eBlock, const WorksetNeeds &needs)
const std::string & getElementBlock(const int block=0) const
Get element block name.
std::map< std::string, WorksetNeeds > ebToNeeds_
How to construct worksets.
void getVolumeWorksetsFromContainer(WorksetContainer &wc, const std::vector< std::string > &elementBlockNames, std::map< std::string, Teuchos::RCP< std::vector< Workset > > > &volumeWksts)