Skip to content

Commit fc0558a

Browse files
authored
Merge pull request #1612 from PrincetonUniversity/1554-medium2
1554 - Move medium to core Part II: refactor medium_container
2 parents 5bcd127 + 2dbbb97 commit fc0558a

File tree

210 files changed

+1262
-1255
lines changed

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

210 files changed

+1262
-1255
lines changed

core/specfem/assembly/assembly/impl/helper.hpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -178,8 +178,8 @@ class helper {
178178
}();
179179

180180
// Call the compute_wavefield function
181-
specfem::medium::compute_wavefield<dimension_tag, MediumTag,
182-
PropertyTag>(
181+
specfem::medium_physics::compute_wavefield<dimension_tag, MediumTag,
182+
PropertyTag>(
183183
chunk_index, assembly, lagrange_derivative, displacement,
184184
velocity, acceleration, wavefield_type, wavefield);
185185
});

include/medium/impl/accessor.hpp renamed to core/specfem/assembly/impl/domain_accessor.hpp

Lines changed: 3 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -3,9 +3,7 @@
33
#include "specfem/point.hpp"
44
#include <Kokkos_Core.hpp>
55

6-
namespace specfem {
7-
namespace medium {
8-
namespace impl {
6+
namespace specfem::assembly::impl {
97

108
/**
119
* @brief CRTP base for accessing data containers in device/host memory.
@@ -23,7 +21,7 @@ namespace impl {
2321
*
2422
*/
2523
template <specfem::dimension::type DimensionTag, typename DataContainer>
26-
class Accessor {
24+
class DomainAccessor {
2725

2826
private:
2927
constexpr static auto dimension = DimensionTag;
@@ -288,6 +286,4 @@ class Accessor {
288286
}
289287
};
290288

291-
} // namespace impl
292-
} // namespace medium
293-
} // namespace specfem
289+
} // namespace specfem::assembly::impl
Lines changed: 164 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,164 @@
1+
#pragma once
2+
3+
#include "domain_accessor.hpp"
4+
#include "specfem/medium_container.hpp"
5+
6+
namespace specfem::assembly::impl {
7+
8+
/**
9+
* @brief Misfit kernel storage container for seismic inversion.
10+
*
11+
* Template container that stores sensitivity kernels (Frechet derivatives)
12+
* which represent the gradient of
13+
* the misfit function with respect to material parameters and are computed
14+
* from the interaction of forward and adjoint wavefields.
15+
*
16+
* Specializes for different dimension/medium/property combinations and provides
17+
* efficient accumulation operations for kernel computation during adjoint
18+
* simulations. Inherits from `kernels::data_container` for storage and
19+
* `specfem::assembly::impl::DomainAccessor` for device/host data access.
20+
*
21+
* @tparam DimensionTag Spatial dimension (dim2/dim3)
22+
* @tparam MediumTag Physical medium type
23+
* @tparam PropertyTag Material property type
24+
*
25+
*/
26+
template <specfem::dimension::type DimensionTag,
27+
specfem::element::medium_tag MediumTag,
28+
specfem::element::property_tag PropertyTag>
29+
struct domain_kernels;
30+
31+
/**
32+
* @brief 2D misfit kernels container specialization.
33+
*
34+
* Stores sensitivity kernels for 2D spectral elements at all quadrature
35+
* points. Data layout: kernels[element][ngllz][ngllx] where ngllz and ngllx
36+
* are quadrature points in vertical and horizontal directions.
37+
*
38+
* Used for accumulating kernel contributions during adjoint simulation:
39+
* - Forward wavefield provides displacement/velocity/acceleration
40+
* - Adjoint wavefield provides sensitivity information
41+
* - Kernels accumulate the correlation between these fields
42+
*
43+
* @tparam MediumTag Physical medium (acoustic, elastic, poroelastic)
44+
* @tparam PropertyTag Material symmetry (isotropic, anisotropic, cosserat)
45+
*/
46+
template <specfem::element::medium_tag MediumTag,
47+
specfem::element::property_tag PropertyTag>
48+
struct domain_kernels<specfem::dimension::type::dim2, MediumTag, PropertyTag>
49+
: public specfem::medium_container::kernels::data_container<
50+
specfem::dimension::type::dim2, MediumTag, PropertyTag>,
51+
public DomainAccessor<specfem::dimension::type::dim2,
52+
domain_kernels<specfem::dimension::type::dim2,
53+
MediumTag, PropertyTag> > {
54+
55+
/// Base kernels data container type
56+
using base_type = specfem::medium_container::kernels::data_container<
57+
specfem::dimension::type::dim2, MediumTag, PropertyTag>;
58+
using base_type::base_type;
59+
60+
constexpr static auto dimension_tag =
61+
base_type::dimension_tag; ///< 2D spatial dimension
62+
constexpr static auto medium_tag =
63+
base_type::medium_tag; ///< Physical medium type
64+
constexpr static auto property_tag =
65+
base_type::property_tag; ///< Material property type
66+
67+
/// Default constructor for empty kernels container
68+
domain_kernels() = default;
69+
70+
/**
71+
* @brief Construct 2D kernels container for specified elements.
72+
*
73+
* Initializes kernels storage for the given spectral elements and sets up
74+
* the mapping from global element indices to local kernel storage indices.
75+
* All kernel values are initialized to zero for accumulation.
76+
*
77+
* @param elements Element indices to initialize kernels for
78+
* @param ngllz Number of vertical quadrature points per element
79+
* @param ngllx Number of horizontal quadrature points per element
80+
* @param property_index_mapping Output mapping from element index to kernel
81+
* storage index
82+
*/
83+
domain_kernels(
84+
const Kokkos::View<int *, Kokkos::DefaultHostExecutionSpace> elements,
85+
const int ngllz, const int ngllx,
86+
const specfem::kokkos::HostView1d<int> property_index_mapping)
87+
: base_type(elements.extent(0), ngllz, ngllx) {
88+
const int nelement = elements.extent(0);
89+
int count = 0;
90+
for (int i = 0; i < nelement; ++i) {
91+
const int ispec = elements(i);
92+
property_index_mapping(ispec) = count;
93+
count++;
94+
}
95+
}
96+
};
97+
98+
/**
99+
* @brief 3D misfit kernels container specialization.
100+
*
101+
* Stores sensitivity kernels for 3D spectral elements at all quadrature
102+
* points. Data layout: kernels[element][ngllz][nglly][ngllx] where ngllz,
103+
* nglly, and ngllx are quadrature points in z, y, and x directions.
104+
*
105+
* Used for accumulating kernel contributions during 3D adjoint simulation.
106+
* Provides efficient storage
107+
* and access patterns for 3D kernel accumulation operations.
108+
*
109+
* @tparam MediumTag Physical medium (acoustic, elastic)
110+
* @tparam PropertyTag Material symmetry (isotropic, anisotropic)
111+
*/
112+
template <specfem::element::medium_tag MediumTag,
113+
specfem::element::property_tag PropertyTag>
114+
struct domain_kernels<specfem::dimension::type::dim3, MediumTag, PropertyTag>
115+
: public specfem::medium_container::kernels::data_container<
116+
specfem::dimension::type::dim3, MediumTag, PropertyTag>,
117+
public DomainAccessor<specfem::dimension::type::dim3,
118+
domain_kernels<specfem::dimension::type::dim3,
119+
MediumTag, PropertyTag> > {
120+
121+
/// Base kernels data container type
122+
using base_type = specfem::medium_container::kernels::data_container<
123+
specfem::dimension::type::dim3, MediumTag, PropertyTag>;
124+
using base_type::base_type;
125+
126+
constexpr static auto dimension_tag =
127+
base_type::dimension_tag; ///< 3D spatial dimension
128+
constexpr static auto medium_tag =
129+
base_type::medium_tag; ///< Physical medium type
130+
constexpr static auto property_tag =
131+
base_type::property_tag; ///< Material property type
132+
133+
/// Default constructor for empty kernels container
134+
domain_kernels() = default;
135+
136+
/**
137+
* @brief Construct 3D kernels container for specified elements.
138+
*
139+
* Initializes kernels storage for the given spectral elements and sets up
140+
* the mapping from global element indices to local kernel storage indices.
141+
* All kernel values are initialized to zero for accumulation.
142+
*
143+
* @param elements Element indices to initialize kernels for
144+
* @param ngllz Number of z-direction quadrature points per element
145+
* @param nglly Number of y-direction quadrature points per element
146+
* @param ngllx Number of x-direction quadrature points per element
147+
* @param property_index_mapping Output mapping from element index to kernel
148+
* storage index
149+
*/
150+
domain_kernels(
151+
const Kokkos::View<int *, Kokkos::DefaultHostExecutionSpace> elements,
152+
const int ngllz, const int nglly, const int ngllx,
153+
const specfem::kokkos::HostView1d<int> property_index_mapping)
154+
: base_type(elements.extent(0), ngllz, nglly, ngllx) {
155+
const int nelement = elements.extent(0);
156+
int count = 0;
157+
for (int i = 0; i < nelement; ++i) {
158+
const int ispec = elements(i);
159+
property_index_mapping(ispec) = count;
160+
count++;
161+
}
162+
}
163+
};
164+
} // namespace specfem::assembly::impl
Lines changed: 178 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,178 @@
1+
#pragma once
2+
3+
#include "domain_accessor.hpp"
4+
#include "enumerations/medium.hpp"
5+
#include "specfem/medium_container.hpp"
6+
#include <Kokkos_Core.hpp>
7+
8+
namespace specfem::assembly {
9+
class mesh_to_compute_mapping;
10+
template <specfem::dimension::type Dimension> struct mesh;
11+
} // namespace specfem::assembly
12+
13+
namespace specfem::mesh {
14+
template <specfem::dimension::type Dimension> struct materials;
15+
}
16+
17+
namespace specfem::assembly::impl {
18+
19+
/**
20+
* @brief Material properties storage for spectral elements at quadrature
21+
* points.
22+
*
23+
* Template container that stores material properties (density, elastic moduli,
24+
* etc.) for all quadrature points within spectral elements. Specializes for
25+
* different dimension/medium/property combinations and provides efficient
26+
* Kokkos-based storage.
27+
*
28+
* Inherits from `properties::data_container` for storage and `impl::Accessor`
29+
* for element-wise access. Material properties are organized as:
30+
* - **2D**: properties[nspec][ngllz][ngllx]
31+
* - **3D**: properties[nspec][ngllz][nglly][ngllx]
32+
*
33+
* Template parameters select appropriate specialization:
34+
* - **DimensionTag**: Spatial dimension (dim2/dim3)
35+
* - **MediumTag**: Physical medium (acoustic, elastic, poroelastic)
36+
* - **PropertyTag**: Material symmetry (isotropic, anisotropic)
37+
*
38+
* @tparam DimensionTag Spatial dimension
39+
* @tparam MediumTag Physical medium type
40+
* @tparam PropertyTag Material property type
41+
*/
42+
template <specfem::dimension::type DimensionTag,
43+
specfem::element::medium_tag MediumTag,
44+
specfem::element::property_tag PropertyTag>
45+
struct domain_properties;
46+
47+
/**
48+
* @brief 2D material properties container specialization.
49+
*
50+
* Stores material properties at quadrature points for 2D spectral elements.
51+
* Data layout: properties[element][ngllz][ngllx] where ngllz and ngllx are
52+
* quadrature points in vertical and horizontal directions.
53+
*
54+
* @tparam MediumTag Physical medium (acoustic, elastic, poroelastic)
55+
* @tparam PropertyTag Material symmetry (isotropic, anisotropic, cosserat)
56+
*/
57+
template <specfem::element::medium_tag MediumTag,
58+
specfem::element::property_tag PropertyTag>
59+
struct domain_properties<specfem::dimension::type::dim2, MediumTag, PropertyTag>
60+
: public specfem::medium_container::properties::data_container<
61+
specfem::dimension::type::dim2, MediumTag, PropertyTag>,
62+
public DomainAccessor<specfem::dimension::type::dim2,
63+
domain_properties<specfem::dimension::type::dim2,
64+
MediumTag, PropertyTag> > {
65+
66+
/// Base data container type for property storage
67+
using base_type = specfem::medium_container::properties::data_container<
68+
specfem::dimension::type::dim2, MediumTag, PropertyTag>;
69+
using base_type::base_type;
70+
71+
constexpr static auto dimension_tag =
72+
base_type::dimension_tag; ///< 2D spatial dimension
73+
constexpr static auto medium_tag =
74+
base_type::medium_tag; ///< Physical medium type
75+
constexpr static auto property_tag =
76+
base_type::property_tag; ///< Material property type
77+
78+
/// Default constructor for empty container
79+
domain_properties() = default;
80+
81+
/**
82+
* @brief Construct 2D properties from mesh and material data.
83+
*
84+
* Extracts material properties from mesh materials and stores them
85+
* at all quadrature points for the specified spectral elements.
86+
*
87+
* @param elements Element indices to initialize
88+
* @param mesh 2D mesh geometry and connectivity
89+
* @param ngllz Number of vertical quadrature points
90+
* @param ngllx Number of horizontal quadrature points
91+
* @param materials Material database indexed by element
92+
* @param has_gll_model Skip material assignment for GLL models
93+
* @param property_index_mapping Element to property mapping
94+
*/
95+
domain_properties(
96+
const Kokkos::View<int *, Kokkos::DefaultHostExecutionSpace> elements,
97+
const specfem::assembly::mesh<dimension_tag> &mesh, const int ngllz,
98+
const int ngllx, const specfem::mesh::materials<dimension_tag> &materials,
99+
const bool has_gll_model,
100+
const specfem::kokkos::HostView1d<int> property_index_mapping);
101+
102+
/// Device value access disabled for this container type
103+
template <typename PointValues, typename IndexType>
104+
KOKKOS_FORCEINLINE_FUNCTION void
105+
add_device_values(const IndexType &index, PointValues &values) const = delete;
106+
107+
/// Host value access disabled for this container type
108+
template <typename PointValues, typename IndexType>
109+
KOKKOS_FORCEINLINE_FUNCTION void
110+
add_host_values(const IndexType &index, PointValues &values) const = delete;
111+
};
112+
113+
/**
114+
* @brief 3D material properties container specialization.
115+
*
116+
* Stores material properties at quadrature points for 3D spectral elements.
117+
* Data layout: properties[element][ngllz][nglly][ngllx] where ngllz, nglly,
118+
* and ngllx are quadrature points in z, y, and x directions respectively.
119+
*
120+
* @tparam MediumTag Physical medium (acoustic, elastic)
121+
* @tparam PropertyTag Material symmetry (isotropic, anisotropic)
122+
*/
123+
template <specfem::element::medium_tag MediumTag,
124+
specfem::element::property_tag PropertyTag>
125+
struct domain_properties<specfem::dimension::type::dim3, MediumTag, PropertyTag>
126+
: public specfem::medium_container::properties::data_container<
127+
specfem::dimension::type::dim3, MediumTag, PropertyTag>,
128+
public DomainAccessor<specfem::dimension::type::dim3,
129+
domain_properties<specfem::dimension::type::dim3,
130+
MediumTag, PropertyTag> > {
131+
132+
/// Base data container type for property storage
133+
using base_type = specfem::medium_container::properties::data_container<
134+
specfem::dimension::type::dim3, MediumTag, PropertyTag>;
135+
using base_type::base_type;
136+
137+
constexpr static auto dimension_tag =
138+
base_type::dimension_tag; ///< 3D spatial dimension
139+
constexpr static auto medium_tag =
140+
base_type::medium_tag; ///< Physical medium type
141+
constexpr static auto property_tag =
142+
base_type::property_tag; ///< Material property type
143+
144+
/// Default constructor for empty container
145+
domain_properties() = default;
146+
147+
/**
148+
* @brief Construct 3D properties from mesh and material data.
149+
*
150+
* Extracts material properties from mesh materials and stores them
151+
* at all quadrature points for the specified spectral elements.
152+
*
153+
* @param elements Element indices to initialize
154+
* @param nspec Total number of spectral elements
155+
* @param ngllz Number of z-direction quadrature points
156+
* @param nglly Number of y-direction quadrature points
157+
* @param ngllx Number of x-direction quadrature points
158+
* @param materials Material database indexed by element
159+
* @param property_index_mapping Element to property mapping
160+
*/
161+
domain_properties(
162+
const Kokkos::View<int *, Kokkos::DefaultHostExecutionSpace> elements,
163+
const int nspec, const int ngllz, const int nglly, const int ngllx,
164+
const specfem::mesh::materials<dimension_tag> &materials,
165+
const specfem::kokkos::HostView1d<int> property_index_mapping);
166+
167+
/// Device value access disabled for this container type
168+
template <typename PointValues, typename IndexType>
169+
KOKKOS_FORCEINLINE_FUNCTION void
170+
add_device_values(const IndexType &index, PointValues &values) const = delete;
171+
172+
/// Host value access disabled for this container type
173+
template <typename PointValues, typename IndexType>
174+
KOKKOS_FORCEINLINE_FUNCTION void
175+
add_host_values(const IndexType &index, PointValues &values) const = delete;
176+
};
177+
178+
} // namespace specfem::assembly::impl

0 commit comments

Comments
 (0)