Skip to content

Commit 51be740

Browse files
Rohit-KakodkarRohit-Kakodkar
authored
Merge pull request #1329 from PrincetonUniversity/issue-1281
Issue 1281 - Compute material properties using MESHFFEM3D mesh
2 parents b39bf84 + 3813efd commit 51be740

File tree

40 files changed

+1472
-111
lines changed

40 files changed

+1472
-111
lines changed

‎CMakeLists.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -362,6 +362,7 @@ add_library(
362362
src/mesh/dim3/generate_database/tags/tags.cpp
363363

364364
src/mesh/dim3/meshfem3d/adjacency_graph/assert_symmetry.cpp
365+
src/mesh/dim3/meshfem3d/tags/tags.cpp
365366
src/mesh/dim3/meshfem3d/mesh.cpp
366367
)
367368

‎core/specfem/assembly/assembly/dim3/assembly.cpp‎

Lines changed: 69 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -70,6 +70,75 @@ specfem::assembly::assembly<specfem::dimension::type::dim3>::assembly(
7070
return;
7171
}
7272

73+
specfem::assembly::assembly<specfem::dimension::type::dim3>::assembly(
74+
const specfem::mesh::meshfem3d::mesh<dimension_tag> &mesh,
75+
const specfem::quadrature::quadratures &quadratures,
76+
std::vector<std::shared_ptr<specfem::sources::source<dimension_tag> > >
77+
&sources,
78+
const std::vector<
79+
std::shared_ptr<specfem::receivers::receiver<dimension_tag> > >
80+
&receivers,
81+
const std::vector<specfem::wavefield::type> &stypes, const type_real t0,
82+
const type_real dt, const int max_timesteps, const int max_sig_step,
83+
const int nsteps_between_samples,
84+
const specfem::simulation::type simulation,
85+
const bool allocate_boundary_values,
86+
const std::shared_ptr<specfem::io::reader> &property_reader) {
87+
88+
const int nspec = mesh.nspec;
89+
const int ngllz = mesh.element_grid.ngllz;
90+
const int nglly = mesh.element_grid.nglly;
91+
const int ngllx = mesh.element_grid.ngllx;
92+
const int ngnod = mesh.control_nodes.ngnod;
93+
94+
this->mesh = { nspec,
95+
ngnod,
96+
ngllz,
97+
nglly,
98+
ngllx,
99+
mesh.adjacency_graph,
100+
mesh.control_nodes,
101+
quadratures };
102+
103+
this->element_types = { nspec, ngllz, nglly, ngllx, this->mesh, mesh.tags };
104+
105+
this->jacobian_matrix = { this->mesh };
106+
107+
this->properties = { nspec, ngllz, nglly,
108+
ngllx, mesh.materials, this->element_types };
109+
110+
this->kernels = { this->mesh.nspec, this->mesh.element_grid.ngllz,
111+
this->mesh.element_grid.nglly,
112+
this->mesh.element_grid.ngllx, this->element_types };
113+
114+
this->sources = {
115+
sources, this->mesh, this->jacobian_matrix, this->element_types,
116+
t0, dt, max_timesteps
117+
};
118+
this->receivers = {
119+
max_sig_step, dt, t0, nsteps_between_samples, receivers,
120+
stypes, this->mesh, mesh.tags, this->element_types
121+
};
122+
// this->boundaries = { this->mesh.nspec,
123+
// this->mesh.element_grid.ngllz,
124+
// this->mesh.element_grid.ngllx,
125+
// mesh,
126+
// this->mesh,
127+
// this->jacobian_matrix };
128+
// this->coupled_interfaces = { mesh, this->mesh, this->jacobian_matrix,
129+
// this->element_types };
130+
this->fields = { this->mesh, this->element_types, simulation };
131+
132+
// if (allocate_boundary_values)
133+
// this->boundary_values = { max_timesteps, this->mesh, this->element_types,
134+
// this->boundaries };
135+
136+
// Currently done in the mesher!
137+
this->check_jacobian_matrix();
138+
139+
return;
140+
}
141+
73142
std::string
74143
specfem::assembly::assembly<specfem::dimension::type::dim3>::print() const {
75144
std::ostringstream message;

‎core/specfem/assembly/assembly/dim3/assembly.hpp‎

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -99,6 +99,21 @@ template <> struct assembly<specfem::dimension::type::dim3> {
9999
const bool allocate_boundary_values,
100100
const std::shared_ptr<specfem::io::reader> &property_reader);
101101

102+
assembly(
103+
const specfem::mesh::meshfem3d::mesh<dimension_tag> &mesh,
104+
const specfem::quadrature::quadratures &quadratures,
105+
std::vector<std::shared_ptr<specfem::sources::source<dimension_tag> > >
106+
&sources,
107+
const std::vector<
108+
std::shared_ptr<specfem::receivers::receiver<dimension_tag> > >
109+
&receivers,
110+
const std::vector<specfem::wavefield::type> &stypes, const type_real t0,
111+
const type_real dt, const int max_timesteps, const int max_sig_step,
112+
const int nsteps_between_samples,
113+
const specfem::simulation::type simulation,
114+
const bool allocate_boundary_values,
115+
const std::shared_ptr<specfem::io::reader> &property_reader);
116+
102117
assembly() = default;
103118

104119
/**

‎core/specfem/assembly/element_types/dim3/element_types.cpp‎

Lines changed: 97 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -16,8 +16,7 @@ specfem::assembly::element_types<specfem::dimension::type::dim3>::element_types(
1616
}
1717

1818
FOR_EACH_IN_PRODUCT(
19-
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC), PROPERTY_TAG(ISOTROPIC),
20-
BOUNDARY_TAG(NONE)),
19+
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC)),
2120
CAPTURE(elements, h_elements) {
2221
int count = 0;
2322
int index = 0;
@@ -38,13 +37,109 @@ specfem::assembly::element_types<specfem::dimension::type::dim3>::element_types(
3837
Kokkos::deep_copy(_elements_, _h_elements_);
3938
})
4039

40+
FOR_EACH_IN_PRODUCT(
41+
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC), PROPERTY_TAG(ISOTROPIC)),
42+
CAPTURE(elements, h_elements) {
43+
int count = 0;
44+
int index = 0;
45+
46+
for (int ispec = 0; ispec < nspec; ispec++) {
47+
if (medium_tags(ispec) == _medium_tag_ &&
48+
property_tags(ispec) == _property_tag_) {
49+
count++;
50+
}
51+
}
52+
53+
_elements_ =
54+
IndexViewType("specfem::assembly::element_types::elements", count);
55+
_h_elements_ = Kokkos::create_mirror_view(_elements_);
56+
57+
for (int ispec = 0; ispec < nspec; ispec++) {
58+
if (medium_tags(ispec) == _medium_tag_ &&
59+
property_tags(ispec) == _property_tag_) {
60+
_h_elements_(index) = ispec;
61+
index++;
62+
}
63+
}
64+
65+
Kokkos::deep_copy(_elements_, _h_elements_);
66+
})
67+
4168
FOR_EACH_IN_PRODUCT(
4269
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC), PROPERTY_TAG(ISOTROPIC),
4370
BOUNDARY_TAG(NONE)),
4471
CAPTURE(elements, h_elements) {
4572
int count = 0;
4673
int index = 0;
4774

75+
for (int ispec = 0; ispec < nspec; ispec++) {
76+
if (medium_tags(ispec) == _medium_tag_ &&
77+
property_tags(ispec) == _property_tag_ &&
78+
boundary_tags(ispec) == _boundary_tag_) {
79+
count++;
80+
}
81+
}
82+
83+
_elements_ =
84+
IndexViewType("specfem::assembly::element_types::elements", count);
85+
_h_elements_ = Kokkos::create_mirror_view(_elements_);
86+
87+
for (int ispec = 0; ispec < nspec; ispec++) {
88+
if (medium_tags(ispec) == _medium_tag_ &&
89+
property_tags(ispec) == _property_tag_ &&
90+
boundary_tags(ispec) == _boundary_tag_) {
91+
_h_elements_(index) = ispec;
92+
index++;
93+
}
94+
}
95+
96+
Kokkos::deep_copy(_elements_, _h_elements_);
97+
})
98+
}
99+
100+
specfem::assembly::element_types<specfem::dimension::type::dim3>::element_types(
101+
const int nspec, const int ngllz, const int nglly, const int ngllx,
102+
const specfem::assembly::mesh<specfem::dimension::type::dim3> &mesh,
103+
const specfem::mesh::meshfem3d::tags<specfem::dimension::type::dim3> &tags)
104+
: nspec(nspec),
105+
medium_tags("specfem::assembly::element_types::medium_tags", nspec),
106+
property_tags("specfem::assembly::element_types::property_tags", nspec),
107+
boundary_tags("specfem::assembly::element_types::boundary_tags", nspec) {
108+
109+
for (int ispec = 0; ispec < nspec; ispec++) {
110+
medium_tags(ispec) = tags.tags_container(ispec).medium_tag;
111+
property_tags(ispec) = tags.tags_container(ispec).property_tag;
112+
boundary_tags(ispec) = tags.tags_container(ispec).boundary_tag;
113+
}
114+
115+
FOR_EACH_IN_PRODUCT(
116+
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC)),
117+
CAPTURE(elements, h_elements) {
118+
int count = 0;
119+
int index = 0;
120+
for (int ispec = 0; ispec < nspec; ispec++) {
121+
if (medium_tags(ispec) == _medium_tag_) {
122+
count++;
123+
}
124+
}
125+
_elements_ =
126+
IndexViewType("specfem::assembly::element_types::elements", count);
127+
_h_elements_ = Kokkos::create_mirror_view(_elements_);
128+
for (int ispec = 0; ispec < nspec; ispec++) {
129+
if (medium_tags(ispec) == _medium_tag_) {
130+
_h_elements_(index) = ispec;
131+
index++;
132+
}
133+
}
134+
Kokkos::deep_copy(_elements_, _h_elements_);
135+
})
136+
137+
FOR_EACH_IN_PRODUCT(
138+
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC), PROPERTY_TAG(ISOTROPIC)),
139+
CAPTURE(elements, h_elements) {
140+
int count = 0;
141+
int index = 0;
142+
48143
for (int ispec = 0; ispec < nspec; ispec++) {
49144
if (medium_tags(ispec) == _medium_tag_ &&
50145
property_tags(ispec) == _property_tag_) {

‎core/specfem/assembly/element_types/dim3/element_types.hpp‎

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -66,6 +66,27 @@ template <> struct element_types<specfem::dimension::type::dim3> {
6666
const specfem::assembly::mesh<dimension_tag> &mesh,
6767
const specfem::mesh::tags<dimension_tag> &tags);
6868

69+
/**
70+
* @brief Constructor for element_types class in 3D assembly
71+
*
72+
* Initializes the element types container with spectral element mesh
73+
* parameters and mesh configuration data for 3D finite element assembly
74+
* operations.
75+
*
76+
* @param nspec Total number of spectral elements in the mesh
77+
* @param ngllz Number of Gauss-Lobatto-Legendre points in the z-direction
78+
* @param nglly Number of Gauss-Lobatto-Legendre points in the y-direction
79+
* @param ngllx Number of Gauss-Lobatto-Legendre points in the x-direction
80+
* @param mesh Reference to the 3D assembly mesh containing geometric and
81+
* topological information
82+
* @param tags Reference to the mesh tags containing element type
83+
* classifications and material properties
84+
*/
85+
element_types(const int nspec, const int ngllz, const int nglly,
86+
const int ngllx,
87+
const specfem::assembly::mesh<dimension_tag> &mesh,
88+
const specfem::mesh::meshfem3d::tags<dimension_tag> &tags);
89+
6990
Kokkos::View<int *, Kokkos::DefaultHostExecutionSpace>
7091
get_elements_on_host(const specfem::element::medium_tag tag) const;
7192

‎core/specfem/assembly/impl/value_containers/dim3/value_containers.hpp‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -64,7 +64,7 @@ struct value_containers<specfem::dimension::type::dim3, containers_type> {
6464
get_container() const {
6565

6666
FOR_EACH_IN_PRODUCT(
67-
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC_PSV), PROPERTY_TAG(ISOTROPIC)),
67+
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC), PROPERTY_TAG(ISOTROPIC)),
6868
CAPTURE(value) {
6969
if constexpr (_medium_tag_ == MediumTag &&
7070
_property_tag_ == PropertyTag) {

‎core/specfem/assembly/properties/dim3/properties.cpp‎

Lines changed: 29 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -39,3 +39,32 @@ specfem::assembly::properties<specfem::dimension::type::dim3>::properties(
3939

4040
return;
4141
}
42+
43+
specfem::assembly::properties<specfem::dimension::type::dim3>::properties(
44+
const int nspec, const int ngllz, const int nglly, const int ngllx,
45+
const specfem::mesh::meshfem3d::Materials<dimension_tag> &materials,
46+
const specfem::assembly::element_types<dimension_tag> &element_types)
47+
: nspec(nspec), ngllz(ngllz), nglly(nglly), ngllx(ngllx) {
48+
49+
this->property_index_mapping =
50+
Kokkos::View<int *, Kokkos::DefaultExecutionSpace>(
51+
"specfem::assembly::properties::property_index_mapping", nspec);
52+
this->h_property_index_mapping =
53+
Kokkos::create_mirror_view(property_index_mapping);
54+
55+
for (int ispec = 0; ispec < nspec; ++ispec) {
56+
h_property_index_mapping(ispec) = -1;
57+
}
58+
59+
FOR_EACH_IN_PRODUCT(
60+
(DIMENSION_TAG(DIM3), MEDIUM_TAG(ELASTIC), PROPERTY_TAG(ISOTROPIC)),
61+
CAPTURE(value) {
62+
_value_ = specfem::medium::properties_container<
63+
_dimension_tag_, _medium_tag_, _property_tag_>(
64+
element_types.get_elements_on_host(_medium_tag_, _property_tag_),
65+
nspec, ngllz, nglly, ngllx, materials, h_property_index_mapping);
66+
})
67+
68+
Kokkos::deep_copy(property_index_mapping, h_property_index_mapping);
69+
return;
70+
}

‎core/specfem/assembly/properties/dim3/properties.hpp‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -48,8 +48,18 @@ struct properties<specfem::dimension::type::dim3>
4848
const specfem::assembly::mesh<dimension_tag> &mesh,
4949
const specfem::mesh::materials<dimension_tag> &materials);
5050

51+
properties(
52+
const int nspec, const int ngllz, const int nglly, const int ngllx,
53+
const specfem::mesh::meshfem3d::Materials<dimension_tag> &materials,
54+
const specfem::assembly::element_types<dimension_tag> &element_types);
55+
5156
///@}
5257

58+
int nspec; ///< total number of spectral elements
59+
int ngllz; ///< number of quadrature points in z dimension
60+
int nglly; ///< number of quadrature points in y dimension
61+
int ngllx; ///< number of quadrature points in x dimension
62+
5363
/**
5464
* @brief Copy misfit kernel data to host
5565
*

0 commit comments

Comments
 (0)