[MS] Implement 3D Simulations

CMakeLists.txt

@@ -197,12 +197,12 @@ add_library(
         src/IO/sources.cpp
         src/IO/receivers.cpp
         src/IO/mesh.cpp
-        src/IO/mesh/impl/fortran/read_boundaries.cpp
-        src/IO/mesh/impl/fortran/read_elements.cpp
-        src/IO/mesh/impl/fortran/read_material_properties.cpp
-        src/IO/mesh/impl/fortran/read_mesh_database.cpp
-        src/IO/mesh/impl/fortran/read_interfaces.cpp
-        src/IO/mesh/impl/fortran/read_parameters.cpp
+        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
 )
 
 if (NOT HDF5_CXX_BUILD)
@@ -631,13 +631,28 @@ add_executable(
         src/specfem2d.cpp
 )
 
-
-
 target_link_libraries(
         specfem2d
         execute
 )
 
+add_executable(
+        specfem3d
+        src/specfem3d.cpp
+)
+
+target_link_libraries(
+        specfem3d
+        specfem_mpi
+        Kokkos::kokkos
+        mesh
+        parameter_reader
+        periodic_tasks
+        IO
+        ${BOOST_LIBS}
+        yaml-cpp
+)
+
 # Include tests
 if (BUILD_TESTS)
         message("-- Including tests.")

docs/api/IO/Libraries/mesh/index.rst

@@ -5,26 +5,26 @@ Mesh Reader
 
 
 The below functions are helper functions for the Fortran binary mesh database
-reader. The functions are called in :cpp:func:`specfem::IO::read_mesh`
+reader. The functions are called in :cpp:func:`specfem::IO::read_2d_mesh`
 in sequential order and faciliate the reading of the mesh database.
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_mesh_database_header
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_mesh_database_header
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_coorg_elements
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_coorg_elements
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_properties
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_properties
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_mesh_database_attenuation
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_mesh_database_attenuation
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_material_properties
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_material_properties
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_boundaries
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_boundaries
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_coupled_interfaces
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_coupled_interfaces
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_tangential_elements
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_tangential_elements
 
-.. doxygenfunction:: specfem::IO::mesh::impl::fortran::read_axial_elements
+.. doxygenfunction:: specfem::IO::mesh::impl::fortran::dim2::read_axial_elements
 
 Finally, we add tags to the :cpp:struct:`specfem::mesh::mesh` using the
 :cpp:struct:`specfem::mesh::tags` struct. The description of which can be found

docs/api/IO/index.rst

@@ -3,7 +3,7 @@ Input/Output modules
 ====================
 
 There are three high-level functions that are essential for the SPECFEM++
-codebase: :cpp:func:`specfem::IO::read_mesh`,
+codebase: :cpp:func:`specfem::IO::read_2d_mesh`,
 :cpp:func:`specfem::IO::read_sources`, and
 :cpp:func:`specfem::IO::read_receivers`. These functions are
 used to read the mesh, sources and receivers from disk. The underlying
@@ -32,7 +32,7 @@ are exposed through the following modules:
 Read Mesh
 ---------
 
-.. doxygenfunction:: specfem::IO::read_mesh
+.. doxygenfunction:: specfem::IO::read_2d_mesh
 
 
 Read Sources

examples/dim3/homogeneous_halfspace/DATA/CMTSOLUTION

@@ -0,0 +1,13 @@
+PDE  1999 01 01 00 00 00.00  67000 67000 -25000 4.2 4.2 homog_test
+event name:       homog_test
+time shift:       0.0000
+half duration:    5.0
+latorUTM:       67000.0
+longorUTM:      67000.0
+depth:          30.0
+Mrr:      -7.600000e+27
+Mtt:       7.700000e+27
+Mpp:      -2.000000e+26
+Mrt:      -2.500000e+28
+Mrp:       4.000000e+26
+Mtp:      -2.500000e+27

examples/dim3/homogeneous_halfspace/DATA/meshfem3D_files/Mesh_Par_file

@@ -0,0 +1,109 @@
+#-----------------------------------------------------------
+#
+# Meshing input parameters
+#
+#-----------------------------------------------------------
+
+MESH_A_CHUNK_OF_THE_EARTH        = .false.
+CHUNK_MESH_PAR_FILE              = dummy.txt
+
+# coordinates of mesh block in latitude/longitude and depth in km
+LATITUDE_MIN                    = 0.0
+LATITUDE_MAX                    = 134000.0
+LONGITUDE_MIN                   = 0.0
+LONGITUDE_MAX                   = 134000.0
+DEPTH_BLOCK_KM                  = 60.d0
+UTM_PROJECTION_ZONE             = 11
+SUPPRESS_UTM_PROJECTION         = .true.
+
+# file that contains the interfaces of the model / mesh
+INTERFACES_FILE                 = DATA/meshfem3D_files/interfaces.txt
+
+# file that contains the cavity
+CAVITY_FILE                     = no_cavity.dat
+
+# Number of nodes for 2D and 3D shape functions for hexahedra.
+# We use either 8-node mesh elements (bricks) or 27-node elements.
+# If you use our internal mesher, the only option is 8-node bricks (27-node elements are not supported).
+NGNOD                           = 8
+
+# number of elements at the surface along edges of the mesh at the surface
+# (must be 8 * multiple of NPROC below if mesh is not regular and contains mesh doublings)
+# (must be multiple of NPROC below if mesh is regular)
+NEX_XI                          = 36
+NEX_ETA                         = 36
+
+# number of MPI processors along xi and eta (can be different)
+NPROC_XI                        = 1
+NPROC_ETA                       = 1
+
+#-----------------------------------------------------------
+#
+# Doubling layers
+#
+#-----------------------------------------------------------
+
+# Regular/irregular mesh
+USE_REGULAR_MESH                = .true.
+# Only for irregular meshes, number of doubling layers and their position
+NDOUBLINGS                      = 0
+# NZ_DOUBLING_1 is the parameter to set up if there is only one doubling layer
+# (more doubling entries can be added if needed to match NDOUBLINGS value)
+NZ_DOUBLING_1                   = 40
+NZ_DOUBLING_2                   = 48
+
+#-----------------------------------------------------------
+#
+# Visualization
+#
+#-----------------------------------------------------------
+
+# create mesh files for visualisation or further checking
+CREATE_ABAQUS_FILES             = .false.
+CREATE_DX_FILES                 = .false.
+CREATE_VTK_FILES                = .true.
+
+# stores mesh files as cubit-exported files into directory MESH/ (for single process run)
+SAVE_MESH_AS_CUBIT              = .false.
+
+# path to store the databases files
+LOCAL_PATH                      = ./OUTPUT_FILES/DATABASES_MPI
+
+#-----------------------------------------------------------
+#
+# CPML
+#
+#-----------------------------------------------------------
+
+# CPML perfectly matched absorbing layers
+THICKNESS_OF_X_PML              = 12.3d0
+THICKNESS_OF_Y_PML              = 12.3d0
+THICKNESS_OF_Z_PML              = 12.3d0
+
+#-----------------------------------------------------------
+#
+# Domain materials
+#
+#-----------------------------------------------------------
+
+# number of materials
+NMATERIALS                      = 1
+# define the different materials in the model as:
+# #material_id  #rho  #vp  #vs  #Q_Kappa  #Q_mu  #anisotropy_flag  #domain_id
+#     Q_Kappa          : Q_Kappa attenuation quality factor
+#     Q_mu             : Q_mu attenuation quality factor
+#     anisotropy_flag  : 0 = no anisotropy / 1,2,... check the implementation in file aniso_model.f90
+#     domain_id        : 1 = acoustic / 2 = elastic
+1   2300.0   2800.0   1500.0   2444.4    300.0 0 2
+
+#-----------------------------------------------------------
+#
+# Domain regions
+#
+#-----------------------------------------------------------
+
+# number of regions
+NREGIONS                        = 1
+# define the different regions of the model as :
+#NEX_XI_BEGIN  #NEX_XI_END  #NEX_ETA_BEGIN  #NEX_ETA_END  #NZ_BEGIN #NZ_END  #material_id
+1             36             1              36        1        16              1

examples/dim3/homogeneous_halfspace/DATA/meshfem3D_files/interface1.txt

@@ -0,0 +1,4 @@
+0
+0
+0
+0

examples/dim3/homogeneous_halfspace/DATA/meshfem3D_files/interfaces.txt

@@ -0,0 +1,16 @@
+# number of interfaces
+ 1
+#
+# We describe each interface below, structured as a 2D-grid, with several parameters :
+# number of points along XI and ETA, minimal XI ETA coordinates
+# and spacing between points which must be constant.
+# Then the records contain the Z coordinates of the NXI x NETA points.
+#
+# interface number 1 (topography, top of the mesh)
+ .true. 2 2 0.0 0.0 1000.0 1000.0
+ DATA/meshfem3D_files/interface1.txt
+#
+# for each layer, we give the number of spectral elements in the vertical direction
+#
+# layer number 1 (top layer)
+ 16