waves
Basic FE playground
WeightsList.h
Go to the documentation of this file.
1 #ifndef KATOPTRON_WEIGHTSLIST_H
2 #define KATOPTRON_WEIGHTSLIST_H
3 
4 #include "katoptron.h"
5 
6 #include "Map.h"
7 #include "wWeight.h"
8 #include "wTag.h"
9 
10 #include <Tpetra_Map.hpp>
11 #include <Tpetra_Vector.hpp>
12 #include <Teuchos_RCP.hpp>
13 #include <Kokkos_ViewFactory.hpp>
14 
15 #include "EnsembleTraits.h"
16 
17 using Teuchos::RCP;
18 
19 namespace katoptron
20 {
21 
26 template <typename scalar>
28 {
31 
32 private:
33  Kokkos::View<int **, Kokkos::LayoutRight> weights_nodesList;
34  Kokkos::View<scalar **, Kokkos::LayoutRight> weights_values_list;
35 
39 
42 
43 public:
44  WeightsList(Problem &pbl, Teuchos::RCP<Map> map, Teuchos::RCP<ElementsList> elementsList);
45 
50  {
51  return weights_number;
52  }
53 
64  {
66  }
67 
79  {
80  return weights_nodesList(i, weights_nodes + n1);
81  }
82 
94  {
95  return (weights_values_list(i, weights_dofs + j) == scalar(1.));
96  }
97 
109  {
110  return weights_values_list(i, weights_values + j);
111  }
112 };
113 
130 template <typename scalar>
131 WeightsList<scalar>::WeightsList(Problem &pbl, Teuchos::RCP<Map> map, Teuchos::RCP<ElementsList> elementsList)
132 {
133  typedef EnsembleTraits<scalar> ET;
134  const int ensemble_size = ET::size;
135 
136  local_ordinal_type numDPN = map->numPrimalDPN;
137 
138  weights_number = pbl.Weights.size();
139  local_ordinal_type maxNumNodesPerWeight = 0;
140 
141  std::vector<std::vector<global_ordinal_type>> WeightNodes = {};
142 
143  for (auto i = 0; i < weights_number; ++i)
144  {
145  std::vector<global_ordinal_type> tmpMyWeightNodes = {};
146  std::vector<global_ordinal_type> myWeightNodes = {};
147 
148  for (auto j = 0; j < pbl.Weights[i]->nodes.size(); ++j)
149  tmpMyWeightNodes.push_back(pbl.Weights[i]->nodes[j]->row);
150 
151  std::sort(tmpMyWeightNodes.begin(), tmpMyWeightNodes.end());
152 
153  if (tmpMyWeightNodes.size() >= 1)
154  myWeightNodes.push_back(tmpMyWeightNodes[0]);
155 
156  for (auto ii = 1; ii < tmpMyWeightNodes.size(); ++ii)
157  if (tmpMyWeightNodes[ii] != tmpMyWeightNodes[ii - 1])
158  myWeightNodes.push_back(tmpMyWeightNodes[ii]);
159 
160  WeightNodes.push_back(myWeightNodes);
161  if (maxNumNodesPerWeight < myWeightNodes.size())
162  maxNumNodesPerWeight = myWeightNodes.size();
163  }
164 
165  weights_size = 0;
166  weights_nodes = 1;
167  weights_nodesList = Kokkos::View<int **, Kokkos::LayoutRight>("R", weights_number, weights_nodes + maxNumNodesPerWeight);
168 
169  for (auto i = 0; i < weights_number; ++i)
170  {
171  weights_nodesList(i, weights_size) = WeightNodes[i].size();
172 
173  for (auto j = 0; j < WeightNodes[i].size(); ++j)
174  weights_nodesList(i, weights_nodes + j) = WeightNodes[i][j];
175  }
176 
177  weights_dofs = 0;
178  weights_values = weights_dofs + numDPN;
179  weights_values_list = Kokkos::View<scalar **, Kokkos::LayoutRight>("R", weights_number, 2 * numDPN);
180  for (auto i = 0; i < weights_number; ++i)
181  {
182  for (auto j = 0; j < numDPN; ++j)
183  if (pbl.Weights[i]->which_dof[j])
184  weights_values_list(i, weights_dofs + j) = 1;
185  else
186  weights_values_list(i, weights_dofs + j) = 0;
187 
188  for (int s = 0; s < ensemble_size; ++s)
189  {
190  ET::coeff(weights_values_list(i, weights_values + 0), s) = pbl.Weights[i]->x_values[s];
191  ET::coeff(weights_values_list(i, weights_values + 1), s) = pbl.Weights[i]->y_values[s];
192  ET::coeff(weights_values_list(i, weights_values + 2), s) = pbl.Weights[i]->z_values[s];
193  if (numDPN == 4)
194  ET::coeff(weights_values_list(i, weights_values + 3), s) = pbl.Weights[i]->T_values[s];
195  }
196  }
197 }
198 }; // namespace katoptron
199 
200 #endif //KATOPTRON_WEIGHTSLIST_H
katoptron::WeightsList::weights_nodesList
Kokkos::View< int **, Kokkos::LayoutRight > weights_nodesList
Definition: WeightsList.h:33
katoptron::Map::local_ordinal_type
int local_ordinal_type
Definition: Map.h:27
katoptron::WeightsList::getWeightDof
bool getWeightDof(local_ordinal_type i, local_ordinal_type j)
Return whether the DOF j is weighted for a given weighted region.
Definition: WeightsList.h:93
katoptron::WeightsList::weights_nodes
local_ordinal_type weights_nodes
Definition: WeightsList.h:38
katoptron::WeightsList::weights_size
local_ordinal_type weights_size
Definition: WeightsList.h:37
katoptron::WeightsList::weights_dofs
local_ordinal_type weights_dofs
Definition: WeightsList.h:40
katoptron::Map::global_ordinal_type
int global_ordinal_type
Definition: Map.h:28
katoptron::WeightsList
Class which is used to store the list of the weights used for the residual computation.
Definition: WeightsList.h:27
katoptron::Problem::Weights
std::vector< Weight * > Weights
Definition: wProblem.h:32
katoptron::WeightsList::getNode
local_ordinal_type getNode(local_ordinal_type i, local_ordinal_type n1)
Return the global ID of the node n1 of a given weighted region.
Definition: WeightsList.h:78
katoptron::WeightsList::global_ordinal_type
Map::global_ordinal_type global_ordinal_type
Definition: WeightsList.h:30
katoptron::WeightsList::getWeightRegionSize
local_ordinal_type getWeightRegionSize(local_ordinal_type i)
Return the number of nodes for a given weighted region.
Definition: WeightsList.h:63
EnsembleTraits
Definition: EnsembleTraits.h:8
katoptron::WeightsList::weights_values
local_ordinal_type weights_values
Definition: WeightsList.h:41
katoptron
katoptron namespace
Definition: Algebraic.h:18
EnsembleTraits.h
katoptron::WeightsList::weights_number
local_ordinal_type weights_number
Definition: WeightsList.h:36
katoptron::WeightsList::getWeightRegionsNumber
local_ordinal_type getWeightRegionsNumber()
Return the number of weighted regions.
Definition: WeightsList.h:49
katoptron::WeightsList::local_ordinal_type
Map::local_ordinal_type local_ordinal_type
Definition: WeightsList.h:29
katoptron::WeightsList::WeightsList
WeightsList(Problem &pbl, Teuchos::RCP< Map > map, Teuchos::RCP< ElementsList > elementsList)
WeightsList constructor.
Definition: WeightsList.h:131
katoptron::Problem
Class which is used to specify in Python the thermomechanical to solve.
Definition: wProblem.h:19
katoptron.h
katoptron::WeightsList::weights_values_list
Kokkos::View< scalar **, Kokkos::LayoutRight > weights_values_list
Definition: WeightsList.h:34
Map.h
katoptron::WeightsList::getWeightValue
scalar getWeightValue(local_ordinal_type i, local_ordinal_type j)
Return the value of weight for the DOF for a given weighted region.
Definition: WeightsList.h:108
wWeight.h