waves
Basic FE playground
SourcesList.h
Go to the documentation of this file.
1 #ifndef KATOPTRON_SOURCESLIST_H
2 #define KATOPTRON_SOURCESLIST_H
3 
4 #include "katoptron.h"
5 
6 #include "Map.h"
7 #include "wSource.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 namespace katoptron
16 {
17 
21 template <typename scalar>
23 {
26 
27 private:
28  Kokkos::View<int **, Kokkos::LayoutRight> sources_elementsList;
29  Kokkos::View<scalar ***, Kokkos::LayoutRight> sources_values_list;
30 
34 
37 
38 public:
39  SourcesList();
40  SourcesList(Problem &pbl, Teuchos::RCP<Map> map, Teuchos::RCP<ElementsList> elementsList);
41 
46 
57  {
59  }
60 
72  {
73  return sources_elementsList(i, sources_elements + ea);
74  }
75 
85  {
86  return sources_values_list(i, 0, sources_values);
87  }
88 };
89 
106 template <typename scalar>
107 SourcesList<scalar>::SourcesList(Problem &pbl, Teuchos::RCP<Map> map, Teuchos::RCP<ElementsList> elementsList)
108 {
111 
112  sources_number = pbl.Sources.size();
113  local_ordinal_type maxNumElemsPerSource = 0;
114 
115  std::vector<std::vector<global_ordinal_type>> SourceElems = {};
116 
117  for (auto i = 0; i < sources_number; ++i)
118  {
119  std::vector<global_ordinal_type> mySourceElems = {};
120 
121  for (auto j = 0; j < pbl.Sources[i]->tag->elems.size(); ++j)
122  {
123  local_ordinal_type e = map->mapElems->getLocalElement(pbl.Sources[i]->tag->elems[j]->no - 1);
124 
125  if (elementsList->getElementType(e) == static_cast<int>(tbox::ELTYPE::TETRA4))
126  ;
127  else if (elementsList->getElementType(e) == static_cast<int>(tbox::ELTYPE::HEX8))
128  ;
129  else
130  continue;
131 
132  mySourceElems.push_back(e);
133  }
134 
135  SourceElems.push_back(mySourceElems);
136  if (maxNumElemsPerSource < mySourceElems.size())
137  maxNumElemsPerSource = mySourceElems.size();
138  }
139 
140  sources_size = 0;
141  sources_elements = 1;
142  sources_elementsList = Kokkos::View<int **, Kokkos::LayoutRight>("R", sources_number, sources_elements + maxNumElemsPerSource);
143 
144  for (auto i = 0; i < sources_number; ++i)
145  {
146  sources_elementsList(i, sources_size) = SourceElems[i].size();
147 
148  for (auto j = 0; j < SourceElems[i].size(); ++j)
149  sources_elementsList(i, sources_elements + j) = SourceElems[i][j];
150  }
151 
152  sources_time = 0;
153  sources_values = 1;
154  sources_values_list = Kokkos::View<scalar ***, Kokkos::LayoutRight>("R", sources_number, 1, 2);
155  for (auto i = 0; i < sources_number; ++i)
156  {
157  sources_values_list(i, 0, sources_time) = 0.; //to change latter
158 
159  sources_values_list(i, 0, sources_values) = pbl.Sources[i]->value;
160  }
161 }
162 }; // namespace katoptron
163 
164 #endif //KATOPTRON_SOURCESLIST_H
katoptron::Map::local_ordinal_type
int local_ordinal_type
Definition: Map.h:27
katoptron::SourcesList::sources_values
local_ordinal_type sources_values
Definition: SourcesList.h:36
katoptron::SourcesList::local_ordinal_type
Map::local_ordinal_type local_ordinal_type
Definition: SourcesList.h:24
katoptron::SourcesList::sources_elements
local_ordinal_type sources_elements
Definition: SourcesList.h:33
katoptron::SourcesList::getSourceElement
local_ordinal_type getSourceElement(local_ordinal_type i, local_ordinal_type ea)
Return the local ID of the element e1 of a source.
Definition: SourcesList.h:71
katoptron::SourcesList::getSourceValue
scalar getSourceValue(local_ordinal_type i)
Return the value of the source.
Definition: SourcesList.h:84
katoptron::Map::global_ordinal_type
int global_ordinal_type
Definition: Map.h:28
katoptron::SourcesList
Class which is used to store the list of the volumetric heat sources.
Definition: SourcesList.h:22
local_ordinal_type
katoptron::Map::local_ordinal_type local_ordinal_type
Definition: ResultsETI.cpp:4
katoptron::SourcesList::getSourceNumber
local_ordinal_type getSourceNumber()
Return the number of sources.
Definition: SourcesList.h:45
katoptron::SourcesList::sources_values_list
Kokkos::View< scalar ***, Kokkos::LayoutRight > sources_values_list
Definition: SourcesList.h:29
katoptron::SourcesList::sources_time
local_ordinal_type sources_time
Definition: SourcesList.h:35
wSource.h
katoptron::Problem::Sources
std::vector< Source * > Sources
Definition: wProblem.h:30
katoptron::SourcesList::global_ordinal_type
Map::global_ordinal_type global_ordinal_type
Definition: SourcesList.h:25
katoptron
katoptron namespace
Definition: Algebraic.h:18
katoptron::SourcesList::sources_size
local_ordinal_type sources_size
Definition: SourcesList.h:32
katoptron::Problem
Class which is used to specify in Python the thermomechanical to solve.
Definition: wProblem.h:19
katoptron::SourcesList::SourcesList
SourcesList()
katoptron::SourcesList::sources_number
local_ordinal_type sources_number
Definition: SourcesList.h:31
katoptron.h
katoptron::SourcesList::sources_elementsList
Kokkos::View< int **, Kokkos::LayoutRight > sources_elementsList
Definition: SourcesList.h:28
Map.h
katoptron::SourcesList::getSourceSize
local_ordinal_type getSourceSize(local_ordinal_type i)
Return the number of elements for a given source.
Definition: SourcesList.h:56
global_ordinal_type
katoptron::Map::global_ordinal_type global_ordinal_type
Definition: ResultsETI.cpp:5