[cig-commits] [commit] devel: framework for the future broadcast of the mesh and model databases in the case of simultaneous runs for the same mesh and model (but different source and/or receivers) (01fb7e7)

cig_noreply at geodynamics.org cig_noreply at geodynamics.org
Tue Sep 23 15:57:37 PDT 2014


Repository : https://github.com/geodynamics/specfem3d

On branch  : devel
Link       : https://github.com/geodynamics/specfem3d/compare/428a00421372cb9e3236b3adf6ab58cb8c3c88b5...f69391023886e2173d4e185197e649e04b9580f0

>---------------------------------------------------------------

commit 01fb7e752006ce2f94e2c6c0eba6c9bcab5806fa
Author: Dimitri Komatitsch <komatitsch at lma.cnrs-mrs.fr>
Date:   Tue Sep 23 23:22:15 2014 +0200

    framework for the future broadcast of the mesh and model databases in the case of simultaneous runs for the same mesh and model (but different source and/or receivers)


>---------------------------------------------------------------

01fb7e752006ce2f94e2c6c0eba6c9bcab5806fa
 src/generate_databases/model_tomography.f90 |   2 +-
 src/shared/constants_mod.f90                |   4 +
 src/shared/parallel.f90                     | 161 +++++++++-
 src/shared/read_parameter_file.f90          |   4 +-
 src/shared/serial.f90                       |  29 ++
 src/specfem3D/read_mesh_databases.F90       | 462 ++++++++++++++++------------
 6 files changed, 455 insertions(+), 207 deletions(-)

diff --git a/src/generate_databases/model_tomography.f90 b/src/generate_databases/model_tomography.f90
index 8ec345c..2c31041 100644
--- a/src/generate_databases/model_tomography.f90
+++ b/src/generate_databases/model_tomography.f90
@@ -101,7 +101,7 @@
   ! only master reads in model file
   !if(myrank == 0) call read_model_tomography()
   ! broadcast the information read on the master to the nodes, e.g.
-  !call bcast_all_i(nrecord,1)
+  !call bcast_all_one_i(nrecord)
 
   !if( myrank /= 0 ) then
   ! allocate( vp_tomography(1:nrecord) ,stat=ier)
diff --git a/src/shared/constants_mod.f90 b/src/shared/constants_mod.f90
index 99ff0de..a277a77 100644
--- a/src/shared/constants_mod.f90
+++ b/src/shared/constants_mod.f90
@@ -38,5 +38,9 @@ module constants
   ! if NUMBER_OF_SIMULTANEOUS_RUNS > 1
   character(len=MAX_STRING_LEN) :: OUTPUT_FILES_PATH = OUTPUT_FILES_PATH_BASE
 
+  ! in the case of simultaneous runs for the same mesh and model, see who reads the mesh and the model and broadcasts it to others
+  ! we put default values here
+  logical :: I_should_read_the_database = .true., I_should_broadcast_the_database = .false.
+
 end module constants
 
diff --git a/src/shared/parallel.f90 b/src/shared/parallel.f90
index c131c69..2286a1b 100644
--- a/src/shared/parallel.f90
+++ b/src/shared/parallel.f90
@@ -57,7 +57,7 @@ module my_mpi
 
   implicit none
 
-  integer :: my_local_mpi_comm_world
+  integer :: my_local_mpi_comm_world, my_local_mpi_comm_for_bcast
 
 end module my_mpi
 
@@ -110,6 +110,24 @@ end module my_mpi
   end subroutine synchronize_all
 
 !
+!---- broadcast using the default communicator for the whole run
+!
+
+  subroutine bcast_all_one_i(buffer)
+
+  use my_mpi
+
+  implicit none
+
+  integer :: buffer
+
+  integer ier
+
+  call MPI_BCAST(buffer,1,MPI_INTEGER,0,my_local_mpi_comm_world,ier)
+
+  end subroutine bcast_all_one_i
+
+!
 !----
 !
 
@@ -188,6 +206,102 @@ end module my_mpi
 
   end subroutine bcast_all_r
 
+!
+!---- broadcast using the communicator to send the mesh and model to other simultaneous runs
+!
+
+  subroutine bcast_all_one_i_for_database(buffer)
+
+  use my_mpi
+
+  implicit none
+
+  integer :: buffer
+
+  integer ier
+
+  call MPI_BCAST(buffer,1,MPI_INTEGER,0,my_local_mpi_comm_for_bcast,ier)
+
+  end subroutine bcast_all_one_i_for_database
+
+!
+!----
+!
+
+  subroutine bcast_all_i_for_database(buffer, countval)
+
+  use my_mpi
+
+  implicit none
+
+  integer countval
+  integer, dimension(countval) :: buffer
+
+  integer ier
+
+  call MPI_BCAST(buffer,countval,MPI_INTEGER,0,my_local_mpi_comm_for_bcast,ier)
+
+  end subroutine bcast_all_i_for_database
+
+!
+!----
+!
+
+  subroutine bcast_all_cr_for_database(buffer, countval)
+
+  use my_mpi
+  use constants,only: CUSTOM_REAL
+
+  implicit none
+
+  include "precision.h"
+
+  integer countval
+  real(kind=CUSTOM_REAL), dimension(countval) :: buffer
+
+  integer ier
+
+  call MPI_BCAST(buffer,countval,CUSTOM_MPI_TYPE,0,my_local_mpi_comm_for_bcast,ier)
+
+  end subroutine bcast_all_cr_for_database
+
+!
+!----
+!
+
+  subroutine bcast_all_dp_for_database(buffer, countval)
+
+  use my_mpi
+
+  implicit none
+
+  integer countval
+  double precision, dimension(countval) :: buffer
+
+  integer ier
+
+  call MPI_BCAST(buffer,countval,MPI_DOUBLE_PRECISION,0,my_local_mpi_comm_for_bcast,ier)
+
+  end subroutine bcast_all_dp_for_database
+
+!
+!----
+!
+
+  subroutine bcast_all_r_for_database(buffer, countval)
+
+  use my_mpi
+
+  implicit none
+
+  integer countval
+  real, dimension(countval) :: buffer
+
+  integer ier
+
+  call MPI_BCAST(buffer,countval,MPI_REAL,0,my_local_mpi_comm_for_bcast,ier)
+
+  end subroutine bcast_all_r_for_database
 
 !
 !----
@@ -1154,12 +1268,11 @@ end module my_mpi
 
   use my_mpi
   use constants,only: MAX_STRING_LEN,NUMBER_OF_SIMULTANEOUS_RUNS,OUTPUT_FILES_PATH, &
-    IMAIN,ISTANDARD_OUTPUT, &
-    mygroup
+    IMAIN,ISTANDARD_OUTPUT,mygroup,BROADCAST_SAME_MESH_AND_MODEL,I_should_read_the_database,I_should_broadcast_the_database
 
   implicit none
 
-  integer :: sizeval,myrank,ier,key
+  integer :: sizeval,myrank,ier,key,my_group_for_bcast,my_local_rank_for_bcast,NPROC
 
   character(len=MAX_STRING_LEN) :: path_to_add
 
@@ -1176,8 +1289,14 @@ end module my_mpi
 
     my_local_mpi_comm_world = MPI_COMM_WORLD
 
+! no broadcast of the mesh and model databases to other runs in that case
+    my_group_for_bcast = 0
+    my_local_mpi_comm_for_bcast = MPI_COMM_NULL
+
   else
 
+!--- create a subcommunicator for each independent run
+
     call MPI_COMM_RANK(MPI_COMM_WORLD,myrank,ier)
 !   create the different groups of processes, one for each independent run
     mygroup = mod(myrank,NUMBER_OF_SIMULTANEOUS_RUNS)
@@ -1192,6 +1311,33 @@ end module my_mpi
     write(path_to_add,"('run',i4.4,'/')") mygroup + 1
     OUTPUT_FILES_PATH = path_to_add(1:len_trim(path_to_add))//OUTPUT_FILES_PATH(1:len_trim(OUTPUT_FILES_PATH))
 
+!--- create a subcommunicator to broadcast the identical mesh and model databases if needed
+    if(BROADCAST_SAME_MESH_AND_MODEL) then
+
+      call MPI_COMM_RANK(MPI_COMM_WORLD,myrank,ier)
+!     to broadcast the model, split along similar ranks per run instead
+      NPROC = sizeval / NUMBER_OF_SIMULTANEOUS_RUNS
+      my_group_for_bcast = mod(myrank,NPROC)
+      key = myrank
+      if(my_group_for_bcast < 0 .or. my_group_for_bcast > NPROC-1) stop 'invalid value of my_group_for_bcast'
+
+!     build the sub-communicators
+      call MPI_COMM_SPLIT(MPI_COMM_WORLD, my_group_for_bcast, key, my_local_mpi_comm_for_bcast, ier)
+      if(ier /= 0) stop 'error while trying to create the sub-communicators'
+
+!     see if that process will need to read the mesh and model database and then broadcast it to others
+      call MPI_COMM_RANK(my_local_mpi_comm_for_bcast,my_local_rank_for_bcast,ier)
+      if(my_local_rank_for_bcast > 0) I_should_read_the_database = .false.
+      if(my_local_rank_for_bcast == 0) I_should_broadcast_the_database = .true.
+
+    else
+
+! no broadcast of the mesh and model databases to other runs in that case
+      my_group_for_bcast = 0
+      my_local_mpi_comm_for_bcast = MPI_COMM_NULL
+
+    endif
+
   endif
 
   end subroutine world_split
@@ -1204,13 +1350,16 @@ end module my_mpi
   subroutine world_unsplit()
 
   use my_mpi
-  use constants,only: NUMBER_OF_SIMULTANEOUS_RUNS
+  use constants,only: NUMBER_OF_SIMULTANEOUS_RUNS,BROADCAST_SAME_MESH_AND_MODEL
 
   implicit none
 
   integer :: ier
 
-  if(NUMBER_OF_SIMULTANEOUS_RUNS > 1) call MPI_COMM_FREE(my_local_mpi_comm_world,ier)
+  if(NUMBER_OF_SIMULTANEOUS_RUNS > 1) then
+    call MPI_COMM_FREE(my_local_mpi_comm_world,ier)
+    if(BROADCAST_SAME_MESH_AND_MODEL) call MPI_COMM_FREE(my_local_mpi_comm_for_bcast,ier)
+  endif
 
   end subroutine world_unsplit
 
diff --git a/src/shared/read_parameter_file.f90 b/src/shared/read_parameter_file.f90
index af2c070..663acb8 100644
--- a/src/shared/read_parameter_file.f90
+++ b/src/shared/read_parameter_file.f90
@@ -441,6 +441,8 @@ subroutine read_adios_parameters(ADIOS_ENABLED, ADIOS_FOR_DATABASES,       &
   endif
   call close_parameter_file()
 
-end subroutine read_adios_parameters
+  if(NUMBER_OF_SIMULTANEOUS_RUNS > 1 .and. BROADCAST_SAME_MESH_AND_MODEL .and. ADIOS_ENABLED) &
+         stop 'ADIOS not yet supported by option BROADCAST_SAME_MESH_AND_MODEL'
 
+end subroutine read_adios_parameters
 
diff --git a/src/shared/serial.f90 b/src/shared/serial.f90
index 48bcdf4..2dc46ba 100644
--- a/src/shared/serial.f90
+++ b/src/shared/serial.f90
@@ -61,6 +61,21 @@
 !----
 !
 
+  subroutine bcast_all_one_i(buffer)
+
+  use unused_mod
+  implicit none
+
+  integer :: buffer
+
+  unused_i4 = buffer
+
+  end subroutine bcast_all_one_i
+
+!
+!----
+!
+
   subroutine bcast_all_i(buffer, count)
 
   use unused_mod
@@ -123,6 +138,20 @@
 
   end subroutine bcast_all_r
 
+!
+!----
+!
+
+  subroutine bcast_all_one_i_for_database(buffer)
+
+  use unused_mod
+  implicit none
+
+  integer :: buffer
+
+  unused_i4 = buffer
+
+  end subroutine bcast_all_one_i_for_database
 
 !
 !----
diff --git a/src/specfem3D/read_mesh_databases.F90 b/src/specfem3D/read_mesh_databases.F90
index 41f30ff..94603b7 100644
--- a/src/specfem3D/read_mesh_databases.F90
+++ b/src/specfem3D/read_mesh_databases.F90
@@ -45,40 +45,44 @@
 
 ! info about external mesh simulation
   call create_name_database(prname,myrank,LOCAL_PATH)
-  open(unit=27,file=prname(1:len_trim(prname))//'external_mesh.bin',status='old',&
+  if(I_should_read_the_database) then
+    open(unit=27,file=prname(1:len_trim(prname))//'external_mesh.bin',status='old',&
        action='read',form='unformatted',iostat=ier)
-  if( ier /= 0 ) then
-    print*,'error: could not open database '
-    print*,'path: ',prname(1:len_trim(prname))//'external_mesh.bin'
-    call exit_mpi(myrank,'error opening database')
+    if( ier /= 0 ) then
+      print*,'error: could not open database '
+      print*,'path: ',prname(1:len_trim(prname))//'external_mesh.bin'
+      call exit_mpi(myrank,'error opening database')
+    endif
   endif
 
-  read(27) NSPEC_AB
-  read(27) NGLOB_AB
-
-  read(27) ibool
-
-  read(27) xstore
-  read(27) ystore
-  read(27) zstore
-
-  read(27) xix
-  read(27) xiy
-  read(27) xiz
-  read(27) etax
-  read(27) etay
-  read(27) etaz
-  read(27) gammax
-  read(27) gammay
-  read(27) gammaz
-  read(27) jacobian
-
-  read(27) kappastore
-  read(27) mustore
-
-  read(27) ispec_is_acoustic
-  read(27) ispec_is_elastic
-  read(27) ispec_is_poroelastic
+  if(I_should_read_the_database) then
+    read(27) NSPEC_AB
+    read(27) NGLOB_AB
+
+    read(27) ibool
+
+    read(27) xstore
+    read(27) ystore
+    read(27) zstore
+
+    read(27) xix
+    read(27) xiy
+    read(27) xiz
+    read(27) etax
+    read(27) etay
+    read(27) etaz
+    read(27) gammax
+    read(27) gammay
+    read(27) gammaz
+    read(27) jacobian
+
+    read(27) kappastore
+    read(27) mustore
+
+    read(27) ispec_is_acoustic
+    read(27) ispec_is_elastic
+    read(27) ispec_is_poroelastic
+  endif
 
   ! acoustic
   ! number of acoustic elements in this partition
@@ -100,7 +104,7 @@
     ! mass matrix, density
     allocate(rmass_acoustic(NGLOB_AB),stat=ier)
     if( ier /= 0 ) stop 'error allocating array rmass_acoustic'
-    read(27) rmass_acoustic
+    if(I_should_read_the_database) read(27) rmass_acoustic
 
     ! initializes mass matrix contribution
     allocate(rmassz_acoustic(NGLOB_AB),stat=ier)
@@ -112,7 +116,7 @@
 ! thus we now allocate it and read it in all cases (whether the simulation is acoustic, elastic, or acoustic/elastic)
   allocate(rhostore(NGLLX,NGLLY,NGLLZ,NSPEC_AB),stat=ier)
   if( ier /= 0 ) stop 'error allocating array rhostore'
-  read(27) rhostore
+  if(I_should_read_the_database) read(27) rhostore
 
   ! elastic
   ! number of elastic elements in this partition
@@ -209,14 +213,14 @@
     if( ier /= 0 ) stop 'error allocating array one_minus_sum_beta_kappa etc.'
 
     ! reads mass matrices
-    read(27,iostat=ier) rmass
+    if(I_should_read_the_database) read(27,iostat=ier) rmass
     if( ier /= 0 ) stop 'error reading in array rmass'
 
     if( APPROXIMATE_OCEAN_LOAD ) then
       ! ocean mass matrix
       allocate(rmass_ocean_load(NGLOB_AB),stat=ier)
       if( ier /= 0 ) stop 'error allocating array rmass_ocean_load'
-      read(27) rmass_ocean_load
+      if(I_should_read_the_database) read(27) rmass_ocean_load
     else
       ! dummy allocation
       allocate(rmass_ocean_load(1),stat=ier)
@@ -224,9 +228,9 @@
     endif
 
     !pll material parameters for stacey conditions
-    read(27,iostat=ier) rho_vp
+    if(I_should_read_the_database) read(27,iostat=ier) rho_vp
     if( ier /= 0 ) stop 'error reading in array rho_vp'
-    read(27,iostat=ier) rho_vs
+    if(I_should_read_the_database) read(27,iostat=ier) rho_vs
     if( ier /= 0 ) stop 'error reading in array rho_vs'
 
   else
@@ -288,35 +292,38 @@
              epsilonw_trace_over_3(NGLLX,NGLLY,NGLLZ,NSPEC_ADJOINT),stat=ier)
     if( ier /= 0 ) stop 'error allocating array epsilons_trace_over_3 etc.'
 
-
-    read(27) rmass_solid_poroelastic
-    read(27) rmass_fluid_poroelastic
-    read(27) rhoarraystore
-    read(27) kappaarraystore
-    read(27) etastore
-    read(27) tortstore
-    read(27) permstore
-    read(27) phistore
-    read(27) rho_vpI
-    read(27) rho_vpII
-    read(27) rho_vsI
+    if(I_should_read_the_database) then
+      read(27) rmass_solid_poroelastic
+      read(27) rmass_fluid_poroelastic
+      read(27) rhoarraystore
+      read(27) kappaarraystore
+      read(27) etastore
+      read(27) tortstore
+      read(27) permstore
+      read(27) phistore
+      read(27) rho_vpI
+      read(27) rho_vpII
+      read(27) rho_vsI
+    endif
   endif
 
   ! checks simulation types are valid
   if( (.not. ACOUSTIC_SIMULATION ) .and. &
      (.not. ELASTIC_SIMULATION ) .and. &
      (.not. POROELASTIC_SIMULATION ) ) then
-    close(27)
+    if(I_should_read_the_database) close(27)
     call exit_mpi(myrank,'error no simulation type defined')
   endif
 
   ! C-PML absorbing boundary conditions
   NSPEC_CPML = 0
   if( PML_CONDITIONS ) then
-    read(27) NSPEC_CPML
-    read(27) CPML_width_x
-    read(27) CPML_width_y
-    read(27) CPML_width_z
+    if(I_should_read_the_database) then
+      read(27) NSPEC_CPML
+      read(27) CPML_width_x
+      read(27) CPML_width_y
+      read(27) CPML_width_z
+    endif
 
     allocate(is_CPML(NSPEC_AB),stat=ier)
     if(ier /= 0) stop 'error allocating array is_CPML'
@@ -349,31 +356,33 @@
       allocate(alpha_store_z(NGLLX,NGLLY,NGLLZ,NSPEC_CPML),stat=ier)
       if(ier /= 0) stop 'error allocating array alpha_store'
 
-      read(27) CPML_regions
-      read(27) CPML_to_spec
-      read(27) is_CPML
-      read(27) d_store_x
-      read(27) d_store_y
-      read(27) d_store_z
-      read(27) k_store_x
-      read(27) k_store_y
-      read(27) k_store_z
-      read(27) alpha_store_x
-      read(27) alpha_store_y
-      read(27) alpha_store_z
+      if(I_should_read_the_database) then
+        read(27) CPML_regions
+        read(27) CPML_to_spec
+        read(27) is_CPML
+        read(27) d_store_x
+        read(27) d_store_y
+        read(27) d_store_z
+        read(27) k_store_x
+        read(27) k_store_y
+        read(27) k_store_z
+        read(27) alpha_store_x
+        read(27) alpha_store_y
+        read(27) alpha_store_z
+      endif
 
       if((SIMULATION_TYPE == 1 .and. SAVE_FORWARD) .or. SIMULATION_TYPE == 3) then
-        read(27) nglob_interface_PML_acoustic
-        read(27) nglob_interface_PML_elastic
+        if(I_should_read_the_database) read(27) nglob_interface_PML_acoustic
+        if(I_should_read_the_database) read(27) nglob_interface_PML_elastic
         if(nglob_interface_PML_acoustic > 0) then
           allocate(points_interface_PML_acoustic(nglob_interface_PML_acoustic),stat=ier)
           if(ier /= 0) stop 'error allocating array points_interface_PML_acoustic'
-          read(27) points_interface_PML_acoustic
+          if(I_should_read_the_database) read(27) points_interface_PML_acoustic
         endif
         if(nglob_interface_PML_elastic > 0) then
           allocate(points_interface_PML_elastic(nglob_interface_PML_elastic),stat=ier)
           if(ier /= 0) stop 'error allocating array points_interface_PML_elastic'
-          read(27) points_interface_PML_elastic
+          if(I_should_read_the_database) read(27) points_interface_PML_elastic
         endif
       endif
     endif
@@ -383,7 +392,7 @@
   endif
 
   ! absorbing boundary surface
-  read(27) num_abs_boundary_faces
+  if(I_should_read_the_database) read(27) num_abs_boundary_faces
 
   ! checks
   if( num_abs_boundary_faces < 0 ) then
@@ -400,9 +409,9 @@
   if (COUPLE_WITH_DSM) then
     allocate(Veloc_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces))
     allocate(Tract_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces))
-    open(unit=IIN_veloc_dsm,file=dsmname(1:len_trim(dsmname))//'vel.bin',status='old', &
+    if(I_should_read_the_database) open(unit=IIN_veloc_dsm,file=dsmname(1:len_trim(dsmname))//'vel.bin',status='old', &
          action='read',form='unformatted',iostat=ier)
-    open(unit=IIN_tract_dsm,file=dsmname(1:len_trim(dsmname))//'tract.bin',status='old', &
+    if(I_should_read_the_database) open(unit=IIN_tract_dsm,file=dsmname(1:len_trim(dsmname))//'tract.bin',status='old', &
          action='read',form='unformatted',iostat=ier)
   else
     allocate(Veloc_dsm_boundary(1,1,1,1))
@@ -412,93 +421,109 @@
 
   if(PML_CONDITIONS)then
     if (num_abs_boundary_faces > 0) then
-      read(27) abs_boundary_ispec
-      read(27) abs_boundary_ijk
-      read(27) abs_boundary_jacobian2Dw
-      read(27) abs_boundary_normal
+      if(I_should_read_the_database) then
+        read(27) abs_boundary_ispec
+        read(27) abs_boundary_ijk
+        read(27) abs_boundary_jacobian2Dw
+        read(27) abs_boundary_normal
+      endif
     endif
   else
     if (num_abs_boundary_faces > 0) then
-      read(27) abs_boundary_ispec
-      read(27) abs_boundary_ijk
-      read(27) abs_boundary_jacobian2Dw
-      read(27) abs_boundary_normal
+      if(I_should_read_the_database) then
+        read(27) abs_boundary_ispec
+        read(27) abs_boundary_ijk
+        read(27) abs_boundary_jacobian2Dw
+        read(27) abs_boundary_normal
+      endif
       if (STACEY_ABSORBING_CONDITIONS) then
         ! store mass matrix contributions
         if (ELASTIC_SIMULATION) then
-          read(27) rmassx
-          read(27) rmassy
-          read(27) rmassz
+          if(I_should_read_the_database) then
+            read(27) rmassx
+            read(27) rmassy
+            read(27) rmassz
+          endif
         endif
         if (ACOUSTIC_SIMULATION) then
-          read(27) rmassz_acoustic
+          if(I_should_read_the_database) read(27) rmassz_acoustic
         endif
       endif
     endif
   endif
 
-  read(27) nspec2D_xmin
-  read(27) nspec2D_xmax
-  read(27) nspec2D_ymin
-  read(27) nspec2D_ymax
-  read(27) NSPEC2D_BOTTOM
-  read(27) NSPEC2D_TOP
+  if(I_should_read_the_database) then
+    read(27) nspec2D_xmin
+    read(27) nspec2D_xmax
+    read(27) nspec2D_ymin
+    read(27) nspec2D_ymax
+    read(27) NSPEC2D_BOTTOM
+    read(27) NSPEC2D_TOP
+  endif
 
   allocate(ibelm_xmin(nspec2D_xmin),ibelm_xmax(nspec2D_xmax), &
            ibelm_ymin(nspec2D_ymin),ibelm_ymax(nspec2D_ymax), &
            ibelm_bottom(NSPEC2D_BOTTOM),ibelm_top(NSPEC2D_TOP),stat=ier)
   if(ier /= 0) stop 'error allocating arrays ibelm_xmin,ibelm_xmax etc.'
-  read(27) ibelm_xmin
-  read(27) ibelm_xmax
-  read(27) ibelm_ymin
-  read(27) ibelm_ymax
-  read(27) ibelm_bottom
-  read(27) ibelm_top
+  if(I_should_read_the_database) then
+    read(27) ibelm_xmin
+    read(27) ibelm_xmax
+    read(27) ibelm_ymin
+    read(27) ibelm_ymax
+    read(27) ibelm_bottom
+    read(27) ibelm_top
+  endif
 
   ! free surface
-  read(27) num_free_surface_faces
+  if(I_should_read_the_database) read(27) num_free_surface_faces
   allocate(free_surface_ispec(num_free_surface_faces), &
            free_surface_ijk(3,NGLLSQUARE,num_free_surface_faces), &
            free_surface_jacobian2Dw(NGLLSQUARE,num_free_surface_faces), &
            free_surface_normal(NDIM,NGLLSQUARE,num_free_surface_faces),stat=ier)
   if(ier /= 0) stop 'error allocating arrays free_surface_ispec etc.'
   if( num_free_surface_faces > 0 ) then
-    read(27) free_surface_ispec
-    read(27) free_surface_ijk
-    read(27) free_surface_jacobian2Dw
-    read(27) free_surface_normal
+    if(I_should_read_the_database) then
+      read(27) free_surface_ispec
+      read(27) free_surface_ijk
+      read(27) free_surface_jacobian2Dw
+      read(27) free_surface_normal
+    endif
   endif
 
   ! acoustic-elastic coupling surface
-  read(27) num_coupling_ac_el_faces
+  if(I_should_read_the_database) read(27) num_coupling_ac_el_faces
   allocate(coupling_ac_el_normal(NDIM,NGLLSQUARE,num_coupling_ac_el_faces), &
            coupling_ac_el_jacobian2Dw(NGLLSQUARE,num_coupling_ac_el_faces), &
            coupling_ac_el_ijk(3,NGLLSQUARE,num_coupling_ac_el_faces), &
            coupling_ac_el_ispec(num_coupling_ac_el_faces),stat=ier)
   if( ier /= 0 ) stop 'error allocating array coupling_ac_el_normal etc.'
   if( num_coupling_ac_el_faces > 0 ) then
-    read(27) coupling_ac_el_ispec
-    read(27) coupling_ac_el_ijk
-    read(27) coupling_ac_el_jacobian2Dw
-    read(27) coupling_ac_el_normal
+    if(I_should_read_the_database) then
+      read(27) coupling_ac_el_ispec
+      read(27) coupling_ac_el_ijk
+      read(27) coupling_ac_el_jacobian2Dw
+      read(27) coupling_ac_el_normal
+    endif
   endif
 
   ! acoustic-poroelastic coupling surface
-  read(27) num_coupling_ac_po_faces
+  if(I_should_read_the_database) read(27) num_coupling_ac_po_faces
   allocate(coupling_ac_po_normal(NDIM,NGLLSQUARE,num_coupling_ac_po_faces), &
            coupling_ac_po_jacobian2Dw(NGLLSQUARE,num_coupling_ac_po_faces), &
            coupling_ac_po_ijk(3,NGLLSQUARE,num_coupling_ac_po_faces), &
            coupling_ac_po_ispec(num_coupling_ac_po_faces),stat=ier)
   if( ier /= 0 ) stop 'error allocating array coupling_ac_po_normal etc.'
   if( num_coupling_ac_po_faces > 0 ) then
-    read(27) coupling_ac_po_ispec
-    read(27) coupling_ac_po_ijk
-    read(27) coupling_ac_po_jacobian2Dw
-    read(27) coupling_ac_po_normal
+    if(I_should_read_the_database) then
+      read(27) coupling_ac_po_ispec
+      read(27) coupling_ac_po_ijk
+      read(27) coupling_ac_po_jacobian2Dw
+      read(27) coupling_ac_po_normal
+    endif
   endif
 
   ! elastic-poroelastic coupling surface
-  read(27) num_coupling_el_po_faces
+  if(I_should_read_the_database) read(27) num_coupling_el_po_faces
   allocate(coupling_el_po_normal(NDIM,NGLLSQUARE,num_coupling_el_po_faces), &
            coupling_el_po_jacobian2Dw(NGLLSQUARE,num_coupling_el_po_faces), &
            coupling_el_po_ijk(3,NGLLSQUARE,num_coupling_el_po_faces), &
@@ -507,106 +532,124 @@
            coupling_po_el_ispec(num_coupling_el_po_faces),stat=ier)
   if( ier /= 0 ) stop 'error allocating array coupling_el_po_normal etc.'
   if( num_coupling_el_po_faces > 0 ) then
-    read(27) coupling_el_po_ispec
-    read(27) coupling_po_el_ispec
-    read(27) coupling_el_po_ijk
-    read(27) coupling_po_el_ijk
-    read(27) coupling_el_po_jacobian2Dw
-    read(27) coupling_el_po_normal
+    if(I_should_read_the_database) then
+      read(27) coupling_el_po_ispec
+      read(27) coupling_po_el_ispec
+      read(27) coupling_el_po_ijk
+      read(27) coupling_po_el_ijk
+      read(27) coupling_el_po_jacobian2Dw
+      read(27) coupling_el_po_normal
+    endif
   endif
 
   ! MPI interfaces
-  read(27) num_interfaces_ext_mesh
+  if(I_should_read_the_database) read(27) num_interfaces_ext_mesh
   allocate(my_neighbours_ext_mesh(num_interfaces_ext_mesh), &
            nibool_interfaces_ext_mesh(num_interfaces_ext_mesh),stat=ier)
   if( ier /= 0 ) stop 'error allocating array my_neighbours_ext_mesh etc.'
   if( num_interfaces_ext_mesh > 0 ) then
-    read(27) max_nibool_interfaces_ext_mesh
+    if(I_should_read_the_database) read(27) max_nibool_interfaces_ext_mesh
     allocate(ibool_interfaces_ext_mesh(max_nibool_interfaces_ext_mesh,num_interfaces_ext_mesh),stat=ier)
     if( ier /= 0 ) stop 'error allocating array ibool_interfaces_ext_mesh'
-    read(27) my_neighbours_ext_mesh
-    read(27) nibool_interfaces_ext_mesh
-    read(27) ibool_interfaces_ext_mesh
+    if(I_should_read_the_database) then
+      read(27) my_neighbours_ext_mesh
+      read(27) nibool_interfaces_ext_mesh
+      read(27) ibool_interfaces_ext_mesh
+    endif
   else
     max_nibool_interfaces_ext_mesh = 0
     allocate(ibool_interfaces_ext_mesh(0,0),stat=ier)
   endif
 
   if( ELASTIC_SIMULATION .and. ANISOTROPY ) then
-    read(27) c11store
-    read(27) c12store
-    read(27) c13store
-    read(27) c14store
-    read(27) c15store
-    read(27) c16store
-    read(27) c22store
-    read(27) c23store
-    read(27) c24store
-    read(27) c25store
-    read(27) c26store
-    read(27) c33store
-    read(27) c34store
-    read(27) c35store
-    read(27) c36store
-    read(27) c44store
-    read(27) c45store
-    read(27) c46store
-    read(27) c55store
-    read(27) c56store
-    read(27) c66store
+    if(I_should_read_the_database) then
+      read(27) c11store
+      read(27) c12store
+      read(27) c13store
+      read(27) c14store
+      read(27) c15store
+      read(27) c16store
+      read(27) c22store
+      read(27) c23store
+      read(27) c24store
+      read(27) c25store
+      read(27) c26store
+      read(27) c33store
+      read(27) c34store
+      read(27) c35store
+      read(27) c36store
+      read(27) c44store
+      read(27) c45store
+      read(27) c46store
+      read(27) c55store
+      read(27) c56store
+      read(27) c66store
+    endif
   endif
 
   ! inner / outer elements
   allocate(ispec_is_inner(NSPEC_AB),stat=ier)
   if( ier /= 0 ) stop 'error allocating array ispec_is_inner'
-  read(27) ispec_is_inner
+  if(I_should_read_the_database) read(27) ispec_is_inner
 
   if( ACOUSTIC_SIMULATION ) then
-    read(27) nspec_inner_acoustic,nspec_outer_acoustic
-    read(27) num_phase_ispec_acoustic
+    if(I_should_read_the_database) then
+      read(27) nspec_inner_acoustic,nspec_outer_acoustic
+      read(27) num_phase_ispec_acoustic
+    endif
     if( num_phase_ispec_acoustic < 0 ) stop 'error acoustic simulation: num_phase_ispec_acoustic is < zero'
     allocate( phase_ispec_inner_acoustic(num_phase_ispec_acoustic,2),stat=ier)
     if( ier /= 0 ) stop 'error allocating array phase_ispec_inner_acoustic'
-    if(num_phase_ispec_acoustic > 0 ) read(27) phase_ispec_inner_acoustic
+    if(num_phase_ispec_acoustic > 0 ) then
+      if(I_should_read_the_database) read(27) phase_ispec_inner_acoustic
+    endif
   endif
 
   if( ELASTIC_SIMULATION ) then
-    read(27) nspec_inner_elastic,nspec_outer_elastic
-    read(27) num_phase_ispec_elastic
+    if(I_should_read_the_database) then
+      read(27) nspec_inner_elastic,nspec_outer_elastic
+      read(27) num_phase_ispec_elastic
+    endif
     if( num_phase_ispec_elastic < 0 ) stop 'error elastic simulation: num_phase_ispec_elastic is < zero'
     allocate( phase_ispec_inner_elastic(num_phase_ispec_elastic,2),stat=ier)
     if( ier /= 0 ) stop 'error allocating array phase_ispec_inner_elastic'
-    if(num_phase_ispec_elastic > 0 ) read(27) phase_ispec_inner_elastic
+    if(num_phase_ispec_elastic > 0 ) then
+      if(I_should_read_the_database) read(27) phase_ispec_inner_elastic
+    endif
   endif
 
   if( POROELASTIC_SIMULATION ) then
-    read(27) nspec_inner_poroelastic,nspec_outer_poroelastic
-    read(27) num_phase_ispec_poroelastic
+    if(I_should_read_the_database) then
+      read(27) nspec_inner_poroelastic,nspec_outer_poroelastic
+      read(27) num_phase_ispec_poroelastic
+    endif
     if( num_phase_ispec_poroelastic < 0 ) stop 'error poroelastic simulation: num_phase_ispec_poroelastic is < zero'
     allocate( phase_ispec_inner_poroelastic(num_phase_ispec_poroelastic,2),stat=ier)
     if( ier /= 0 ) stop 'error allocating array phase_ispec_inner_poroelastic'
-    if(num_phase_ispec_poroelastic > 0 ) read(27) phase_ispec_inner_poroelastic
+    if(num_phase_ispec_poroelastic > 0 ) then
+      if(I_should_read_the_database) read(27) phase_ispec_inner_poroelastic
+    endif
   endif
 
 ! mesh coloring for GPUs
   if( USE_MESH_COLORING_GPU ) then
     ! acoustic domain colors
     if( ACOUSTIC_SIMULATION ) then
-      read(27) num_colors_outer_acoustic,num_colors_inner_acoustic
+      if(I_should_read_the_database) read(27) num_colors_outer_acoustic,num_colors_inner_acoustic
 
       allocate(num_elem_colors_acoustic(num_colors_outer_acoustic + num_colors_inner_acoustic),stat=ier)
       if( ier /= 0 ) stop 'error allocating num_elem_colors_acoustic array'
 
-      read(27) num_elem_colors_acoustic
+      if(I_should_read_the_database) read(27) num_elem_colors_acoustic
     endif
     ! elastic domain colors
     if( ELASTIC_SIMULATION ) then
-      read(27) num_colors_outer_elastic,num_colors_inner_elastic
+      if(I_should_read_the_database) read(27) num_colors_outer_elastic,num_colors_inner_elastic
 
       allocate(num_elem_colors_elastic(num_colors_outer_elastic + num_colors_inner_elastic),stat=ier)
       if( ier /= 0 ) stop 'error allocating num_elem_colors_elastic array'
 
-      read(27) num_elem_colors_elastic
+      if(I_should_read_the_database) read(27) num_elem_colors_elastic
     endif
   else
     ! allocates dummy arrays
@@ -623,7 +666,7 @@
       if( ier /= 0 ) stop 'error allocating num_elem_colors_elastic array'
     endif
   endif
-  close(27)
+  if(I_should_read_the_database) close(27)
 
   ! outputs total element numbers
   call sum_all_i(count(ispec_is_acoustic(:)),inum)
@@ -904,16 +947,17 @@
         call read_moho_mesh_adjoint_adios()
       else
         ! boundary elements
-        !open(unit=27,file=prname(1:len_trim(prname))//'ibelm_moho.bin',status='unknown',form='unformatted')
-        open(unit=27,file=prname(1:len_trim(prname))//'ibelm_moho.bin',status='old',&
-             form='unformatted',iostat=ier)
-        if( ier /= 0 ) then
-          print*,'error: could not open ibelm_moho '
-          print*,'path: ',prname(1:len_trim(prname))//'ibelm_moho.bin'
-          call exit_mpi(myrank,'error opening ibelm_moho')
+        if(I_should_read_the_database) then
+          open(unit=27,file=prname(1:len_trim(prname))//'ibelm_moho.bin',status='old',&
+             action='read',form='unformatted',iostat=ier)
+          if( ier /= 0 ) then
+            print*,'error: could not open ibelm_moho '
+            print*,'path: ',prname(1:len_trim(prname))//'ibelm_moho.bin'
+            call exit_mpi(myrank,'error opening ibelm_moho')
+          endif
         endif
 
-        read(27) NSPEC2D_MOHO
+        if(I_should_read_the_database) read(27) NSPEC2D_MOHO
 
         ! allocates arrays for moho mesh
         allocate(ibelm_moho_bot(NSPEC2D_MOHO), &
@@ -924,39 +968,49 @@
                  ijk_moho_top(3,NGLLSQUARE,NSPEC2D_MOHO),stat=ier)
         if( ier /= 0 ) stop 'error allocating array ibelm_moho_bot etc.'
 
-        read(27) ibelm_moho_top
-        read(27) ibelm_moho_bot
-        read(27) ijk_moho_top
-        read(27) ijk_moho_bot
+        if(I_should_read_the_database) then
+          read(27) ibelm_moho_top
+          read(27) ibelm_moho_bot
+          read(27) ijk_moho_top
+          read(27) ijk_moho_bot
+        endif
 
-        close(27)
+        if(I_should_read_the_database) close(27)
 
         ! normals
-        open(unit=27,file=prname(1:len_trim(prname))//'normal_moho.bin',status='old',&
-             form='unformatted',iostat=ier)
-        if( ier /= 0 ) then
-          print*,'error: could not open normal_moho '
-          print*,'path: ',prname(1:len_trim(prname))//'normal_moho.bin'
-          call exit_mpi(myrank,'error opening normal_moho')
+        if(I_should_read_the_database) then
+          open(unit=27,file=prname(1:len_trim(prname))//'normal_moho.bin',status='old',&
+             action='read',form='unformatted',iostat=ier)
+          if( ier /= 0 ) then
+            print*,'error: could not open normal_moho '
+            print*,'path: ',prname(1:len_trim(prname))//'normal_moho.bin'
+            call exit_mpi(myrank,'error opening normal_moho')
+          endif
         endif
 
-        read(27) normal_moho_top
-        read(27) normal_moho_bot
-        close(27)
+        if(I_should_read_the_database) then
+          read(27) normal_moho_top
+          read(27) normal_moho_bot
+        endif
+        if(I_should_read_the_database) close(27)
 
         ! flags
-        open(unit=27,file=prname(1:len_trim(prname))//'is_moho.bin',status='old',&
-             form='unformatted',iostat=ier)
-        if( ier /= 0 ) then
-          print*,'error: could not open is_moho '
-          print*,'path: ',prname(1:len_trim(prname))//'is_moho.bin'
-          call exit_mpi(myrank,'error opening is_moho')
+        if(I_should_read_the_database) then
+          open(unit=27,file=prname(1:len_trim(prname))//'is_moho.bin',status='old',&
+             action='read',form='unformatted',iostat=ier)
+          if( ier /= 0 ) then
+            print*,'error: could not open is_moho '
+            print*,'path: ',prname(1:len_trim(prname))//'is_moho.bin'
+            call exit_mpi(myrank,'error opening is_moho')
+          endif
         endif
 
-        read(27) is_moho_top
-        read(27) is_moho_bot
+        if(I_should_read_the_database) then
+          read(27) is_moho_top
+          read(27) is_moho_bot
+        endif
 
-        close(27)
+        if(I_should_read_the_database) close(27)
       endif
 
       ! moho kernel
@@ -1135,16 +1189,26 @@
 
   integer :: ier
 
-  open(unit=IIN,file=prname(1:len_trim(prname))//'external_mesh.bin',status='old',&
+  if(I_should_read_the_database) then
+    open(unit=IIN,file=prname(1:len_trim(prname))//'external_mesh.bin',status='old',&
        action='read',form='unformatted',iostat=ier)
-  if( ier /= 0 ) then
-    print*,'error: could not open database '
-    print*,'path: ',prname(1:len_trim(prname))//'external_mesh.bin'
-    call exit_mpi(myrank,'error opening database')
+    if( ier /= 0 ) then
+      print*,'error: could not open database '
+      print*,'path: ',prname(1:len_trim(prname))//'external_mesh.bin'
+      call exit_mpi(myrank,'error opening database')
+    endif
   endif
-  read(IIN) NSPEC_AB
-  read(IIN) NGLOB_AB
-  close(IIN)
+
+  if(I_should_read_the_database) then
+    read(IIN) NSPEC_AB
+    read(IIN) NGLOB_AB
+  endif
+  if(I_should_broadcast_the_database) then
+    call bcast_all_one_i_for_database(NSPEC_AB)
+    call bcast_all_one_i_for_database(NGLOB_AB)
+  endif
+
+  if(I_should_read_the_database) close(IIN)
 
   end subroutine read_mesh_for_init
 



More information about the CIG-COMMITS mailing list