Merge pull request #462 from PrincetonUniversity/issue-461

Issue 461 -- Implement Parameter Reader

CMakeLists.txt

@@ -197,12 +197,15 @@ add_library(
         src/IO/sources.cpp
         src/IO/receivers.cpp
         src/IO/mesh.cpp
+        # Fortran 2D part
         src/IO/mesh/impl/fortran/dim2/read_boundaries.cpp
         src/IO/mesh/impl/fortran/dim2/read_elements.cpp
         src/IO/mesh/impl/fortran/dim2/read_material_properties.cpp
         src/IO/mesh/impl/fortran/dim2/read_mesh_database.cpp
         src/IO/mesh/impl/fortran/dim2/read_interfaces.cpp
         src/IO/mesh/impl/fortran/dim2/read_parameters.cpp
+        # Fortran 3D part
+        src/IO/mesh/impl/fortran/dim3/read_parameters.cpp
 )
 
 if (NOT HDF5_CXX_BUILD)

CMakePresets.json

@@ -0,0 +1,26 @@
+{
+  "version": 6,
+  "configurePresets": [
+    {
+      "name": "debug",
+      "displayName": "Debug (Serial)",
+      "binaryDir": "build/debug",
+      "cacheVariables": {
+        "CMAKE_BUILD_TYPE": "Debug",
+        "BUILD_TESTS": "ON",
+        "ENABLE_SIMD": "OFF",
+        "Kokkos_ARCH_NATIVE": "ON",
+        "Kokkos_ENABLE_AGGRESSIVE_VECTORIZATION": "ON",
+        "Kokkos_ENABLE_ATOMICS_BYPASS": "ON"
+      }
+    }
+  ],
+  "buildPresets": [
+    {
+      "name": "debug",
+      "configurePreset": "debug",
+      "targets": ["all"],
+      "jobs": 8
+    }
+  ]
+}

examples/dim3/homogeneous_halfspace/specfem_config.yaml

@@ -45,6 +45,7 @@ parameters:
   ## databases
   databases:
     mesh-database: "OUTPUT_FILES/DATABASES_MPI/proc000000_external_mesh.bin"
+    mesh-parameters: "OUTPUT_FILES/DATABASES_MPI/mesh_parameters.bin"
 
   ## sources
   sources: "/Users/lsawade/SPECFEMPP/examples/homogeneous-medium-flat-topography/source.yaml"

fortran/meshfem3d/generate_databases/save_parameters.f90

@@ -58,7 +58,6 @@ subroutine save_parameters()
       WRITE(IOUT) APPROXIMATE_OCEAN_LOAD
       WRITE(IOUT) USE_MESH_COLORING_GPU
 
-
       ! Write test parameter
       itest = 9998
       WRITE(IOUT) itest
@@ -72,32 +71,32 @@ subroutine save_parameters()
       WRITE(IOUT) NGLLSQUARE
 
       ! Write test parameter
-      WRITE(IOUT) itest
       itest = 9997
+      WRITE(IOUT) itest
 
       WRITE(IOUT) NSPEC_AB ! nspec
       WRITE(IOUT) NSPEC_PORO ! nspec_poro == nspec, if POROELASTIC_SIMULATION
       WRITE(IOUT) NGLOB_AB ! nglob
       WRITE(IOUT) NGLOB_OCEAN
-      WRITE(IOUT) NSPEC2D_BOTTOM
-      WRITE(IOUT) NSPEC2D_TOP
-      WRITE(IOUT) nspec_irregular
 
       ! Write test parameter
-      WRITE(IOUT) itest
       itest = 9996
+      WRITE(IOUT) itest
 
-      WRITE(IOUT) num_neighbors_all
-      WRITE(IOUT) nfaces_surface
-      WRITE(IOUT) nspec2D_xmin
+      WRITE(IOUT) NSPEC2D_BOTTOM
+      WRITE(IOUT) NSPEC2D_TOP
       WRITE(IOUT) nspec2D_xmin
+      WRITE(IOUT) nspec2D_xmax
       WRITE(IOUT) nspec2D_ymin
-      WRITE(IOUT) nspec2D_ymin
+      WRITE(IOUT) nspec2D_ymax
+      WRITE(IOUT) nspec_irregular
 
       ! Write test parameter
-      WRITE(IOUT) itest
       itest = 9995
+      WRITE(IOUT) itest
 
+      WRITE(IOUT) num_neighbors_all
+      WRITE(IOUT) nfaces_surface
       WRITE(IOUT) num_abs_boundary_faces
       WRITE(IOUT) num_free_surface_faces
       WRITE(IOUT) num_coupling_ac_el_faces
@@ -108,8 +107,8 @@ subroutine save_parameters()
       WRITE(IOUT) max_nibool_interfaces_ext_mesh
 
       ! Write test parameter
-      WRITE(IOUT) itest
       itest = 9994
+      WRITE(IOUT) itest
 
       WRITE(IOUT) nspec_inner_acoustic
       WRITE(IOUT) nspec_outer_acoustic
@@ -119,8 +118,8 @@ subroutine save_parameters()
       WRITE(IOUT) nspec_outer_poroelastic
 
       ! Write test parameter
-      WRITE(IOUT) itest
       itest = 9993
+      WRITE(IOUT) itest
 
       WRITE(IOUT) num_phase_ispec_acoustic
       WRITE(IOUT) num_phase_ispec_elastic
@@ -131,8 +130,8 @@ subroutine save_parameters()
       WRITE(IOUT) num_colors_outer_elastic
 
       ! Write test parameter
-      WRITE(IOUT) itest
       itest = 9992
+      WRITE(IOUT) itest
   endif
 
 end subroutine save_parameters

include/IO/interface.hpp

@@ -25,7 +25,9 @@ specfem::mesh::mesh<specfem::dimension::type::dim2>
 read_2d_mesh(const std::string filename, const specfem::MPI::MPI *mpi);
 
 specfem::mesh::mesh<specfem::dimension::type::dim3>
-read_3d_mesh(const std::string filename, const specfem::MPI::MPI *mpi);
+read_3d_mesh(const std::string mesh_parameters_file,
+             const std::string mesh_databases_file,
+             const specfem::MPI::MPI *mpi);
 /**
  * @brief Read station file
  *

include/IO/mesh/impl/fortran/dim3/read_parameters.hpp

@@ -0,0 +1,30 @@
+#pragma once
+
+#include "IO/fortranio/interface.hpp"
+#include "mesh/parameters/parameters.hpp"
+#include "specfem_mpi/interface.hpp"
+
+namespace specfem {
+namespace IO {
+namespace mesh {
+namespace impl {
+namespace fortran {
+namespace dim3 {
+
+/*
+ * @brief Read paramters from 3D mesh database
+ *
+ * @param stream Input stream
+ * @param mpi MPI object
+ * @return specfem::mesh::parameters<specfem::dimension::type::dim2> Mesh
+ * parameters
+ */
+specfem::mesh::parameters<specfem::dimension::type::dim3>
+read_mesh_parameters(std::ifstream &stream, const specfem::MPI::MPI *mpi);
+
+} // namespace dim3
+} // namespace fortran
+} // namespace impl
+} // namespace mesh
+} // namespace IO
+} // namespace specfem

include/mesh/mesh.hpp

@@ -5,6 +5,7 @@
 #include "coupled_interfaces/coupled_interfaces.hpp"
 #include "elements/axial_elements.hpp"
 #include "elements/tangential_elements.hpp"
+#include "enumerations/dimension.hpp"
 #include "enumerations/interface.hpp"
 #include "materials/materials.hpp"
 #include "mesh/tags/tags.hpp"
@@ -104,6 +105,9 @@ template <> struct mesh<specfem::dimension::type::dim3> {
   constexpr static auto dimension =
       specfem::dimension::type::dim3; ///< Dimension
 
+  // Struct to store all the mesh parameter
+  specfem::mesh::parameters<dimension> parameters;
+
   // int npgeo; ///< Total number of spectral element control nodes
   // int nspec; ///< Total number of spectral elements
   // int nproc; ///< Total number of processors

include/mesh/parameters/parameters.hpp

@@ -76,6 +76,186 @@ template <> struct parameters<specfem::dimension::type::dim2> {
         nelem_on_the_axis(nelem_on_the_axis),
         plot_lowerleft_corner_only(plot_lowerleft_corner_only){};
 };
+
+/**
+ * @brief Template specialization for 3D mesh parameters
+ */
+template <> struct parameters<specfem::dimension::type::dim3> {
+  constexpr static auto dimension =
+      specfem::dimension::type::dim3; ///< Dimension
+                                      ///< type
+
+  // Flags
+  bool acoustic_simulation;    ///< Flag for acoustic simulation
+  bool elastic_simulation;     ///< Flag for elastic simulation
+  bool poroelastic_simulation; ///< Flag for poroelastic simulation
+  bool anisotropy;             ///< Flag for anisotropy
+  bool stacey_abc;             ///< Flag for Stacey absorbing boundary c.
+  bool pml_abc;                ///< Flag for PML absorbing boundary c.
+  bool approximate_ocean_load; ///< Flag for approx. ocean load
+  bool use_mesh_coloring;      ///< Flag for mesh coloring
+
+  // Integer Parameters: Dimensions and GLL layouts
+  int ndim;       ///< Number of dimensions
+  int ngllx;      ///< Number of GLL points in x
+  int nglly;      ///< Number of GLL points in y
+  int ngllz;      ///< Number of GLL points in z
+  int ngllsquare; ///< Number of GLL points in square
+
+  // Integer Parameters: Elements/Nodes
+  int nspec;           ///< Number of spectral elements (SEs)
+  int nspec_poro;      ///< Number of poroelastic SEs
+  int nglob;           ///< Number of global nodes
+  int nglob_ocean;     ///< Number of global ocean nodes
+  int nspec2D_bottom;  ///< Number of 2D SEs at the bottom
+  int nspec2D_top;     ///< Number of 2D SEs at the top
+  int nspec2D_xmin;    ///< Number of 2D SEs at the left
+  int nspec2D_xmax;    ///< Number of 2D SEs at the right
+  int nspec2D_ymin;    ///< Number of 2D SEs at the front
+  int nspec2D_ymax;    ///< Number of 2D SEs at the back
+  int nspec_irregular; ///< Number of irregular SEs
+
+  // Integer Parameters: Mesh
+  int num_neighbors;                  ///< Number of neighbors
+  int nfaces_surface;                 ///< Number of faces on the surface
+  int num_abs_boundary_faces;         ///< Number of absorbing boundary faces
+  int num_free_surface_faces;         ///< Number of free surface faces
+  int num_coupling_ac_el_faces;       ///< Number of acoustic-elastic faces
+  int num_coupling_ac_po_faces;       ///< Number of acoustic-poroelastic faces
+  int num_coupling_el_po_faces;       ///< Number of elastic-poroelastic faces
+  int num_coupling_po_el_faces;       ///< Number of poroelastic-elastic faces
+  int num_interfaces_ext_mesh;        ///< Number of external mesh interfaces
+  int max_nibool_interfaces_ext_mesh; ///< Maximum number of interfaces
+
+  // Integer Parameters: Nspec Inner/Outer
+  int nspec_inner_acoustic;    ///< Number of inner acoustic SEs
+  int nspec_outer_acoustic;    ///< Number of outer acoustic SEs
+  int nspec_inner_elastic;     ///< Number of inner elastic SEs
+  int nspec_outer_elastic;     ///< Number of outer elastic SEs
+  int nspec_inner_poroelastic; ///< Number of inner poroelastic SEs
+  int nspec_outer_poroelastic; ///< Number of outer poroelastic SEs
+
+  // Integer Parameters: coloring
+  int num_phase_ispec_acoustic;
+  int num_phase_ispec_elastic;
+  int num_phase_ispec_poroelastic;
+  int num_colors_inner_acoustic;
+  int num_colors_outer_acoustic;
+  int num_colors_inner_elastic;
+  int num_colors_outer_elastic;
+
+  /**
+   * @brief Default constructor
+   *
+   */
+  parameters(){};
+
+  /** Constructor
+   * @param acoustic_simulation Flag for acoustic simulation
+   * @param elastic_simulation Flag for elastic simulation
+   * @param poroelastic_simulation Flag for poroelastic simulation
+   * @param anisotropy Flag for anisotropy
+   * @param stacey_abc Flag for Stacey absorbing boundary c.
+   * @param pml_abc Flag for PML absorbing boundary c.
+   * @param approximate_ocean_load Flag for approx. ocean load
+   * @param use_mesh_coloring Flag for mesh coloring
+   *
+   *
+   *
+   * @param nspec Number of spectral elements (SEs)
+   * @param nspec_poro Number of poroelastic SEs
+   * @param nglob Number of global nodes
+   * @param nglob_ocean Number of global ocean nodes
+   * @param nspec2D_bottom Number of 2D SEs at the bottom
+   * @param nspec2D_top Number of 2D SEs at the top
+   * @param nspec2D_xmin Number of 2D SEs at the left
+   * @param nspec2D_xmax Number of 2D SEs at the right
+   * @param nspec2D_ymin Number of 2D SEs at the front
+   * @param nspec2D_ymax Number of 2D SEs at the back
+   * @param nspec_irregular Number of irregular SEs
+   *
+   * @param num_neighbors Number of neighbors
+   * @param nfaces_surface Number of faces on the surface
+   * @param num_abs_boundary_faces Number of absorbing boundary faces
+   * @param num_free_surface_faces Number of free surface faces
+   * @param num_coupling_ac_el_faces Number of acoustic-elastic faces
+   * @param num_coupling_ac_po_faces Number of acoustic-poroelastic faces
+   * @param num_coupling_el_po_faces Number of elastic-poroelastic faces
+   * @param num_coupling_po_el_faces Number of poroelastic-elastic faces
+   * @param num_interfaces_ext_mesh Number of external mesh interfaces
+   * @param max_nibool_interfaces_ext_mesh Maximum number of interfaces
+   *
+   * @param nspec_inner_acoustic Number of inner acoustic SEs
+   * @param nspec_outer_acoustic Number of outer acoustic SEs
+   * @param nspec_inner_elastic Number of inner elastic SEs
+   * @param nspec_outer_elastic Number of outer elastic SEs
+   * @param nspec_inner_poroelastic Number of inner poroelastic SEs
+   * @param nspec_outer_poroelastic Number of outer poroelastic SEs
+   *
+   * @param num_phase_ispec_acoustic
+   * @param num_phase_ispec_elastic
+   * @param num_phase_ispec_poroelastic
+   * @param num_colors_inner_acoustic
+   * @param num_colors_outer_acoustic
+   * @param num_colors_inner_elastic
+   * @param num_colors_outer_elastic
+   *
+   */
+  parameters(
+      const bool acoustic_simulation, const bool elastic_simulation,
+      const bool poroelastic_simulation, const bool anisotropy,
+      const bool stacey_abc, const bool pml_abc,
+      const bool approximate_ocean_load, const bool use_mesh_coloring,
+      const int nspec, const int nspec_poro, const int nglob,
+      const int nglob_ocean, const int nspec2D_bottom, const int nspec2D_top,
+      const int nspec2D_xmin, const int nspec2D_xmax, const int nspec2D_ymin,
+      const int nspec2D_ymax, const int nspec_irregular,
+      const int num_neighbors, const int nfaces_surface,
+      const int num_abs_boundary_faces, const int num_free_surface_faces,
+      const int num_coupling_ac_el_faces, const int num_coupling_ac_po_faces,
+      const int num_coupling_el_po_faces, const int num_coupling_po_el_faces,
+      const int num_interfaces_ext_mesh,
+      const int max_nibool_interfaces_ext_mesh, const int nspec_inner_acoustic,
+      const int nspec_outer_acoustic, const int nspec_inner_elastic,
+      const int nspec_outer_elastic, const int nspec_inner_poroelastic,
+      const int nspec_outer_poroelastic, const int num_phase_ispec_acoustic,
+      const int num_phase_ispec_elastic, const int num_phase_ispec_poroelastic,
+      const int num_colors_inner_acoustic, const int num_colors_outer_acoustic,
+      const int num_colors_inner_elastic, const int num_colors_outer_elastic)
+      : acoustic_simulation(acoustic_simulation),
+        elastic_simulation(elastic_simulation),
+        poroelastic_simulation(poroelastic_simulation), anisotropy(anisotropy),
+        stacey_abc(stacey_abc), pml_abc(pml_abc),
+        approximate_ocean_load(approximate_ocean_load),
+        use_mesh_coloring(use_mesh_coloring), nspec(nspec),
+        nspec_poro(nspec_poro), nglob(nglob), nglob_ocean(nglob_ocean),
+        nspec2D_bottom(nspec2D_bottom), nspec2D_top(nspec2D_top),
+        nspec2D_xmin(nspec2D_xmin), nspec2D_xmax(nspec2D_xmax),
+        nspec2D_ymin(nspec2D_ymin), nspec2D_ymax(nspec2D_ymax),
+        nspec_irregular(nspec_irregular), num_neighbors(num_neighbors),
+        nfaces_surface(nfaces_surface),
+        num_abs_boundary_faces(num_abs_boundary_faces),
+        num_free_surface_faces(num_free_surface_faces),
+        num_coupling_ac_el_faces(num_coupling_ac_el_faces),
+        num_coupling_ac_po_faces(num_coupling_ac_po_faces),
+        num_coupling_el_po_faces(num_coupling_el_po_faces),
+        num_coupling_po_el_faces(num_coupling_po_el_faces),
+        num_interfaces_ext_mesh(num_interfaces_ext_mesh),
+        max_nibool_interfaces_ext_mesh(max_nibool_interfaces_ext_mesh),
+        nspec_inner_acoustic(nspec_inner_acoustic),
+        nspec_outer_acoustic(nspec_outer_acoustic),
+        nspec_inner_elastic(nspec_inner_elastic),
+        nspec_outer_elastic(nspec_outer_elastic),
+        nspec_inner_poroelastic(nspec_inner_poroelastic),
+        nspec_outer_poroelastic(nspec_outer_poroelastic),
+        num_phase_ispec_acoustic(num_phase_ispec_acoustic),
+        num_phase_ispec_elastic(num_phase_ispec_elastic),
+        num_phase_ispec_poroelastic(num_phase_ispec_poroelastic),
+        num_colors_inner_acoustic(num_colors_inner_acoustic),
+        num_colors_outer_acoustic(num_colors_outer_acoustic),
+        num_colors_inner_elastic(num_colors_inner_elastic),
+        num_colors_outer_elastic(num_colors_outer_elastic){};
+};
 } // namespace mesh
 } // namespace specfem