Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_Workset_Builder.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 
43 #include "PanzerDiscFE_config.hpp"
44 
46 
49 
50 template
52 panzer::buildWorksets(const WorksetNeeds & needs,
53  const std::string & elementBlock,
54  const std::vector<std::size_t>& local_cell_ids,
55  const Kokkos::DynRankView<double,PHX::Device>& vertex_coordinates);
56 
57 template
59 panzer::buildBCWorkset(const WorksetNeeds& needs,
60  const std::string& elementBlock,
61  const std::vector<std::size_t>& local_cell_ids,
62  const std::vector<std::size_t>& local_side_ids,
63  const Kokkos::DynRankView<double,PHX::Device>& vertex_coordinates,
64  const bool populate_value_arrays);
65 
66 template
68 panzer::buildBCWorkset(const WorksetNeeds & needs_a,
69  const std::string & blockid_a,
70  const std::vector<std::size_t>& local_cell_ids_a,
71  const std::vector<std::size_t>& local_side_ids_a,
72  const Kokkos::DynRankView<double,PHX::Device> & vertex_coordinates_a,
73  const panzer::WorksetNeeds & needs_b,
74  const std::string & blockid_b,
75  const std::vector<std::size_t>& local_cell_ids_b,
76  const std::vector<std::size_t>& local_side_ids_b,
77  const Kokkos::DynRankView<double,PHX::Device> & vertex_coordinates_b);
78 
79 template
82  const std::string & eblock_a,
83  const std::vector<std::size_t>& local_cell_ids_a,
84  const std::vector<std::size_t>& local_side_ids_a,
85  const Kokkos::DynRankView<double,PHX::Device> & vertex_coordinates_a,
86  const panzer::WorksetNeeds & needs_b,
87  const std::string & eblock_b,
88  const std::vector<std::size_t>& local_cell_ids_b,
89  const std::vector<std::size_t>& local_side_ids_b,
90  const Kokkos::DynRankView<double,PHX::Device> & vertex_coordinates_b);
91 
92 namespace panzer {
93 
94 void populateValueArrays(std::size_t num_cells,bool isSide,const WorksetNeeds & needs,
95  WorksetDetails & details,const Teuchos::RCP<WorksetDetails> other_details)
96 {
97  using Teuchos::RCP;
98  using Teuchos::rcp;
99 
100  // check to see if the number of owned, ghosted and virtual cells are needed
101  // this is required for backwards compatibility
103  if(details.numOwnedCells()==-1)
104  details.setNumberOfCells(num_cells,0,0);
105  if(other_details!=Teuchos::null and other_details->numOwnedCells()==-1)
106  other_details->setNumberOfCells(num_cells,0,0);
107 
108  panzer::MDFieldArrayFactory mdArrayFactory("",true);
109 
110  // setup the integration rules and bases
111 
112  std::vector<RCP<const panzer::PureBasis> > bases;
113  std::vector<RCP<const panzer::IntegrationRule> > int_rules;
114  if(isSide) {
115  const panzer::CellData side_cell_data(num_cells,
116  details.subcell_index,
117  needs.cellData.getCellTopology());
118 
119  for(std::size_t i=0;i<needs.int_rules.size();i++)
120  int_rules.push_back(rcp(new IntegrationRule(needs.int_rules[i]->cubature_degree,side_cell_data)));
121 
122  for(std::size_t i=0;i<needs.bases.size();i++)
123  bases.push_back(rcp(new PureBasis(needs.bases[i]->type(),needs.bases[i]->order(),side_cell_data)));
124  }
125  else {
126  int_rules = needs.int_rules;
127  bases = needs.bases;
128  }
129 
130  details.ir_degrees = rcp(new std::vector<int>(0));
131  details.basis_names = rcp(new std::vector<std::string>(0));
132 
133  for(std::size_t i=0;i<int_rules.size();i++) {
134 
135  details.ir_degrees->push_back(int_rules[i]->cubature_degree);
136 
139  iv2->setupArrays(int_rules[i]);
140  if (Teuchos::nonnull(other_details))
141  iv2->evaluateValues(details.cell_vertex_coordinates, other_details->int_rules[i]->ip_coordinates,num_cells);
142  else
143  iv2->evaluateValues(details.cell_vertex_coordinates,num_cells);
144 
145  details.int_rules.push_back(iv2);
146 
147  // Need to create all combinations of basis/ir pairings
148  for(std::size_t b=0;b<bases.size();b++) {
149  RCP<panzer::BasisIRLayout> b_layout = rcp(new panzer::BasisIRLayout(bases[b],*int_rules[i]));
150  details.basis_names->push_back(b_layout->name());
151 
152  std::size_t int_degree_index = std::distance(details.ir_degrees->begin(),
153  std::find(details.ir_degrees->begin(),
154  details.ir_degrees->end(),
155  int_rules[i]->order()));
157  rcp(new panzer::BasisValues2<double>("",true,true));
158  bv2->setupArrays(b_layout);
159  bv2->evaluateValues(details.int_rules[int_degree_index]->cub_points,
160  details.int_rules[int_degree_index]->jac,
161  details.int_rules[int_degree_index]->jac_det,
162  details.int_rules[int_degree_index]->jac_inv,
163  details.int_rules[int_degree_index]->weighted_measure,
164  details.cell_vertex_coordinates,
165  true,
166  num_cells);
167 
168  details.bases.push_back(bv2);
169  }
170  }
171 
172 }
173 
174 }
std::vector< Teuchos::RCP< const PureBasis > > bases
Teuchos::RCP< std::vector< Workset > > buildEdgeWorksets(const WorksetNeeds &needs_a, const std::string &eblock_a, const std::vector< std::size_t > &local_cell_ids_a, const std::vector< std::size_t > &local_side_ids_a, const ArrayT &vertex_coordinates_a, const WorksetNeeds &needs_b, const std::string &eblock_b, const std::vector< std::size_t > &local_cell_ids_b, const std::vector< std::size_t > &local_side_ids_b, const ArrayT &vertex_coordinates_b)
bool nonnull(const std::shared_ptr< T > &p)
Teuchos::RCP< std::vector< int > > ir_degrees
If workset corresponds to a sub cell, what is the index?
std::vector< Teuchos::RCP< panzer::BasisValues2< double > > > bases
Static basis function data, key is basis name, value is index in the static_bases vector...
Teuchos::RCP< std::vector< std::string > > basis_names
Value corresponds to basis type. Use the offest for indexing.
int numOwnedCells() const
Number of cells owned by this workset.
Teuchos::RCP< std::vector< Workset > > buildWorksets(const WorksetNeeds &needs, const std::string &elementBlock, const std::vector< std::size_t > &local_cell_ids, const ArrayT &vertex_coordinates)
Teuchos::RCP< const shards::CellTopology > getCellTopology() const
Get CellTopology for the base cell.
std::vector< Teuchos::RCP< const IntegrationRule > > int_rules
CellCoordArray cell_vertex_coordinates
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< std::map< unsigned, Workset > > buildBCWorkset(const WorksetNeeds &needs, const std::string &elementBlock, const std::vector< std::size_t > &local_cell_ids, const std::vector< std::size_t > &local_side_ids, const ArrayT &vertex_coordinates, const bool populate_value_arrays=true)
Data for determining cell topology and dimensionality.
void setNumberOfCells(int o_cells, int g_cells, int v_cells)
Provides access to set numbers of cells (required for backwards compatibility)
std::vector< Teuchos::RCP< panzer::IntegrationValues2< double > > > int_rules
std::string name() const
Unique key for workset indexing composed of basis name and point rule name.
Description and data layouts associated with a particular basis.
void populateValueArrays(std::size_t num_cells, bool isSide, const WorksetNeeds &needs, WorksetDetails &details, const Teuchos::RCP< WorksetDetails > other_details)