Skip to content

Commit f0c75ab

Browse files
lsawadelsawade
authored
Merge pull request #1630 from PrincetonUniversity/issue-1601-squash-merge
Issue 1601 - Implement check mesh resolution
2 parents 6091976 + a8a08cd commit f0c75ab

File tree

40 files changed

+2224
-17
lines changed

40 files changed

+2224
-17
lines changed

‎core/specfem/assembly/CMakeLists.txt‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
file(GLOB_RECURSE COMPUTE_SOURCE_FILES "./*.cpp" "./**/*.cpp")
1+
file(GLOB_RECURSE COMPUTE_SOURCE_FILES "*.cpp")
22

33
add_library(assembly ${COMPUTE_SOURCE_FILES})
44
add_library(specfem::assembly ALIAS assembly)

‎core/specfem/assembly/assembly/dim2/assembly.cpp‎

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,7 @@
22
#include "enumerations/interface.hpp"
33
#include "specfem/io.hpp"
44
#include "specfem/mesh.hpp"
5+
#include "specfem/assembly/info.hpp"
56

67
specfem::assembly::assembly<specfem::element::dimension_tag::dim2>::assembly(
78
const specfem::mesh::mesh<dimension_tag> &mesh,
@@ -66,6 +67,8 @@ specfem::assembly::assembly<specfem::element::dimension_tag::dim2>::assembly(
6667
this->edge_types, this->mesh };
6768
this->fields = { this->mesh, this->element_types, simulation };
6869

70+
this->info = { this->mesh, this->properties, this->element_types };
71+
6972
if (allocate_boundary_values)
7073
this->boundary_values = { max_timesteps, this->mesh, this->element_types,
7174
this->boundaries };
@@ -102,7 +105,8 @@ specfem::assembly::assembly<specfem::element::dimension_tag::dim2>::print()
102105
<< "Total number of geometric points : "
103106
<< this->mesh.element_grid.ngllz << "\n"
104107
<< "Total number of distinct quadrature points : " << this->mesh.nglob
105-
<< "\n";
108+
<< "\n"
109+
<< this->info.string() << "\n";
106110

107111
int total_elements = 0;
108112

‎core/specfem/assembly/assembly/dim2/assembly.hpp‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -7,6 +7,7 @@
77
#include "specfem/assembly/conforming_interfaces.hpp"
88
#include "specfem/assembly/edge_types.hpp"
99
#include "specfem/assembly/element_types.hpp"
10+
#include "specfem/assembly/info.hpp"
1011
#include "specfem/assembly/fields.hpp"
1112
#include "specfem/assembly/jacobian_matrix.hpp"
1213
#include "specfem/assembly/kernels.hpp"
@@ -155,6 +156,12 @@ template <> struct assembly<specfem::element::dimension_tag::dim2> {
155156
///< the
156157
///< boundaries
157158

159+
/**
160+
* @brief Info
161+
*
162+
*/
163+
specfem::assembly::Info<dimension_tag> info;
164+
158165
///@}
159166

160167
/**
@@ -242,6 +249,7 @@ template <> struct assembly<specfem::element::dimension_tag::dim2> {
242249
*
243250
*/
244251
void check_jacobian_matrix() const;
252+
245253
};
246254

247255
} // namespace specfem::assembly

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

Lines changed: 5 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -70,6 +70,9 @@ specfem::assembly::assembly<specfem::element::dimension_tag::dim3>::assembly(
7070
// Currently done in the mesher!
7171
this->check_jacobian_matrix();
7272

73+
74+
this->info = { this->mesh, this->properties, this->element_types };
75+
7376
return;
7477
}
7578

@@ -82,9 +85,8 @@ specfem::assembly::assembly<specfem::element::dimension_tag::dim3>::print()
8285
<< " Total number of spectral elements : "
8386
<< this->mesh.nspec << "\n"
8487
<< " Total number of quadrature points per element : "
85-
<< this->mesh.element_grid.ngllz << "\n";
86-
// << "Total number of distinct quadrature points : "
87-
// << this->mesh.nglob << "\n";
88+
<< this->mesh.element_grid.ngllz << "\n"
89+
<< this->info.string() << "\n";
8890

8991
int total_elements = 0;
9092

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

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,7 @@
66
#include "specfem/assembly/compute_source_array.hpp"
77
#include "specfem/assembly/conforming_interfaces.hpp"
88
#include "specfem/assembly/fields.hpp"
9+
#include "specfem/assembly/info.hpp"
910
#include "specfem/assembly/jacobian_matrix.hpp"
1011
#include "specfem/assembly/kernels.hpp"
1112
#include "specfem/assembly/mesh.hpp"
@@ -113,6 +114,10 @@ template <> struct assembly<specfem::element::dimension_tag::dim3> {
113114
///< the
114115
///< boundaries
115116

117+
specfem::assembly::Info<dimension_tag>
118+
info; ///< Information about the mesh and
119+
///< simulation
120+
116121
///@}
117122

118123
/**

‎core/specfem/assembly/info.cpp‎

Lines changed: 37 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,37 @@
1+
#include "specfem/assembly/info.hpp"
2+
3+
template <>
4+
std::string specfem::assembly::Info<specfem::element::dimension_tag::dim2>::string() const {
5+
std::ostringstream oss;
6+
oss << "Mesh Information (2D):\n";
7+
oss << " Domain X: ............... [" << domain_bounds.x().min << ", " << domain_bounds.x().max << "]\n";
8+
oss << " Domain Z: ............... [" << domain_bounds.z().min << ", " << domain_bounds.z().max << "]\n";
9+
oss << " VP: ..................... [" << vp.min << ", " << vp.max << "]\n";
10+
oss << " VS: ..................... [" << vs.min << ", " << vs.max << "]\n";
11+
oss << " V: .......................[" << v.min << ", " << v.max << "]\n";
12+
oss << " Rho: .................... [" << rho.min << ", " << rho.max << "]\n";
13+
oss << " Element Size: ........... [" << element_size.min << ", " << element_size.max << "]\n";
14+
oss << " GLL Distance: ........... [" << gll_distance.min << ", " << gll_distance.max << "]\n";
15+
oss << " Largest Minimum Period: . " << largest_minimum_period << "\n";
16+
oss << " Suggested Time Step: .... " << suggested_time_step << "\n";
17+
return oss.str();
18+
}
19+
20+
21+
template <>
22+
std::string specfem::assembly::Info<specfem::element::dimension_tag::dim3>::string() const {
23+
std::ostringstream oss;
24+
oss << "Mesh Information (3D):\n";
25+
oss << " Domain X: ............... [" << domain_bounds.x().min << ", " << domain_bounds.x().max << "]\n";
26+
oss << " Domain Y: ............... [" << domain_bounds.y().min << ", " << domain_bounds.y().max << "]\n";
27+
oss << " Domain Z: ............... [" << domain_bounds.z().min << ", " << domain_bounds.z().max << "]\n";
28+
oss << " VP: ..................... [" << vp.min << ", " << vp.max << "]\n";
29+
oss << " VS: ..................... [" << vs.min << ", " << vs.max << "]\n";
30+
oss << " V: .......................[" << v.min << ", " << v.max << "]\n";
31+
oss << " Rho: .................... [" << rho.min << ", " << rho.max << "]\n";
32+
oss << " Element Size: ........... [" << element_size.min << ", " << element_size.max << "]\n";
33+
oss << " GLL Distance: ........... [" << gll_distance.min << ", " << gll_distance.max << "]\n";
34+
oss << " Largest Minimum Period: . " << largest_minimum_period << "\n";
35+
oss << " Suggested Time Step: .... " << suggested_time_step << "\n";
36+
return oss.str();
37+
}

‎core/specfem/assembly/info.hpp‎

Lines changed: 70 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,70 @@
1+
#pragma once
2+
3+
#include "specfem/assembly/element_types.hpp"
4+
#include "specfem/assembly/info/impl/bounding_box.hpp"
5+
#include "specfem/assembly/info/impl/bounds.hpp"
6+
#include "specfem/assembly/mesh.hpp"
7+
#include "specfem/assembly/properties.hpp"
8+
#include "specfem/element.hpp"
9+
#include "specfem_setup.hpp"
10+
11+
namespace specfem::assembly {
12+
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 {
31+
constexpr static auto dimension_tag = DimensionTag; ///< Dimension tag
32+
33+
Info() = default;
34+
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+
*/
43+
Info(const specfem::assembly::mesh<dimension_tag> &mesh,
44+
const specfem::assembly::properties<dimension_tag> &properties,
45+
const specfem::assembly::element_types<dimension_tag> &element_types);
46+
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+
*/
65+
std::string string() const;
66+
};
67+
68+
} // namespace specfem::assembly
69+
70+
#include "specfem/assembly/info.tpp"

0 commit comments

Comments
 (0)