1 #ifndef KATOPTRON_MATERIALSLIST_H
2 #define KATOPTRON_MATERIALSLIST_H
4 #include "Stokhos_Sacado_Kokkos_MP_Vector.hpp"
18 #include <Kokkos_Core.hpp>
21 #include <Tpetra_Map.hpp>
22 #include <Tpetra_Vector.hpp>
23 #include <Teuchos_RCP.hpp>
24 #include <Kokkos_ViewFactory.hpp>
35 template <
typename scalar>
39 Kokkos::View<scalar ***, Kokkos::LayoutRight>
materials;
66 template <
typename scalar>
69 materials_number = pbl.
media.size();
78 const int ensemble_size = ET::size;
80 materials = Kokkos::View<scalar ***, Kokkos::LayoutRight>(
"R", materials_number, 1, 5);
81 for (
auto i = 0; i < materials_number; ++i)
83 materials(i, 0, materials_T) = (scalar)0;
84 for (
int j = 0; j < ensemble_size; ++j)
86 ET::coeff(materials(i, 0, materials_E), j) = (double)pbl.
media[i]->E_vector[j];
87 ET::coeff(materials(i, 0, materials_nu), j) = (
double)pbl.
media[i]->nu_vector[j];
88 ET::coeff(materials(i, 0, materials_k), j) = (double)pbl.
media[i]->k_vector[j];
89 ET::coeff(materials(i, 0, materials_beta), j) = (
double)pbl.
media[i]->beta_vector[j];
91 materials(i, 0, materials_k) = (scalar)pbl.
media[i]->k;
92 materials(i, 0, materials_beta) = (scalar)pbl.
media[i]->beta;
105 template <
typename scalar>
113 size_t materials_T = 0;
116 const int ensemble_size = ET::size;
118 for (
int j = 0; j < ensemble_size; ++j)
120 typename ET::value_type Tj = 0.;
123 size_t index_1 = materials.extent(1) - 1;
125 if (Tj <= ET::coeff(materials(current_material, index_0, materials_T), j))
127 ET::coeff(E, j) = ET::coeff(materials(current_material, index_0, 1), j);
128 ET::coeff(nu, j) = ET::coeff(materials(current_material, index_0, 2), j);
129 ET::coeff(k, j) = ET::coeff(materials(current_material, index_0, 3), j);
130 ET::coeff(beta, j) = ET::coeff(materials(current_material, index_0, 4), j);
132 else if (Tj > ET::coeff(materials(current_material, index_1, materials_T), j))
134 ET::coeff(E, j) = ET::coeff(materials(current_material, index_1, 1), j);
135 ET::coeff(nu, j) = ET::coeff(materials(current_material, index_1, 2), j);
136 ET::coeff(k, j) = ET::coeff(materials(current_material, index_1, 3), j);
137 ET::coeff(beta, j) = ET::coeff(materials(current_material, index_1, 4), j);
141 while (index_1 != index_0 + 1)
143 size_t index_c = index_0 + (index_1 - index_0) / 2;
145 if (Tj == ET::coeff(materials(current_material, index_c, materials_T), j))
148 index_0 = index_1 - 1;
150 else if (Tj < ET::coeff(materials(current_material, index_c, materials_T), j))
156 double coeff = (Tj - ET::coeff(materials(current_material, index_0, materials_T), j)) / (ET::coeff(materials(current_material, index_1, materials_T), j) - ET::coeff(materials(current_material, index_0, materials_T), j));
158 ET::coeff(E, j) = ET::coeff(materials(current_material, index_0, 1), j) + coeff * (ET::coeff(materials(current_material, index_1, 1), j) - ET::coeff(materials(current_material, index_0, 1), j));
159 ET::coeff(nu, j) = ET::coeff(materials(current_material, index_0, 2), j) + coeff * (ET::coeff(materials(current_material, index_1, 2), j) - ET::coeff(materials(current_material, index_0, 2), j));
160 ET::coeff(k, j) = ET::coeff(materials(current_material, index_0, 3), j) + coeff * (ET::coeff(materials(current_material, index_1, 3), j) - ET::coeff(materials(current_material, index_0, 3), j));
161 ET::coeff(beta, j) = ET::coeff(materials(current_material, index_0, 4), j) + coeff * (ET::coeff(materials(current_material, index_1, 4), j) - ET::coeff(materials(current_material, index_0, 4), j));
170 #endif //KATOPTRON_MATERIALSLIST_H