Skip to content

Commit a8a08cd

Browse files
lsawadelsawade
committed
Updated the documentation
1 parent 6f87eae commit a8a08cd

File tree

5 files changed

+207
-86
lines changed

5 files changed

+207
-86
lines changed

‎core/specfem/assembly/info.hpp‎

Lines changed: 48 additions & 43 deletions
Original file line numberDiff line numberDiff line change
@@ -1,65 +1,70 @@
11
#pragma once
22

3-
#include "specfem/element.hpp"
4-
#include "specfem_setup.hpp"
5-
#include "specfem/assembly/info/impl/bounds.hpp"
3+
#include "specfem/assembly/element_types.hpp"
64
#include "specfem/assembly/info/impl/bounding_box.hpp"
5+
#include "specfem/assembly/info/impl/bounds.hpp"
76
#include "specfem/assembly/mesh.hpp"
87
#include "specfem/assembly/properties.hpp"
9-
#include "specfem/assembly/element_types.hpp"
10-
8+
#include "specfem/element.hpp"
9+
#include "specfem_setup.hpp"
1110

1211
namespace specfem::assembly {
1312

14-
/*
15-
! estimation of minimum period resolved
16-
! based on average GLL distance within element and minimum velocity
17-
!
18-
! rule of thumb (Komatitsch et al. 2005):
19-
! "average number of points per minimum wavelength in an element should be around 5."
20-
21-
! average distance between GLL points within this element
22-
avg_distance = elemsize_max / ( NGLLX - 1) ! since NGLLX = NGLLY = NGLLZ
23-
24-
! largest possible minimum period such that number of points per minimum wavelength
25-
! npts = ( min(vpmin,vsmin) * pmax ) / avg_distance is about ~ NPTS_PER_WAVELENGTH
26-
!
27-
! note: obviously, this estimation depends on the choice of points per wavelength
28-
! which is empirical at the moment.
29-
! also, keep in mind that the minimum period is just an estimation and
30-
! there is no such sharp cut-off period for valid synthetics.
31-
! seismograms become just more and more inaccurate for periods shorter than this estimate.
32-
*/
33-
34-
template <specfem::element::dimension_tag DimensionTag>
35-
struct Info {
13+
/**
14+
* @brief Computes and stores mesh statistics and numerical stability
15+
* parameters.
16+
*
17+
* Analyzes the assembled mesh to extract spatial bounds, material property
18+
* ranges, and element geometry statistics. Estimates the minimum resolvable
19+
* period and suggests a time step based on the CFL condition.
20+
*
21+
* The minimum period estimation follows Komatitsch et al. (2005):
22+
* "average number of points per minimum wavelength in an element should be
23+
* around 5."
24+
*
25+
* @tparam DimensionTag Spatial dimension (dim2 or dim3)
26+
*
27+
* @note The minimum period is an empirical estimate, not a sharp cutoff.
28+
* Synthetics become progressively less accurate for shorter periods.
29+
*/
30+
template <specfem::element::dimension_tag DimensionTag> struct Info {
3631
constexpr static auto dimension_tag = DimensionTag; ///< Dimension tag
3732

3833
Info() = default;
3934

35+
/**
36+
* @brief Construct mesh info by analyzing mesh geometry and material
37+
* properties.
38+
*
39+
* @param mesh Assembled mesh containing element geometry and GLL points
40+
* @param properties Material properties at all mesh points
41+
* @param element_types Element classification by medium and property type
42+
*/
4043
Info(const specfem::assembly::mesh<dimension_tag> &mesh,
4144
const specfem::assembly::properties<dimension_tag> &properties,
4245
const specfem::assembly::element_types<dimension_tag> &element_types);
4346

44-
// Computed mesh properties
45-
info::impl::BoundingBox<dimension_tag> domain_bounds;
46-
info::impl::Bounds element_size;
47-
info::impl::Bounds gll_distance;
48-
info::impl::Bounds vp;
49-
info::impl::Bounds vs;
50-
info::impl::Bounds v;
51-
info::impl::Bounds rho;
52-
info::impl::Bounds vp_vs_ratio;
53-
54-
// Suggested time step based on CFL condition
55-
type_real suggested_time_step;
56-
type_real largest_minimum_period;
57-
58-
47+
info::impl::BoundingBox<dimension_tag> domain_bounds; ///< Spatial extent of
48+
///< the mesh domain
49+
info::impl::Bounds element_size; ///< Element corner-to-corner distances
50+
info::impl::Bounds gll_distance; ///< Distances between adjacent GLL points
51+
info::impl::Bounds vp; ///< P-wave velocity range
52+
info::impl::Bounds vs; ///< S-wave velocity range
53+
info::impl::Bounds v; ///< Combined wave velocity range
54+
info::impl::Bounds rho; ///< Density range
55+
info::impl::Bounds vp_vs_ratio; ///< Vp/Vs ratio range
56+
57+
type_real suggested_time_step; ///< Time step satisfying CFL condition
58+
type_real largest_minimum_period; ///< Maximum of minimum resolvable periods
59+
///< across elements
60+
61+
/**
62+
* @brief Generate formatted string representation of mesh statistics.
63+
* @return Multi-line string with labeled mesh properties
64+
*/
5965
std::string string() const;
6066
};
6167

62-
6368
} // namespace specfem::assembly
6469

6570
#include "specfem/assembly/info.tpp"

‎core/specfem/assembly/info/impl/bounding_box.hpp‎

Lines changed: 70 additions & 33 deletions
Original file line numberDiff line numberDiff line change
@@ -8,37 +8,67 @@
88
#include <vector>
99

1010
namespace specfem::assembly::info::impl {
11-
template <specfem::element::dimension_tag DimensionTag>
12-
struct BoundingBox {
11+
12+
/**
13+
* @brief Axis-aligned bounding box for mesh domain.
14+
*
15+
* Stores min/max bounds for each spatial dimension (x, y, z).
16+
* Provides dimension-specific accessors that handle 2D vs 3D cases.
17+
*
18+
* @tparam DimensionTag Spatial dimension (dim2 or dim3)
19+
*/
20+
template <specfem::element::dimension_tag DimensionTag> struct BoundingBox {
1321

1422
constexpr static auto dimension_tag = DimensionTag;
1523
constexpr static int ndim =
1624
specfem::element::dimension<dimension_tag>::dim; ///< Number of dimensions
1725

18-
specfem::datatype::RegisterArray<
19-
Bounds,
20-
Kokkos::extents<size_t, ndim>,
21-
Kokkos::layout_left> bounds_array;
26+
specfem::datatype::RegisterArray<Bounds, Kokkos::extents<size_t, ndim>,
27+
Kokkos::layout_left>
28+
bounds_array; ///< Array of bounds per dimension
2229

2330
BoundingBox() = default;
2431

32+
/**
33+
* @brief Construct 2D bounding box from explicit min/max values.
34+
* @param x_min Minimum x coordinate
35+
* @param x_max Maximum x coordinate
36+
* @param z_min Minimum z coordinate
37+
* @param z_max Maximum z coordinate
38+
*/
2539
template <specfem::element::dimension_tag U = dimension_tag,
26-
typename std::enable_if<U == specfem::element::dimension_tag::dim2>::type
27-
* = nullptr>
28-
BoundingBox(const type_real x_min, const type_real x_max, const type_real z_min,
29-
const type_real z_max)
30-
: bounds_array(Bounds(x_min, x_max), Bounds(z_min, z_max)) {};
31-
40+
typename std::enable_if<
41+
U == specfem::element::dimension_tag::dim2>::type * = nullptr>
42+
BoundingBox(const type_real x_min, const type_real x_max,
43+
const type_real z_min, const type_real z_max)
44+
: bounds_array(Bounds(x_min, x_max), Bounds(z_min, z_max)){};
45+
46+
/**
47+
* @brief Construct 3D bounding box from explicit min/max values.
48+
* @param x_min Minimum x coordinate
49+
* @param x_max Maximum x coordinate
50+
* @param y_min Minimum y coordinate
51+
* @param y_max Maximum y coordinate
52+
* @param z_min Minimum z coordinate
53+
* @param z_max Maximum z coordinate
54+
*/
3255
template <specfem::element::dimension_tag U = dimension_tag,
33-
typename std::enable_if<U == specfem::element::dimension_tag::dim3>::type
34-
* = nullptr>
35-
BoundingBox(const type_real x_min, const type_real x_max, const type_real y_min,
36-
const type_real y_max, const type_real z_min, const type_real z_max)
37-
: bounds_array(Bounds(x_min, x_max), Bounds(y_min, y_max), Bounds(z_min, z_max)) {};
38-
56+
typename std::enable_if<
57+
U == specfem::element::dimension_tag::dim3>::type * = nullptr>
58+
BoundingBox(const type_real x_min, const type_real x_max,
59+
const type_real y_min, const type_real y_max,
60+
const type_real z_min, const type_real z_max)
61+
: bounds_array(Bounds(x_min, x_max), Bounds(y_min, y_max),
62+
Bounds(z_min, z_max)){};
63+
64+
/**
65+
* @brief Construct 2D bounding box from vector of Bounds.
66+
* @param bounds Vector of 2 Bounds objects (x, z)
67+
* @throws std::invalid_argument if bounds.size() != 2
68+
*/
3969
template <specfem::element::dimension_tag U = dimension_tag,
40-
typename std::enable_if<U == specfem::element::dimension_tag::dim2>::type
41-
* = nullptr>
70+
typename std::enable_if<
71+
U == specfem::element::dimension_tag::dim2>::type * = nullptr>
4272
BoundingBox(const std::vector<Bounds> &bounds) {
4373
if (bounds.size() != 2) {
4474
throw std::invalid_argument(
@@ -47,9 +77,14 @@ struct BoundingBox {
4777
bounds_array = decltype(bounds_array)(bounds[0], bounds[1]);
4878
}
4979

80+
/**
81+
* @brief Construct 3D bounding box from vector of Bounds.
82+
* @param bounds Vector of 3 Bounds objects (x, y, z)
83+
* @throws std::invalid_argument if bounds.size() != 3
84+
*/
5085
template <specfem::element::dimension_tag U = dimension_tag,
51-
typename std::enable_if<U == specfem::element::dimension_tag::dim3>::type
52-
* = nullptr>
86+
typename std::enable_if<
87+
U == specfem::element::dimension_tag::dim3>::type * = nullptr>
5388
BoundingBox(const std::vector<Bounds> &bounds) {
5489
if (bounds.size() != 3) {
5590
throw std::invalid_argument(
@@ -58,28 +93,29 @@ struct BoundingBox {
5893
bounds_array = decltype(bounds_array)(bounds[0], bounds[1], bounds[2]);
5994
}
6095

61-
Bounds &x() {
62-
return bounds_array(0);
63-
}
96+
/** @brief Access x-direction bounds. */
97+
Bounds &x() { return bounds_array(0); }
6498

65-
const Bounds &x() const {
66-
return bounds_array(0);
67-
}
99+
/** @brief Access x-direction bounds (const). */
100+
const Bounds &x() const { return bounds_array(0); }
68101

102+
/** @brief Access y-direction bounds (3D only). */
69103
template <specfem::element::dimension_tag U = dimension_tag,
70-
typename std::enable_if<U == specfem::element::dimension_tag::dim3>::type
71-
* = nullptr>
104+
typename std::enable_if<
105+
U == specfem::element::dimension_tag::dim3>::type * = nullptr>
72106
Bounds &y() {
73107
return bounds_array(1);
74108
}
75109

110+
/** @brief Access y-direction bounds (3D only, const). */
76111
template <specfem::element::dimension_tag U = dimension_tag,
77-
typename std::enable_if<U == specfem::element::dimension_tag::dim3>::type
78-
* = nullptr>
112+
typename std::enable_if<
113+
U == specfem::element::dimension_tag::dim3>::type * = nullptr>
79114
const Bounds &y() const {
80115
return bounds_array(1);
81116
}
82117

118+
/** @brief Access z-direction bounds (index 1 for 2D, index 2 for 3D). */
83119
Bounds &z() {
84120
if constexpr (dimension_tag == specfem::element::dimension_tag::dim2) {
85121
return bounds_array(1);
@@ -88,14 +124,15 @@ struct BoundingBox {
88124
}
89125
}
90126

127+
/** @brief Access z-direction bounds (const, index 1 for 2D, index 2 for 3D).
128+
*/
91129
const Bounds &z() const {
92130
if constexpr (dimension_tag == specfem::element::dimension_tag::dim2) {
93131
return bounds_array(1);
94132
} else {
95133
return bounds_array(2);
96134
}
97135
}
98-
99136
};
100137

101138
} // namespace specfem::assembly::info::impl

‎core/specfem/assembly/info/impl/bounds.hpp‎

Lines changed: 35 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -4,30 +4,55 @@
44

55
namespace specfem::assembly::info::impl {
66

7+
/**
8+
* @brief Min/max range for a scalar quantity.
9+
*
10+
* Stores the minimum and maximum values observed for a property
11+
* (e.g., velocity, density, distance). Provides utility methods
12+
* for computing derived quantities.
13+
*/
714
struct Bounds {
815
public:
9-
type_real min;
10-
type_real max;
16+
type_real min; ///< Minimum value
17+
type_real max; ///< Maximum value
1118

19+
/** @brief Default constructor initializes to zero bounds. */
1220
Bounds() : min(0), max(0) {}
1321

14-
Bounds(type_real min_in, type_real max_in)
15-
: min(min_in),
16-
max(max_in) {}
22+
/**
23+
* @brief Construct bounds from explicit min/max values.
24+
* @param min_in Minimum value
25+
* @param max_in Maximum value
26+
*/
27+
Bounds(type_real min_in, type_real max_in) : min(min_in), max(max_in) {}
1728

29+
/** @brief Compute the range (max - min). */
1830
type_real length() const { return this->max - this->min; }
19-
type_real ratio() const {
31+
32+
/**
33+
* @brief Compute the ratio (max / min).
34+
* @throws std::runtime_error if min is zero
35+
*/
36+
type_real ratio() const {
2037
if (this->min == 0) {
21-
throw std::runtime_error("Bounds::ratio(): min is zero, cannot compute ratio.");
38+
throw std::runtime_error(
39+
"Bounds::ratio(): min is zero, cannot compute ratio.");
2240
}
23-
return this->max / this->min; }
41+
return this->max / this->min;
42+
}
43+
44+
/** @brief Compute the midpoint of the range. */
2445
type_real center() const { return 0.5 * (this->max + this->min); }
2546

26-
Bounds&operator=(const type_real value) {
47+
/**
48+
* @brief Set both min and max to the same value.
49+
* @param value Value to assign to both bounds
50+
*/
51+
Bounds &operator=(const type_real value) {
2752
this->min = value;
2853
this->max = value;
2954
return *this;
3055
}
3156
};
3257

33-
} // namespace specfem::assembly::info::impl
58+
} // namespace specfem::assembly::info::impl

‎docs/sections/api/specfem/assembly/index.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,7 @@
1919
edge_types/edge_types
2020
face_types/face_types
2121
fields/fields
22+
info/info
2223
jacobian_matrix/jacobian_matrix
2324
kernels/kernels
2425
mesh/index

0 commit comments

Comments
 (0)