 |
waves
Basic FE playground
|
Go to the documentation of this file. 1 #ifndef KATOPTRON_NEUMANNLIST_H
2 #define KATOPTRON_NEUMANNLIST_H
10 #include <Tpetra_Map.hpp>
11 #include <Tpetra_Vector.hpp>
12 #include <Teuchos_RCP.hpp>
13 #include <Kokkos_ViewFactory.hpp>
23 template <
typename scalar>
42 NeumannList(
Problem &pbl, Teuchos::RCP<Map> map, Teuchos::RCP<ElementsList> elementsList);
125 template <
typename scalar>
132 const int ensemble_size = ET::size;
136 neumann_number = pbl.
nBCs.size();
139 std::vector<std::vector<global_ordinal_type>> nbcElems = {};
141 for (
auto i = 0; i < neumann_number; ++i)
143 std::vector<global_ordinal_type> mynBcElems = {};
144 for (
auto j = 0; j < pbl.
nBCs[i]->tag->elems.size(); ++j)
148 if (elementsList->getElementType(e) ==
static_cast<int>(tbox::ELTYPE::QUAD4))
150 else if (elementsList->getElementType(e) ==
static_cast<int>(tbox::ELTYPE::TRI3))
155 mynBcElems.push_back(e);
157 nbcElems.push_back(mynBcElems);
158 if (maxNumElemsPernBC < mynBcElems.size())
159 maxNumElemsPernBC = mynBcElems.size();
163 neumann_elements = 1;
164 neumann_elementsList = Kokkos::View<int **, Kokkos::LayoutRight>(
"R", neumann_number, neumann_elements + maxNumElemsPernBC);
166 for (
auto i = 0; i < neumann_number; ++i)
168 neumann_elementsList(i, neumann_size) = nbcElems[i].size();
170 for (
auto j = 0; j < nbcElems[i].size(); ++j)
171 neumann_elementsList(i, neumann_elements + j) = nbcElems[i][j];
176 neumann_values = 1 + numDPN;
177 neumann_values_list = Kokkos::View<scalar ***, Kokkos::LayoutRight>(
"R", neumann_number, 1, 1 + 2 * numDPN);
178 for (
auto i = 0; i < neumann_number; ++i)
180 neumann_values_list(i, 0, neumann_time) = 0.;
181 for (
auto j = 0; j < numDPN; ++j)
182 if (pbl.
nBCs[i]->which_dof[j])
183 neumann_values_list(i, 0, neumann_dofs + j) = 1;
185 neumann_values_list(i, 0, neumann_dofs + j) = 0;
187 for (
auto j = 0; j < numDPN; ++j)
188 for (
int s = 0; s < ensemble_size; ++s)
190 ET::coeff(neumann_values_list(i, 0, neumann_values + j), s) = pbl.
nBCs[i]->values[j * ensemble_size + s];
196 #endif //KATOPTRON_NEUMANNLIST_H
int local_ordinal_type
Definition: Map.h:27
Class which is used to store the list of the Neumann boundary conditions (BC).
Definition: NeumannList.h:24
int global_ordinal_type
Definition: Map.h:28
katoptron::Map::local_ordinal_type local_ordinal_type
Definition: ResultsETI.cpp:4
local_ordinal_type getNeumannNumber()
Return the number of Neumann BC.
Definition: NeumannList.h:47
local_ordinal_type neumann_size
Definition: NeumannList.h:34
scalar getNeumannDof(local_ordinal_type i, local_ordinal_type j)
Return whether there is a surface load for the given BC for a given DOF index.
Definition: NeumannList.h:88
Definition: EnsembleTraits.h:8
local_ordinal_type neumann_dofs
Definition: NeumannList.h:38
Kokkos::View< scalar ***, Kokkos::LayoutRight > neumann_values_list
Definition: NeumannList.h:31
katoptron namespace
Definition: Algebraic.h:18
local_ordinal_type neumann_number
Definition: NeumannList.h:33
local_ordinal_type neumann_time
Definition: NeumannList.h:37
Map::global_ordinal_type global_ordinal_type
Definition: NeumannList.h:27
Map::local_ordinal_type local_ordinal_type
Definition: NeumannList.h:26
NeumannList(Problem &pbl, Teuchos::RCP< Map > map, Teuchos::RCP< ElementsList > elementsList)
NeumannList constructor.
Definition: NeumannList.h:126
local_ordinal_type neumann_elements
Definition: NeumannList.h:35
scalar getNeumannValue(local_ordinal_type i, local_ordinal_type j)
Return the value of the surface load for the given BC for a given DOF index.
Definition: NeumannList.h:103
Class which is used to specify in Python the thermomechanical to solve.
Definition: wProblem.h:19
std::vector< Neumann * > nBCs
Definition: wProblem.h:28
local_ordinal_type getNeumannSize(local_ordinal_type i)
Return the number of elements for a given Neumann BC.
Definition: NeumannList.h:58
katoptron::Map::global_ordinal_type global_ordinal_type
Definition: ResultsETI.cpp:5
local_ordinal_type getNeumannElement(local_ordinal_type i, local_ordinal_type e1)
Return the local ID of the element e1 of a given BC.
Definition: NeumannList.h:73
local_ordinal_type neumann_values
Definition: NeumannList.h:39
Kokkos::View< int **, Kokkos::LayoutRight > neumann_elementsList
Definition: NeumannList.h:30