[cig-commits] r12863 - in seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta: DATA/crust2.0 setup src
dkomati1 at geodynamics.org
dkomati1 at geodynamics.org
Thu Sep 11 08:03:35 PDT 2008
Author: dkomati1
Date: 2008-09-11 08:03:34 -0700 (Thu, 11 Sep 2008)
New Revision: 12863
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/crust2.0/smooth_crust2.0_once_and_for_all.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/assemble_MPI_vector.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/attenuation_model.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_chunk_buffers.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/crustal_model.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/extract_all_seismos_from_large_file.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_event_info.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_receivers.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_sources.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/main_program.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.F90
Log:
fixed the MPI_BCAST problem, drastically reduced the total number of calls to MPI_BCAST in the code, and permanently smoothed the crust in the data file instead of smoothing it in the code, which is costly at very high resolution
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/crust2.0/smooth_crust2.0_once_and_for_all.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/crust2.0/smooth_crust2.0_once_and_for_all.f90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/crust2.0/smooth_crust2.0_once_and_for_all.f90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -37,6 +37,11 @@
include "../../setup/constants.h"
+! crustal model parameters for crust2.0
+ integer, parameter :: NKEYS_CRUST = 359
+ integer, parameter :: NLAYERS_CRUST = 8
+ integer, parameter :: NCAP_CRUST = 180
+
! crustal_model_variables
type crustal_model_variables
sequence
@@ -117,6 +122,11 @@
implicit none
include "../../setup/constants.h"
+! crustal model parameters for crust2.0
+ integer, parameter :: NKEYS_CRUST = 359
+ integer, parameter :: NLAYERS_CRUST = 8
+ integer, parameter :: NCAP_CRUST = 180
+
! crustal_model_variables
type crustal_model_variables
sequence
@@ -175,6 +185,11 @@
implicit none
include "../../setup/constants.h"
+! crustal model parameters for crust2.0
+ integer, parameter :: NKEYS_CRUST = 359
+ integer, parameter :: NLAYERS_CRUST = 8
+ integer, parameter :: NCAP_CRUST = 180
+
integer, parameter :: NTHETA = 2
integer, parameter :: NPHI = 10
double precision, parameter :: CAP = 2.0d0*PI/180.0d0
@@ -337,6 +352,11 @@
implicit none
include "../../setup/constants.h"
+! crustal model parameters for crust2.0
+ integer, parameter :: NKEYS_CRUST = 359
+ integer, parameter :: NLAYERS_CRUST = 8
+ integer, parameter :: NCAP_CRUST = 180
+
! argument variables
integer ierr
double precision rhtyp(NLAYERS_CRUST),thtp(NLAYERS_CRUST)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h 2008-09-11 15:03:34 UTC (rev 12863)
@@ -35,16 +35,15 @@
!! DK DK to false in all cases etc
logical, parameter :: PATCH_FOR_GORDON_BELL = .true.
+! save partial seismograms every 10,000 time steps or not
+ logical, parameter :: SAVE_PARTIAL_SEISMOGRAMS = .false.
+
! (much) faster detection of receivers at high resolution: use grid points only
logical, parameter :: FASTER_RECEIVERS_POINTS_ONLY = .true.
! suppress calculation and storage of seismograms if needed
logical, parameter :: COMPUTE_STORE_SEISMOGRAMS = .true.
-! decrease the number of MPI messages by 3 but increase the size
-! of several MPI buffers by 3 in order to do that
- logical, parameter :: FEWER_MESSAGES_LARGER_BUFFERS = .true.
-
!! DK DK for Gordon Bell
! integer, parameter :: SEA99_VS_DIM1 = 100, SEA99_VS_DIM2 = 100, SEA99_VS_DIM3 = 100
integer, parameter :: SEA99_VS_DIM1 = 1, SEA99_VS_DIM2 = 1, SEA99_VS_DIM3 = 1
@@ -79,6 +78,11 @@
! uncomment this to write messages to the screen (slows down the code)
! integer, parameter :: IMAIN = ISTANDARD_OUTPUT
+! number of values read in Par_file that we need to broadcast
+ integer, parameter :: NVALUES_bcast_integer = 41
+ integer, parameter :: NVALUES_bcast_double_precision = 30
+ integer, parameter :: NVALUES_bcast_logical = 33
+
! R_EARTH is the radius of the bottom of the oceans (radius of Earth in m)
double precision, parameter :: R_EARTH = 6371000.d0
! uncomment line below for PREM with oceans
@@ -453,11 +457,6 @@
double precision,parameter :: LON_MIN = 130.d0
double precision,parameter :: DEP_MAX = 500.d0
-! crustal model parameters for crust2.0
- integer, parameter :: NKEYS_CRUST = 359
- integer, parameter :: NLAYERS_CRUST = 8
- integer, parameter :: NCAP_CRUST = 180
-
! use sedimentary layers of crust 2.0
logical, parameter :: INCLUDE_SEDIMENTS_CRUST = .true.
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/assemble_MPI_vector.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/assemble_MPI_vector.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/assemble_MPI_vector.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -45,7 +45,7 @@
buffer_send_chunkcorners_vector,buffer_recv_chunkcorners_vector, &
NUMMSGS_FACES,NUM_MSG_TYPES,NCORNERSCHUNKS, &
NPROC_XI,NPROC_ETA,NGLOB1D_RADIAL_crust_mantle, &
- NGLOB1D_RADIAL_inner_core,NCHUNKS,NDIM_smaller_buffers)
+ NGLOB1D_RADIAL_inner_core,NCHUNKS)
implicit none
@@ -75,7 +75,7 @@
integer npoin2D_faces_inner_core(NUMFACES_SHARED)
integer NGLOB1D_RADIAL_crust_mantle,NGLOB1D_RADIAL_inner_core,NPROC_XI,NPROC_ETA
- integer NUMMSGS_FACES,NUM_MSG_TYPES,NCORNERSCHUNKS,NDIM_smaller_buffers
+ integer NUMMSGS_FACES,NUM_MSG_TYPES,NCORNERSCHUNKS
! for addressing of the slices
integer, dimension(NCHUNKS,0:NPROC_XI-1,0:NPROC_ETA-1) :: addressing
@@ -94,7 +94,7 @@
integer :: npoin2D_max_all
integer, dimension(NGLOB2DMAX_XY_VAL_CM,NUMFACES_SHARED) :: iboolfaces_crust_mantle
integer, dimension(NGLOB2DMAX_XY_VAL_IC,NUMFACES_SHARED) :: iboolfaces_inner_core
- real(kind=CUSTOM_REAL), dimension(NDIM_smaller_buffers,npoin2D_max_all) :: buffer_send_faces_vector,buffer_received_faces_vector
+ real(kind=CUSTOM_REAL), dimension(NDIM,npoin2D_max_all) :: buffer_send_faces_vector,buffer_received_faces_vector
! buffers for send and receive between corners of the chunks
! size of buffers is the sum of two sizes because we handle two regions in the same MPI call
@@ -116,7 +116,7 @@
integer :: ipoin,ipoin2D,ipoin1D
integer :: sender,receiver
- integer :: imsg,imsg_loop,iloop
+ integer :: imsg,imsg_loop
integer :: icount_faces,npoin2D_chunks_all
integer :: npoin2D_xi_all,npoin2D_eta_all,NGLOB1D_RADIAL_all,ioffset
@@ -140,9 +140,6 @@
!---- assemble the contributions between slices using MPI
!----
-! loop three times if using smaller buffers, and only once if using larger buffers
- do iloop = 1,NDIM + 1 - NDIM_smaller_buffers
-
!----
!---- first assemble along xi using the 2-D topology
!----
@@ -155,19 +152,15 @@
! slices copy the right face into the buffer
do ipoin = 1,npoin2D_xi_crust_mantle
- buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(iloop,iboolright_xi_crust_mantle(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolright_xi_crust_mantle(ipoin))
- buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolright_xi_crust_mantle(ipoin))
- endif
+ buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(1,iboolright_xi_crust_mantle(ipoin))
+ buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolright_xi_crust_mantle(ipoin))
+ buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolright_xi_crust_mantle(ipoin))
enddo
do ipoin = 1,npoin2D_xi_inner_core
- buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(iloop,iboolright_xi_inner_core(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolright_xi_inner_core(ipoin))
- buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolright_xi_inner_core(ipoin))
- endif
+ buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(1,iboolright_xi_inner_core(ipoin))
+ buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolright_xi_inner_core(ipoin))
+ buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolright_xi_inner_core(ipoin))
enddo
! send messages forward along each row
@@ -186,8 +179,8 @@
receiver = addressing(ichunk,iproc_xi + 1,iproc_eta)
endif
#ifdef USE_MPI
- call MPI_SENDRECV(buffer_send_faces_vector,NDIM_smaller_buffers*npoin2D_xi_all,CUSTOM_MPI_TYPE,receiver, &
- itag2,buffer_received_faces_vector,NDIM_smaller_buffers*npoin2D_xi_all,CUSTOM_MPI_TYPE,sender, &
+ call MPI_SENDRECV(buffer_send_faces_vector,NDIM*npoin2D_xi_all,CUSTOM_MPI_TYPE,receiver, &
+ itag2,buffer_received_faces_vector,NDIM*npoin2D_xi_all,CUSTOM_MPI_TYPE,sender, &
itag,MPI_COMM_WORLD,msg_status,ier)
#endif
@@ -195,25 +188,21 @@
if(iproc_xi > 0) then
do ipoin = 1,npoin2D_xi_crust_mantle
- accel_crust_mantle(iloop,iboolleft_xi_crust_mantle(ipoin)) = accel_crust_mantle(iloop,iboolleft_xi_crust_mantle(ipoin)) + &
+ accel_crust_mantle(1,iboolleft_xi_crust_mantle(ipoin)) = accel_crust_mantle(1,iboolleft_xi_crust_mantle(ipoin)) + &
buffer_received_faces_vector(1,ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_crust_mantle(2,iboolleft_xi_crust_mantle(ipoin)) = accel_crust_mantle(2,iboolleft_xi_crust_mantle(ipoin)) + &
+ accel_crust_mantle(2,iboolleft_xi_crust_mantle(ipoin)) = accel_crust_mantle(2,iboolleft_xi_crust_mantle(ipoin)) + &
buffer_received_faces_vector(2,ipoin)
- accel_crust_mantle(3,iboolleft_xi_crust_mantle(ipoin)) = accel_crust_mantle(3,iboolleft_xi_crust_mantle(ipoin)) + &
+ accel_crust_mantle(3,iboolleft_xi_crust_mantle(ipoin)) = accel_crust_mantle(3,iboolleft_xi_crust_mantle(ipoin)) + &
buffer_received_faces_vector(3,ipoin)
- endif
enddo
do ipoin = 1,npoin2D_xi_inner_core
- accel_inner_core(iloop,iboolleft_xi_inner_core(ipoin)) = accel_inner_core(iloop,iboolleft_xi_inner_core(ipoin)) + &
+ accel_inner_core(1,iboolleft_xi_inner_core(ipoin)) = accel_inner_core(1,iboolleft_xi_inner_core(ipoin)) + &
buffer_received_faces_vector(1,ioffset + ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_inner_core(2,iboolleft_xi_inner_core(ipoin)) = accel_inner_core(2,iboolleft_xi_inner_core(ipoin)) + &
+ accel_inner_core(2,iboolleft_xi_inner_core(ipoin)) = accel_inner_core(2,iboolleft_xi_inner_core(ipoin)) + &
buffer_received_faces_vector(2,ioffset + ipoin)
- accel_inner_core(3,iboolleft_xi_inner_core(ipoin)) = accel_inner_core(3,iboolleft_xi_inner_core(ipoin)) + &
+ accel_inner_core(3,iboolleft_xi_inner_core(ipoin)) = accel_inner_core(3,iboolleft_xi_inner_core(ipoin)) + &
buffer_received_faces_vector(3,ioffset + ipoin)
- endif
enddo
endif
@@ -222,19 +211,15 @@
! now we have to send the result back to the sender
! all slices copy the left face into the buffer
do ipoin = 1,npoin2D_xi_crust_mantle
- buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(iloop,iboolleft_xi_crust_mantle(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolleft_xi_crust_mantle(ipoin))
- buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolleft_xi_crust_mantle(ipoin))
- endif
+ buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(1,iboolleft_xi_crust_mantle(ipoin))
+ buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolleft_xi_crust_mantle(ipoin))
+ buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolleft_xi_crust_mantle(ipoin))
enddo
do ipoin = 1,npoin2D_xi_inner_core
- buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(iloop,iboolleft_xi_inner_core(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolleft_xi_inner_core(ipoin))
- buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolleft_xi_inner_core(ipoin))
- endif
+ buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(1,iboolleft_xi_inner_core(ipoin))
+ buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolleft_xi_inner_core(ipoin))
+ buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolleft_xi_inner_core(ipoin))
enddo
! send messages backward along each row
@@ -253,8 +238,8 @@
receiver = addressing(ichunk,iproc_xi - 1,iproc_eta)
endif
#ifdef USE_MPI
- call MPI_SENDRECV(buffer_send_faces_vector,NDIM_smaller_buffers*npoin2D_xi_all,CUSTOM_MPI_TYPE,receiver, &
- itag2,buffer_received_faces_vector,NDIM_smaller_buffers*npoin2D_xi_all,CUSTOM_MPI_TYPE,sender, &
+ call MPI_SENDRECV(buffer_send_faces_vector,NDIM*npoin2D_xi_all,CUSTOM_MPI_TYPE,receiver, &
+ itag2,buffer_received_faces_vector,NDIM*npoin2D_xi_all,CUSTOM_MPI_TYPE,sender, &
itag,MPI_COMM_WORLD,msg_status,ier)
#endif
@@ -262,19 +247,15 @@
if(iproc_xi < NPROC_XI-1) then
do ipoin = 1,npoin2D_xi_crust_mantle
- accel_crust_mantle(iloop,iboolright_xi_crust_mantle(ipoin)) = buffer_received_faces_vector(1,ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_crust_mantle(2,iboolright_xi_crust_mantle(ipoin)) = buffer_received_faces_vector(2,ipoin)
- accel_crust_mantle(3,iboolright_xi_crust_mantle(ipoin)) = buffer_received_faces_vector(3,ipoin)
- endif
+ accel_crust_mantle(1,iboolright_xi_crust_mantle(ipoin)) = buffer_received_faces_vector(1,ipoin)
+ accel_crust_mantle(2,iboolright_xi_crust_mantle(ipoin)) = buffer_received_faces_vector(2,ipoin)
+ accel_crust_mantle(3,iboolright_xi_crust_mantle(ipoin)) = buffer_received_faces_vector(3,ipoin)
enddo
do ipoin = 1,npoin2D_xi_inner_core
- accel_inner_core(iloop,iboolright_xi_inner_core(ipoin)) = buffer_received_faces_vector(1,ioffset + ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_inner_core(2,iboolright_xi_inner_core(ipoin)) = buffer_received_faces_vector(2,ioffset + ipoin)
- accel_inner_core(3,iboolright_xi_inner_core(ipoin)) = buffer_received_faces_vector(3,ioffset + ipoin)
- endif
+ accel_inner_core(1,iboolright_xi_inner_core(ipoin)) = buffer_received_faces_vector(1,ioffset + ipoin)
+ accel_inner_core(2,iboolright_xi_inner_core(ipoin)) = buffer_received_faces_vector(2,ioffset + ipoin)
+ accel_inner_core(3,iboolright_xi_inner_core(ipoin)) = buffer_received_faces_vector(3,ioffset + ipoin)
enddo
endif
@@ -293,19 +274,15 @@
! slices copy the right face into the buffer
do ipoin = 1,npoin2D_eta_crust_mantle
- buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(iloop,iboolright_eta_crust_mantle(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolright_eta_crust_mantle(ipoin))
- buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolright_eta_crust_mantle(ipoin))
- endif
+ buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(1,iboolright_eta_crust_mantle(ipoin))
+ buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolright_eta_crust_mantle(ipoin))
+ buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolright_eta_crust_mantle(ipoin))
enddo
do ipoin = 1,npoin2D_eta_inner_core
- buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(iloop,iboolright_eta_inner_core(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolright_eta_inner_core(ipoin))
- buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolright_eta_inner_core(ipoin))
- endif
+ buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(1,iboolright_eta_inner_core(ipoin))
+ buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolright_eta_inner_core(ipoin))
+ buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolright_eta_inner_core(ipoin))
enddo
! send messages forward along each row
@@ -324,8 +301,8 @@
receiver = addressing(ichunk,iproc_xi,iproc_eta + 1)
endif
#ifdef USE_MPI
- call MPI_SENDRECV(buffer_send_faces_vector,NDIM_smaller_buffers*npoin2D_eta_all,CUSTOM_MPI_TYPE,receiver, &
- itag2,buffer_received_faces_vector,NDIM_smaller_buffers*npoin2D_eta_all,CUSTOM_MPI_TYPE,sender, &
+ call MPI_SENDRECV(buffer_send_faces_vector,NDIM*npoin2D_eta_all,CUSTOM_MPI_TYPE,receiver, &
+ itag2,buffer_received_faces_vector,NDIM*npoin2D_eta_all,CUSTOM_MPI_TYPE,sender, &
itag,MPI_COMM_WORLD,msg_status,ier)
#endif
@@ -333,25 +310,21 @@
if(iproc_eta > 0) then
do ipoin = 1,npoin2D_eta_crust_mantle
- accel_crust_mantle(iloop,iboolleft_eta_crust_mantle(ipoin)) = accel_crust_mantle(iloop,iboolleft_eta_crust_mantle(ipoin)) + &
+ accel_crust_mantle(1,iboolleft_eta_crust_mantle(ipoin)) = accel_crust_mantle(1,iboolleft_eta_crust_mantle(ipoin)) + &
buffer_received_faces_vector(1,ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_crust_mantle(2,iboolleft_eta_crust_mantle(ipoin)) = accel_crust_mantle(2,iboolleft_eta_crust_mantle(ipoin)) + &
+ accel_crust_mantle(2,iboolleft_eta_crust_mantle(ipoin)) = accel_crust_mantle(2,iboolleft_eta_crust_mantle(ipoin)) + &
buffer_received_faces_vector(2,ipoin)
- accel_crust_mantle(3,iboolleft_eta_crust_mantle(ipoin)) = accel_crust_mantle(3,iboolleft_eta_crust_mantle(ipoin)) + &
+ accel_crust_mantle(3,iboolleft_eta_crust_mantle(ipoin)) = accel_crust_mantle(3,iboolleft_eta_crust_mantle(ipoin)) + &
buffer_received_faces_vector(3,ipoin)
- endif
enddo
do ipoin = 1,npoin2D_eta_inner_core
- accel_inner_core(iloop,iboolleft_eta_inner_core(ipoin)) = accel_inner_core(iloop,iboolleft_eta_inner_core(ipoin)) + &
+ accel_inner_core(1,iboolleft_eta_inner_core(ipoin)) = accel_inner_core(1,iboolleft_eta_inner_core(ipoin)) + &
buffer_received_faces_vector(1,ioffset + ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_inner_core(2,iboolleft_eta_inner_core(ipoin)) = accel_inner_core(2,iboolleft_eta_inner_core(ipoin)) + &
+ accel_inner_core(2,iboolleft_eta_inner_core(ipoin)) = accel_inner_core(2,iboolleft_eta_inner_core(ipoin)) + &
buffer_received_faces_vector(2,ioffset + ipoin)
- accel_inner_core(3,iboolleft_eta_inner_core(ipoin)) = accel_inner_core(3,iboolleft_eta_inner_core(ipoin)) + &
+ accel_inner_core(3,iboolleft_eta_inner_core(ipoin)) = accel_inner_core(3,iboolleft_eta_inner_core(ipoin)) + &
buffer_received_faces_vector(3,ioffset + ipoin)
- endif
enddo
endif
@@ -360,19 +333,15 @@
! now we have to send the result back to the sender
! all slices copy the left face into the buffer
do ipoin = 1,npoin2D_eta_crust_mantle
- buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(iloop,iboolleft_eta_crust_mantle(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolleft_eta_crust_mantle(ipoin))
- buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolleft_eta_crust_mantle(ipoin))
- endif
+ buffer_send_faces_vector(1,ipoin) = accel_crust_mantle(1,iboolleft_eta_crust_mantle(ipoin))
+ buffer_send_faces_vector(2,ipoin) = accel_crust_mantle(2,iboolleft_eta_crust_mantle(ipoin))
+ buffer_send_faces_vector(3,ipoin) = accel_crust_mantle(3,iboolleft_eta_crust_mantle(ipoin))
enddo
do ipoin = 1,npoin2D_eta_inner_core
- buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(iloop,iboolleft_eta_inner_core(ipoin))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolleft_eta_inner_core(ipoin))
- buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolleft_eta_inner_core(ipoin))
- endif
+ buffer_send_faces_vector(1,ioffset + ipoin) = accel_inner_core(1,iboolleft_eta_inner_core(ipoin))
+ buffer_send_faces_vector(2,ioffset + ipoin) = accel_inner_core(2,iboolleft_eta_inner_core(ipoin))
+ buffer_send_faces_vector(3,ioffset + ipoin) = accel_inner_core(3,iboolleft_eta_inner_core(ipoin))
enddo
! send messages backward along each row
@@ -391,8 +360,8 @@
receiver = addressing(ichunk,iproc_xi,iproc_eta - 1)
endif
#ifdef USE_MPI
- call MPI_SENDRECV(buffer_send_faces_vector,NDIM_smaller_buffers*npoin2D_eta_all,CUSTOM_MPI_TYPE,receiver, &
- itag2,buffer_received_faces_vector,NDIM_smaller_buffers*npoin2D_eta_all,CUSTOM_MPI_TYPE,sender, &
+ call MPI_SENDRECV(buffer_send_faces_vector,NDIM*npoin2D_eta_all,CUSTOM_MPI_TYPE,receiver, &
+ itag2,buffer_received_faces_vector,NDIM*npoin2D_eta_all,CUSTOM_MPI_TYPE,sender, &
itag,MPI_COMM_WORLD,msg_status,ier)
#endif
@@ -400,19 +369,15 @@
if(iproc_eta < NPROC_ETA-1) then
do ipoin = 1,npoin2D_eta_crust_mantle
- accel_crust_mantle(iloop,iboolright_eta_crust_mantle(ipoin)) = buffer_received_faces_vector(1,ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_crust_mantle(2,iboolright_eta_crust_mantle(ipoin)) = buffer_received_faces_vector(2,ipoin)
- accel_crust_mantle(3,iboolright_eta_crust_mantle(ipoin)) = buffer_received_faces_vector(3,ipoin)
- endif
+ accel_crust_mantle(1,iboolright_eta_crust_mantle(ipoin)) = buffer_received_faces_vector(1,ipoin)
+ accel_crust_mantle(2,iboolright_eta_crust_mantle(ipoin)) = buffer_received_faces_vector(2,ipoin)
+ accel_crust_mantle(3,iboolright_eta_crust_mantle(ipoin)) = buffer_received_faces_vector(3,ipoin)
enddo
do ipoin = 1,npoin2D_eta_inner_core
- accel_inner_core(iloop,iboolright_eta_inner_core(ipoin)) = buffer_received_faces_vector(1,ioffset + ipoin)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_inner_core(2,iboolright_eta_inner_core(ipoin)) = buffer_received_faces_vector(2,ioffset + ipoin)
- accel_inner_core(3,iboolright_eta_inner_core(ipoin)) = buffer_received_faces_vector(3,ioffset + ipoin)
- endif
+ accel_inner_core(1,iboolright_eta_inner_core(ipoin)) = buffer_received_faces_vector(1,ioffset + ipoin)
+ accel_inner_core(2,iboolright_eta_inner_core(ipoin)) = buffer_received_faces_vector(2,ioffset + ipoin)
+ accel_inner_core(3,iboolright_eta_inner_core(ipoin)) = buffer_received_faces_vector(3,ioffset + ipoin)
enddo
endif
@@ -452,30 +417,26 @@
ioffset = npoin2D_faces_crust_mantle(icount_faces)
#ifdef USE_MPI
- call MPI_RECV(buffer_received_faces_vector,NDIM_smaller_buffers*npoin2D_chunks_all,CUSTOM_MPI_TYPE,sender, &
+ call MPI_RECV(buffer_received_faces_vector,NDIM*npoin2D_chunks_all,CUSTOM_MPI_TYPE,sender, &
itag,MPI_COMM_WORLD,msg_status,ier)
#endif
do ipoin2D = 1,npoin2D_faces_crust_mantle(icount_faces)
- accel_crust_mantle(iloop,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = &
- accel_crust_mantle(iloop,iboolfaces_crust_mantle(ipoin2D,icount_faces)) + buffer_received_faces_vector(1,ipoin2D)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = &
- accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces)) + buffer_received_faces_vector(2,ipoin2D)
- accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = &
- accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces)) + buffer_received_faces_vector(3,ipoin2D)
- endif
+ accel_crust_mantle(1,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = &
+ accel_crust_mantle(1,iboolfaces_crust_mantle(ipoin2D,icount_faces)) + buffer_received_faces_vector(1,ipoin2D)
+ accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = &
+ accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces)) + buffer_received_faces_vector(2,ipoin2D)
+ accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = &
+ accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces)) + buffer_received_faces_vector(3,ipoin2D)
enddo
do ipoin2D = 1,npoin2D_faces_inner_core(icount_faces)
- accel_inner_core(iloop,iboolfaces_inner_core(ipoin2D,icount_faces)) = &
- accel_inner_core(iloop,iboolfaces_inner_core(ipoin2D,icount_faces)) + buffer_received_faces_vector(1,ioffset + ipoin2D)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces)) = &
- accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces)) + buffer_received_faces_vector(2,ioffset + ipoin2D)
- accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces)) = &
- accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces)) + buffer_received_faces_vector(3,ioffset + ipoin2D)
- endif
+ accel_inner_core(1,iboolfaces_inner_core(ipoin2D,icount_faces)) = &
+ accel_inner_core(1,iboolfaces_inner_core(ipoin2D,icount_faces)) + buffer_received_faces_vector(1,ioffset + ipoin2D)
+ accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces)) = &
+ accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces)) + buffer_received_faces_vector(2,ioffset + ipoin2D)
+ accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces)) = &
+ accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces)) + buffer_received_faces_vector(3,ioffset + ipoin2D)
enddo
endif
@@ -497,23 +458,19 @@
ioffset = npoin2D_faces_crust_mantle(icount_faces)
do ipoin2D = 1,npoin2D_faces_crust_mantle(icount_faces)
- buffer_send_faces_vector(1,ipoin2D) = accel_crust_mantle(iloop,iboolfaces_crust_mantle(ipoin2D,icount_faces))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ipoin2D) = accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces))
- buffer_send_faces_vector(3,ipoin2D) = accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces))
- endif
+ buffer_send_faces_vector(1,ipoin2D) = accel_crust_mantle(1,iboolfaces_crust_mantle(ipoin2D,icount_faces))
+ buffer_send_faces_vector(2,ipoin2D) = accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces))
+ buffer_send_faces_vector(3,ipoin2D) = accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces))
enddo
do ipoin2D = 1,npoin2D_faces_inner_core(icount_faces)
- buffer_send_faces_vector(1,ioffset + ipoin2D) = accel_inner_core(iloop,iboolfaces_inner_core(ipoin2D,icount_faces))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ioffset + ipoin2D) = accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces))
- buffer_send_faces_vector(3,ioffset + ipoin2D) = accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces))
- endif
+ buffer_send_faces_vector(1,ioffset + ipoin2D) = accel_inner_core(1,iboolfaces_inner_core(ipoin2D,icount_faces))
+ buffer_send_faces_vector(2,ioffset + ipoin2D) = accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces))
+ buffer_send_faces_vector(3,ioffset + ipoin2D) = accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces))
enddo
#ifdef USE_MPI
- call MPI_SEND(buffer_send_faces_vector,NDIM_smaller_buffers*npoin2D_chunks_all,CUSTOM_MPI_TYPE,receiver,itag,MPI_COMM_WORLD,ier)
+ call MPI_SEND(buffer_send_faces_vector,NDIM*npoin2D_chunks_all,CUSTOM_MPI_TYPE,receiver,itag,MPI_COMM_WORLD,ier)
#endif
endif
@@ -541,24 +498,20 @@
ioffset = npoin2D_faces_crust_mantle(icount_faces)
#ifdef USE_MPI
- call MPI_RECV(buffer_received_faces_vector,NDIM_smaller_buffers*npoin2D_chunks_all,CUSTOM_MPI_TYPE,sender, &
+ call MPI_RECV(buffer_received_faces_vector,NDIM*npoin2D_chunks_all,CUSTOM_MPI_TYPE,sender, &
itag,MPI_COMM_WORLD,msg_status,ier)
#endif
do ipoin2D = 1,npoin2D_faces_crust_mantle(icount_faces)
- accel_crust_mantle(iloop,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = buffer_received_faces_vector(1,ipoin2D)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = buffer_received_faces_vector(2,ipoin2D)
- accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = buffer_received_faces_vector(3,ipoin2D)
- endif
+ accel_crust_mantle(1,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = buffer_received_faces_vector(1,ipoin2D)
+ accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = buffer_received_faces_vector(2,ipoin2D)
+ accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces)) = buffer_received_faces_vector(3,ipoin2D)
enddo
do ipoin2D = 1,npoin2D_faces_inner_core(icount_faces)
- accel_inner_core(iloop,iboolfaces_inner_core(ipoin2D,icount_faces)) = buffer_received_faces_vector(1,ioffset + ipoin2D)
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces)) = buffer_received_faces_vector(2,ioffset + ipoin2D)
- accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces)) = buffer_received_faces_vector(3,ioffset + ipoin2D)
- endif
+ accel_inner_core(1,iboolfaces_inner_core(ipoin2D,icount_faces)) = buffer_received_faces_vector(1,ioffset + ipoin2D)
+ accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces)) = buffer_received_faces_vector(2,ioffset + ipoin2D)
+ accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces)) = buffer_received_faces_vector(3,ioffset + ipoin2D)
enddo
endif
@@ -580,23 +533,19 @@
ioffset = npoin2D_faces_crust_mantle(icount_faces)
do ipoin2D = 1,npoin2D_faces_crust_mantle(icount_faces)
- buffer_send_faces_vector(1,ipoin2D) = accel_crust_mantle(iloop,iboolfaces_crust_mantle(ipoin2D,icount_faces))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ipoin2D) = accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces))
- buffer_send_faces_vector(3,ipoin2D) = accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces))
- endif
+ buffer_send_faces_vector(1,ipoin2D) = accel_crust_mantle(1,iboolfaces_crust_mantle(ipoin2D,icount_faces))
+ buffer_send_faces_vector(2,ipoin2D) = accel_crust_mantle(2,iboolfaces_crust_mantle(ipoin2D,icount_faces))
+ buffer_send_faces_vector(3,ipoin2D) = accel_crust_mantle(3,iboolfaces_crust_mantle(ipoin2D,icount_faces))
enddo
do ipoin2D = 1,npoin2D_faces_inner_core(icount_faces)
- buffer_send_faces_vector(1,ioffset + ipoin2D) = accel_inner_core(iloop,iboolfaces_inner_core(ipoin2D,icount_faces))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- buffer_send_faces_vector(2,ioffset + ipoin2D) = accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces))
- buffer_send_faces_vector(3,ioffset + ipoin2D) = accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces))
- endif
+ buffer_send_faces_vector(1,ioffset + ipoin2D) = accel_inner_core(1,iboolfaces_inner_core(ipoin2D,icount_faces))
+ buffer_send_faces_vector(2,ioffset + ipoin2D) = accel_inner_core(2,iboolfaces_inner_core(ipoin2D,icount_faces))
+ buffer_send_faces_vector(3,ioffset + ipoin2D) = accel_inner_core(3,iboolfaces_inner_core(ipoin2D,icount_faces))
enddo
#ifdef USE_MPI
- call MPI_SEND(buffer_send_faces_vector,NDIM_smaller_buffers*npoin2D_chunks_all,CUSTOM_MPI_TYPE,receiver,itag,MPI_COMM_WORLD,ier)
+ call MPI_SEND(buffer_send_faces_vector,NDIM*npoin2D_chunks_all,CUSTOM_MPI_TYPE,receiver,itag,MPI_COMM_WORLD,ier)
#endif
endif
@@ -605,8 +554,6 @@
! end of anti-deadlocking loop
enddo
- enddo ! of loop on iloop depending on NDIM_smaller_buffers
-
!----
!---- start MPI assembling corners
!----
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/attenuation_model.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/attenuation_model.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/attenuation_model.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -686,6 +686,7 @@
integer i,j,rmax
#ifdef USE_MPI
+ double precision, dimension(N_SLS+2) :: array_to_broadcast
integer :: ier
#endif
double precision scale_t
@@ -717,9 +718,15 @@
! Synch up after the Read
#ifdef USE_MPI
- call MPI_BCAST(AM_V%QT_c_source,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(tau_s,N_SLS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(AM_V%Qn,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ array_to_broadcast(1:N_SLS) = tau_s(:)
+ array_to_broadcast(N_SLS+1) = AM_V%QT_c_source
+ array_to_broadcast(N_SLS+2) = AM_V%Qn
+
+ call MPI_BCAST(array_to_broadcast,N_SLS+2,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+ tau_s(:) = array_to_broadcast(1:N_SLS)
+ AM_V%QT_c_source = array_to_broadcast(N_SLS+1)
+ AM_V%Qn = nint(array_to_broadcast(N_SLS+2))
#endif
if(myrank /= 0) then
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -247,12 +247,10 @@
! crustal_model_variables
type crustal_model_variables
sequence
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: thlr
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocp
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocs
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: dens
- character(len=2) abbreviation(NCAP_CRUST/2,NCAP_CRUST)
- character(len=2) code(NKEYS_CRUST)
+ real(kind=4) velocp(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) velocs(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) dens(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) thlr(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
end type crustal_model_variables
type (crustal_model_variables) CM_V
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_chunk_buffers.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_chunk_buffers.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_chunk_buffers.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -140,7 +140,7 @@
integer ichunk_receive,iproc_xi_receive,iproc_eta_receive
integer iproc_loop,iproc_xi_loop,iproc_eta_loop
integer iproc_xi_loop_inv,iproc_eta_loop_inv
- integer imember_corner
+ integer imember_corner,npoin2D_value_received
integer iregion_code
@@ -152,11 +152,13 @@
integer ibelm_xmin(NSPEC2DMAX_XMIN_XMAX),ibelm_xmax(NSPEC2DMAX_XMIN_XMAX)
integer ibelm_ymin(NSPEC2DMAX_YMIN_YMAX),ibelm_ymax(NSPEC2DMAX_YMIN_YMAX)
- integer npoin2D,npoin2D_send_local,npoin2D_receive_local
+ integer npoin2D
integer i,j,k,ispec,ispec2D,ipoin2D
#ifdef USE_MPI
+ integer :: sender,receiver
+ integer, dimension(MPI_STATUS_SIZE) :: msg_status
integer :: ier
#endif
@@ -747,29 +749,24 @@
!
! gather information about all the messages on all processes
+#ifdef USE_MPI
do imsg = 1,NUMMSGS_FACES
+! send or receive number of points for sender
+ if(myrank == iprocto_faces(imsg)) then
+ sender = iprocfrom_faces(imsg)
+ call MPI_RECV(npoin2D_value_received,1,MPI_INTEGER,sender,itag,MPI_COMM_WORLD,msg_status,ier)
+ if(npoin2D_value_received /= npoin2D_receive(imsg)) &
+ call exit_MPI(myrank,'incorrect number of points for sender/receiver pair detected')
+ endif
-! gather number of points for sender
- npoin2D_send_local = npoin2D_send(imsg)
-#ifdef USE_MPI
- call MPI_BCAST(npoin2D_send_local,1,MPI_INTEGER,iprocfrom_faces(imsg),MPI_COMM_WORLD,ier)
+ if(myrank == iprocfrom_faces(imsg)) then
+ receiver = iprocto_faces(imsg)
+ call MPI_SEND(npoin2D_send(imsg),1,MPI_INTEGER,receiver,itag,MPI_COMM_WORLD,ier)
+ endif
+ enddo
#endif
- if(myrank /= iprocfrom_faces(imsg)) npoin2D_send(imsg) = npoin2D_send_local
-! gather number of points for receiver
- npoin2D_receive_local = npoin2D_receive(imsg)
-#ifdef USE_MPI
- call MPI_BCAST(npoin2D_receive_local,1,MPI_INTEGER,iprocto_faces(imsg),MPI_COMM_WORLD,ier)
-#endif
- if(myrank /= iprocto_faces(imsg)) npoin2D_receive(imsg) = npoin2D_receive_local
-
- enddo
-
! check the number of points
- do imsg = 1,NUMMSGS_FACES
- if(npoin2D_send(imsg) /= npoin2D_receive(imsg)) &
- call exit_MPI(myrank,'incorrect number of points for sender/receiver pair detected')
- enddo
if(myrank == 0) then
write(IMAIN,*)
write(IMAIN,*) 'all the messages for chunk faces have the right size'
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -67,7 +67,7 @@
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!! DK DK to debug the second sorting routine
- logical, parameter :: DEBUG = .false.
+ logical, parameter :: DEBUG = .true.
!! DK DK added this for merged version
#ifdef USE_MPI
@@ -326,12 +326,10 @@
! crustal_model_variables
type crustal_model_variables
sequence
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: thlr
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocp
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocs
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: dens
- character(len=2) abbreviation(NCAP_CRUST/2,NCAP_CRUST)
- character(len=2) code(NKEYS_CRUST)
+ real(kind=4) velocp(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) velocs(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) dens(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) thlr(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
end type crustal_model_variables
type (crustal_model_variables) CM_V
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/crustal_model.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/crustal_model.f90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/crustal_model.f90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -26,38 +26,118 @@
!=====================================================================
!
-! read and smooth crust2.0 model
-! based on software routines provided with the crust2.0 model by Bassin et al.
+! read or evaluate smoothed crust2.0 model
!
- subroutine crustal_model(lat,lon,x,vp,vs,rho,moho,found_crust,CM_V)
+ subroutine get_smoothed_crust(lat,lon,x,vp,vs,rho,moho,found_crust,CM_V)
implicit none
+
include "constants.h"
+ logical found_crust
+ double precision lat,lon,x,vp,vs,rho,moho
+
+ double precision h_sed,h_uc
+ double precision x3,x4,x5,x6,x7,scaleval
+ double precision vps(3:7),vss(3:7),rhos(3:7),thicks(3:7)
+
+ double precision vps1(3:7),vss1(3:7),rhos1(3:7),thicks1(3:7)
+ double precision vps2(3:7),vss2(3:7),rhos2(3:7),thicks2(3:7)
+ double precision vps3(3:7),vss3(3:7),rhos3(3:7),thicks3(3:7)
+ double precision vps4(3:7),vss4(3:7),rhos4(3:7),thicks4(3:7)
+
! crustal_model_variables
type crustal_model_variables
sequence
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: thlr
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocp
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocs
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: dens
- character(len=2) abbreviation(NCAP_CRUST/2,NCAP_CRUST)
- character(len=2) code(NKEYS_CRUST)
+ real(kind=4) velocp(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) velocs(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) dens(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) thlr(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
end type crustal_model_variables
type (crustal_model_variables) CM_V
! crustal_model_variables
- double precision lat,lon,x,vp,vs,rho,moho
- logical found_crust
+ integer ilon_scaled,ilat_scaled
+ double precision lon_scaled,lat_scaled
+ double precision gamma_lon,gamma_lat
- double precision h_sed,h_uc
- double precision x3,x4,x5,x6,x7,scaleval
- double precision vps(NLAYERS_CRUST),vss(NLAYERS_CRUST),rhos(NLAYERS_CRUST),thicks(NLAYERS_CRUST)
+! scale interval and make all values positive
+ lon_scaled = lon * NFACTOR_CRUST + NLON_CRUST
+ lat_scaled = lat * NFACTOR_CRUST + NLAT_CRUST
- call crust(lat,lon,vps,vss,rhos,thicks,CM_V%abbreviation,CM_V%code,CM_V%thlr,CM_V%velocp,CM_V%velocs,CM_V%dens)
+! select the right cell to use in smoothed file
+ ilon_scaled = int(lon_scaled)
+ ilat_scaled = int(lat_scaled)
+! compute coordinates in interpolation cell
+ gamma_lon = lon_scaled - dble(ilon_scaled)
+ gamma_lat = lat_scaled - dble(ilat_scaled)
+
+! prevent edge effects
+ if(ilon_scaled < 0) then
+ ilon_scaled = 0
+ gamma_lon = 0.d0
+ endif
+
+ if(ilat_scaled < 0) then
+ ilat_scaled = 0
+ gamma_lat = 0.d0
+ endif
+
+ if(ilon_scaled >= 2*NLON_CRUST) then
+ ilon_scaled = 2*NLON_CRUST - 1
+ gamma_lon = 1.d0
+ endif
+
+ if(ilat_scaled >= 2*NLAT_CRUST) then
+ ilat_scaled = 2*NLAT_CRUST - 1
+ gamma_lat = 1.d0
+ endif
+
+! copy four corners of interpolation cell
+ vps1(:) = dble(CM_V%velocp(ilon_scaled,ilat_scaled,:))
+ vss1(:) = dble(CM_V%velocs(ilon_scaled,ilat_scaled,:))
+ rhos1(:) = dble(CM_V%dens(ilon_scaled,ilat_scaled,:))
+ thicks1(:) = dble(CM_V%thlr(ilon_scaled,ilat_scaled,:))
+
+ vps2(:) = dble(CM_V%velocp(ilon_scaled+1,ilat_scaled,:))
+ vss2(:) = dble(CM_V%velocs(ilon_scaled+1,ilat_scaled,:))
+ rhos2(:) = dble(CM_V%dens(ilon_scaled+1,ilat_scaled,:))
+ thicks2(:) = dble(CM_V%thlr(ilon_scaled+1,ilat_scaled,:))
+
+ vps3(:) = dble(CM_V%velocp(ilon_scaled+1,ilat_scaled+1,:))
+ vss3(:) = dble(CM_V%velocs(ilon_scaled+1,ilat_scaled+1,:))
+ rhos3(:) = dble(CM_V%dens(ilon_scaled+1,ilat_scaled+1,:))
+ thicks3(:) = dble(CM_V%thlr(ilon_scaled+1,ilat_scaled+1,:))
+
+ vps4(:) = dble(CM_V%velocp(ilon_scaled,ilat_scaled+1,:))
+ vss4(:) = dble(CM_V%velocs(ilon_scaled,ilat_scaled+1,:))
+ rhos4(:) = dble(CM_V%dens(ilon_scaled,ilat_scaled+1,:))
+ thicks4(:) = dble(CM_V%thlr(ilon_scaled,ilat_scaled+1,:))
+
+! perform interpolation from the four corners
+ vps(:) = (1.d0-gamma_lon)*(1.d0-gamma_lat)*vps1(:) + &
+ gamma_lon*(1.d0-gamma_lat)*vps2(:) + &
+ gamma_lon*gamma_lat*vps3(:) + &
+ (1.d0-gamma_lon)*gamma_lat*vps4(:)
+
+ vss(:) = (1.d0-gamma_lon)*(1.d0-gamma_lat)*vss1(:) + &
+ gamma_lon*(1.d0-gamma_lat)*vss2(:) + &
+ gamma_lon*gamma_lat*vss3(:) + &
+ (1.d0-gamma_lon)*gamma_lat*vss4(:)
+
+ rhos(:) = (1.d0-gamma_lon)*(1.d0-gamma_lat)*rhos1(:) + &
+ gamma_lon*(1.d0-gamma_lat)*rhos2(:) + &
+ gamma_lon*gamma_lat*rhos3(:) + &
+ (1.d0-gamma_lon)*gamma_lat*rhos4(:)
+
+ thicks(:) = (1.d0-gamma_lon)*(1.d0-gamma_lat)*thicks1(:) + &
+ gamma_lon*(1.d0-gamma_lat)*thicks2(:) + &
+ gamma_lon*gamma_lat*thicks3(:) + &
+ (1.d0-gamma_lon)*gamma_lat*thicks4(:)
+
x3 = (R_EARTH-thicks(3)*1000.0d0)/R_EARTH
h_sed = thicks(3) + thicks(4)
x4 = (R_EARTH-h_sed*1000.0d0)/R_EARTH
@@ -67,6 +147,7 @@
x7 = (R_EARTH-(h_uc+thicks(6)+thicks(7))*1000.0d0)/R_EARTH
found_crust = .true.
+
if(x > x3 .and. INCLUDE_SEDIMENTS_CRUST) then
vp = vps(3)
vs = vss(3)
@@ -91,276 +172,77 @@
found_crust = .false.
endif
- if (found_crust) then
+ if(found_crust) then
+
! non-dimensionalize
scaleval = dsqrt(PI*GRAV*RHOAV)
vp = vp*1000.0d0/(R_EARTH*scaleval)
vs = vs*1000.0d0/(R_EARTH*scaleval)
rho = rho*1000.0d0/RHOAV
moho = (h_uc+thicks(6)+thicks(7))*1000.0d0/R_EARTH
- endif
- end subroutine crustal_model
+ endif
+ end subroutine get_smoothed_crust
+
!---------------------------
- subroutine read_crustal_model(CM_V)
+ subroutine read_smoothed_3D_crustal_model(CM_V)
implicit none
+
include "constants.h"
! crustal_model_variables
type crustal_model_variables
sequence
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: thlr
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocp
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocs
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: dens
- character(len=2) abbreviation(NCAP_CRUST/2,NCAP_CRUST)
- character(len=2) code(NKEYS_CRUST)
- end type crustal_model_variables
+ real(kind=4) velocp(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) velocs(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) dens(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) thlr(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ end type crustal_model_variables
type (crustal_model_variables) CM_V
! crustal_model_variables
-! local variables
- integer i
- integer ila,icolat
- integer ikey
+ integer ilon,ilat
- double precision h_moho_min,h_moho_max
+ open(unit=1,file='DATA/crust2.0/smoothed_crust2.0.dat',status='old',action='read')
- character(len=150) CNtype2, CNtype2_key_modif
+! loop on latitude and longitude
+ do ilat = 0,2*NLAT_CRUST
+ do ilon = 0,2*NLON_CRUST
- call get_value_string(CNtype2, 'model.CNtype2', 'DATA/crust2.0/CNtype2.txt')
- call get_value_string(CNtype2_key_modif, 'model.CNtype2_key_modif', 'DATA/crust2.0/CNtype2_key_modif.txt')
+! we only saved layers 3 to 7 because 1, 2 and 8 are never used in the mesher
- open(unit=1,file=CNtype2,status='old',action='read')
- do ila=1,NCAP_CRUST/2
- read(1,*) icolat,(CM_V%abbreviation(ila,i),i=1,NCAP_CRUST)
- enddo
- close(1)
+ read(1,*) CM_V%velocp(ilon,ilat,3)
+ read(1,*) CM_V%velocp(ilon,ilat,4)
+ read(1,*) CM_V%velocp(ilon,ilat,5)
+ read(1,*) CM_V%velocp(ilon,ilat,6)
+ read(1,*) CM_V%velocp(ilon,ilat,7)
- open(unit=1,file=CNtype2_key_modif,status='old',action='read')
- h_moho_min=HUGEVAL
- h_moho_max=-HUGEVAL
- do ikey=1,NKEYS_CRUST
- read (1,"(a2)") CM_V%code(ikey)
- read (1,*) (CM_V%velocp(ikey,i),i=1,NLAYERS_CRUST)
- read (1,*) (CM_V%velocs(ikey,i),i=1,NLAYERS_CRUST)
- read (1,*) (CM_V%dens(ikey,i),i=1,NLAYERS_CRUST)
- read (1,*) (CM_V%thlr(ikey,i),i=1,NLAYERS_CRUST-1),CM_V%thlr(ikey,NLAYERS_CRUST)
- if(CM_V%thlr(ikey,NLAYERS_CRUST) > h_moho_max) h_moho_max=CM_V%thlr(ikey,NLAYERS_CRUST)
- if(CM_V%thlr(ikey,NLAYERS_CRUST) < h_moho_min) h_moho_min=CM_V%thlr(ikey,NLAYERS_CRUST)
- enddo
- close(1)
+ read(1,*) CM_V%velocs(ilon,ilat,3)
+ read(1,*) CM_V%velocs(ilon,ilat,4)
+ read(1,*) CM_V%velocs(ilon,ilat,5)
+ read(1,*) CM_V%velocs(ilon,ilat,6)
+ read(1,*) CM_V%velocs(ilon,ilat,7)
- if(h_moho_min == HUGEVAL .or. h_moho_max == -HUGEVAL) &
- stop 'incorrect moho depths in read_3D_crustal_model'
+ read(1,*) CM_V%dens(ilon,ilat,3)
+ read(1,*) CM_V%dens(ilon,ilat,4)
+ read(1,*) CM_V%dens(ilon,ilat,5)
+ read(1,*) CM_V%dens(ilon,ilat,6)
+ read(1,*) CM_V%dens(ilon,ilat,7)
- end subroutine read_crustal_model
+ read(1,*) CM_V%thlr(ilon,ilat,3)
+ read(1,*) CM_V%thlr(ilon,ilat,4)
+ read(1,*) CM_V%thlr(ilon,ilat,5)
+ read(1,*) CM_V%thlr(ilon,ilat,6)
+ read(1,*) CM_V%thlr(ilon,ilat,7)
-!---------------------------
-
- subroutine crust(lat,lon,velp,vels,rho,thick,abbreviation,code,thlr,velocp,velocs,dens)
-
-! crustal vp and vs in km/s, layer thickness in km
-! crust2.0 is smoothed with a cap of size CAP using NTHETA points
-! in the theta direction and NPHI in the phi direction.
-! The cap is rotated to the North Pole.
-
- implicit none
- include "constants.h"
-
- integer, parameter :: NTHETA = 2
- integer, parameter :: NPHI = 10
- double precision, parameter :: CAP = 2.0d0*PI/180.0d0
-
-! argument variables
- double precision lat,lon
- double precision rho(NLAYERS_CRUST),thick(NLAYERS_CRUST),velp(NLAYERS_CRUST),vels(NLAYERS_CRUST)
- double precision thlr(NKEYS_CRUST,NLAYERS_CRUST),velocp(NKEYS_CRUST,NLAYERS_CRUST)
- double precision velocs(NKEYS_CRUST,NLAYERS_CRUST),dens(NKEYS_CRUST,NLAYERS_CRUST)
- character(len=2) code(NKEYS_CRUST),abbreviation(NCAP_CRUST/2,NCAP_CRUST)
-
-! local variables
- integer i,j,k,icolat,ilon,ierr
- integer itheta,iphi,npoints
- double precision theta,phi,sint,cost,sinp,cosp,dtheta,dphi,cap_area,wght,total
- double precision r_rot,theta_rot,phi_rot
- double precision rotation_matrix(3,3),x(3),xc(3)
- double precision xlon(NTHETA*NPHI),xlat(NTHETA*NPHI),weight(NTHETA*NPHI)
- double precision rhol(NLAYERS_CRUST),thickl(NLAYERS_CRUST),velpl(NLAYERS_CRUST),velsl(NLAYERS_CRUST)
- character(len=2) crustaltype
-
-! get integer colatitude and longitude of crustal cap
-! -90<lat<90 -180<lon<180
- if(lat > 90.0d0 .or. lat < -90.0d0 .or. lon > 180.0d0 .or. lon < -180.0d0) &
- stop 'error in latitude/longitude range in crust'
- if(lat==90.0d0) lat=89.9999d0
- if(lat==-90.0d0) lat=-89.9999d0
- if(lon==180.0d0) lon=179.9999d0
- if(lon==-180.0d0) lon=-179.9999d0
-
- call icolat_ilon(lat,lon,icolat,ilon)
- crustaltype=abbreviation(icolat,ilon)
- call get_crust_structure(crustaltype,velp,vels,rho,thick, &
- code,thlr,velocp,velocs,dens,ierr)
-
-! uncomment the following line to use crust2.0 as is, without smoothing
-!
-! return
-
- theta = (90.0-lat)*PI/180.0
- phi = lon*PI/180.0
-
- sint = sin(theta)
- cost = cos(theta)
- sinp = sin(phi)
- cosp = cos(phi)
-
-! set up rotation matrix to go from cap at North pole
-! to cap around point of interest
- rotation_matrix(1,1) = cosp*cost
- rotation_matrix(1,2) = -sinp
- rotation_matrix(1,3) = cosp*sint
- rotation_matrix(2,1) = sinp*cost
- rotation_matrix(2,2) = cosp
- rotation_matrix(2,3) = sinp*sint
- rotation_matrix(3,1) = -sint
- rotation_matrix(3,2) = 0.0
- rotation_matrix(3,3) = cost
-
- dtheta = CAP/dble(NTHETA)
- dphi = 2.0*PI/dble(NPHI)
- cap_area = 2.0*PI*(1.0-cos(CAP))
-
-! integrate over a cap at the North pole
- i = 0
- total = 0.0
- do itheta = 1,NTHETA
-
- theta = 0.5*dble(2*itheta-1)*CAP/dble(NTHETA)
- cost = cos(theta)
- sint = sin(theta)
- wght = sint*dtheta*dphi/cap_area
-
- do iphi = 1,NPHI
-
- i = i+1
-! get the weight associated with this integration point (same for all phi)
- weight(i) = wght
- total = total + weight(i)
- phi = dble(2*iphi-1)*PI/dble(NPHI)
- cosp = cos(phi)
- sinp = sin(phi)
-! x,y,z coordinates of integration point in cap at North pole
- xc(1) = sint*cosp
- xc(2) = sint*sinp
- xc(3) = cost
-! get x,y,z coordinates in cap around point of interest
- do j=1,3
- x(j) = 0.0
- do k=1,3
- x(j) = x(j)+rotation_matrix(j,k)*xc(k)
- enddo
- enddo
-! get latitude and longitude (degrees) of integration point
- call xyz_2_rthetaphi_dble(x(1),x(2),x(3),r_rot,theta_rot,phi_rot)
- call reduce(theta_rot,phi_rot)
- xlat(i) = (PI/2.0-theta_rot)*180.0/PI
- xlon(i) = phi_rot*180.0/PI
- if(xlon(i) > 180.0) xlon(i) = xlon(i)-360.0
-
- enddo
-
enddo
-
- if(abs(total-1.0) > 0.001) stop 'error in cap integration for crust2.0'
-
- npoints = i
-
- do j=1,NLAYERS_CRUST
- rho(j)=0.0d0
- thick(j)=0.0d0
- velp(j)=0.0d0
- vels(j)=0.0d0
enddo
- do i=1,npoints
- call icolat_ilon(xlat(i),xlon(i),icolat,ilon)
- crustaltype=abbreviation(icolat,ilon)
- call get_crust_structure(crustaltype,velpl,velsl,rhol,thickl, &
- code,thlr,velocp,velocs,dens,ierr)
- if(ierr /= 0) stop 'error in routine get_crust_structure'
- do j=1,NLAYERS_CRUST
- rho(j)=rho(j)+weight(i)*rhol(j)
- thick(j)=thick(j)+weight(i)*thickl(j)
- velp(j)=velp(j)+weight(i)*velpl(j)
- vels(j)=vels(j)+weight(i)*velsl(j)
- enddo
- enddo
+ close(1)
- end subroutine crust
+ end subroutine read_smoothed_3D_crustal_model
-!------------------------------------------------------
-
- subroutine icolat_ilon(xlat,xlon,icolat,ilon)
-
- implicit none
-
-
-! argument variables
- double precision xlat,xlon
- integer icolat,ilon
-
- if(xlat > 90.0d0 .or. xlat < -90.0d0 .or. xlon > 180.0d0 .or. xlon < -180.0d0) &
- stop 'error in latitude/longitude range in icolat_ilon'
- icolat=int(1+((90.d0-xlat)/2.d0))
- if(icolat == 91) icolat=90
- ilon=int(1+((180.d0+xlon)/2.d0))
- if(ilon == 181) ilon=1
-
- if(icolat>90 .or. icolat<1) stop 'error in routine icolat_ilon'
- if(ilon<1 .or. ilon>180) stop 'error in routine icolat_ilon'
-
- end subroutine icolat_ilon
-
-!---------------------------------------------------------------------
-
- subroutine get_crust_structure(type,vptyp,vstyp,rhtyp,thtp, &
- code,thlr,velocp,velocs,dens,ierr)
-
- implicit none
- include "constants.h"
-
-! argument variables
- integer ierr
- double precision rhtyp(NLAYERS_CRUST),thtp(NLAYERS_CRUST)
- double precision vptyp(NLAYERS_CRUST),vstyp(NLAYERS_CRUST)
- character(len=2) type,code(NKEYS_CRUST)
- double precision thlr(NKEYS_CRUST,NLAYERS_CRUST),velocp(NKEYS_CRUST,NLAYERS_CRUST)
- double precision velocs(NKEYS_CRUST,NLAYERS_CRUST),dens(NKEYS_CRUST,NLAYERS_CRUST)
-
-! local variables
- integer i,ikey
-
- ierr=1
- do ikey=1,NKEYS_CRUST
- if (code(ikey) == type) then
- do i=1,NLAYERS_CRUST
- vptyp(i)=velocp(ikey,i)
- vstyp(i)=velocs(ikey,i)
- rhtyp(i)=dens(ikey,i)
- enddo
- do i=1,NLAYERS_CRUST-1
- thtp(i)=thlr(ikey,i)
- enddo
-! get distance to Moho from the bottom of the ocean or the ice
- thtp(NLAYERS_CRUST)=thlr(ikey,NLAYERS_CRUST)-thtp(1)-thtp(2)
- ierr=0
- endif
- enddo
-
- end subroutine get_crust_structure
-
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/extract_all_seismos_from_large_file.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/extract_all_seismos_from_large_file.f90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/extract_all_seismos_from_large_file.f90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -33,15 +33,17 @@
implicit none
- include "constants.h"
+ include "setup/constants.h"
! number of seismogram files stored in the unique large file
integer, parameter :: N_COMPONENTS = 3
- integer, parameter :: NREC = 344 * N_COMPONENTS
+ integer, parameter :: NREC = 35 * N_COMPONENTS
! number of time steps in each seismogram file
- integer, parameter :: NSTEP = 66800
+ integer, parameter :: NSTEP = 10100
+ logical, parameter :: USE_BINARY_FOR_LARGE_FILE = .true.
+
integer :: irec,istep,irepeat
real :: time,U_value
@@ -49,9 +51,9 @@
! open the large seismogram file
if(USE_BINARY_FOR_LARGE_FILE) then
- open(unit=30,file='OUTPUT_FILES/all_seismograms_d.bin',status='old',form='unformatted',action='read')
+ open(unit=30,file='OUTPUT_FILES/all_seismograms.bin',status='old',form='unformatted',action='read')
else
- open(unit=30,file='all_seismograms.ascii',status='old',action='read')
+ open(unit=30,file='OUTPUT_FILES/all_seismograms.ascii',status='old',action='read')
endif
! loop on all the seismogram files
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_event_info.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_event_info.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_event_info.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -53,6 +53,7 @@
integer :: i
#ifdef USE_MPI
+ double precision, dimension(14) :: array_to_broadcast
integer :: ier
#endif
@@ -77,22 +78,44 @@
! broadcast the information read on the master to the nodes
#ifdef USE_MPI
- call MPI_BCAST(yr,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(jda,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ho,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(mi,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSOURCES,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+! integer
+ array_to_broadcast( 1) = yr
+ array_to_broadcast( 2) = jda
+ array_to_broadcast( 3) = ho
+ array_to_broadcast( 4) = mi
+ array_to_broadcast( 5) = NSOURCES
- call MPI_BCAST(sec,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(t_cmt,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(elat,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(elon,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(depth,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(cmt_lat,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(cmt_lon,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(cmt_depth,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(cmt_hdur,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+! double precision
+ array_to_broadcast( 6) = sec
+ array_to_broadcast( 7) = t_cmt
+ array_to_broadcast( 8) = elat
+ array_to_broadcast( 9) = elon
+ array_to_broadcast(10) = depth
+ array_to_broadcast(11) = cmt_lat
+ array_to_broadcast(12) = cmt_lon
+ array_to_broadcast(13) = cmt_depth
+ array_to_broadcast(14) = cmt_hdur
+ call MPI_BCAST(array_to_broadcast,14,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+! integer
+ yr = nint(array_to_broadcast( 1))
+ jda = nint(array_to_broadcast( 2))
+ ho = nint(array_to_broadcast( 3))
+ mi = nint(array_to_broadcast( 4))
+ NSOURCES = nint(array_to_broadcast( 5))
+
+! double precision
+ sec = array_to_broadcast( 6)
+ t_cmt = array_to_broadcast( 7)
+ elat = array_to_broadcast( 8)
+ elon = array_to_broadcast( 9)
+ depth = array_to_broadcast(10)
+ cmt_lat = array_to_broadcast(11)
+ cmt_lon = array_to_broadcast(12)
+ cmt_depth = array_to_broadcast(13)
+ cmt_hdur = array_to_broadcast(14)
+
call MPI_BCAST(ename,12,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
#endif
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -252,12 +252,10 @@
! crustal_model_variables
type crustal_model_variables
sequence
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: thlr
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocp
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocs
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: dens
- character(len=2) abbreviation(NCAP_CRUST/2,NCAP_CRUST)
- character(len=2) code(NKEYS_CRUST)
+ real(kind=4) velocp(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) velocs(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) dens(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) thlr(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
end type crustal_model_variables
type (crustal_model_variables) CM_V
@@ -787,7 +785,7 @@
lat=(PI/2.0d0-theta)*180.0d0/PI
lon=phi*180.0d0/PI
if(lon>180.0d0) lon=lon-360.0d0
- call crustal_model(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,CM_V)
+ call get_smoothed_crust(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,CM_V)
if (found_crust) then
vpv=vpc
vph=vpc
@@ -824,7 +822,7 @@
lat=(PI/2.0d0-theta)*180.0d0/PI
lon=phi*180.0d0/PI
if(lon>180.0d0) lon=lon-360.0d0
- call crustal_model(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,CM_V)
+ call get_smoothed_crust(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,CM_V)
if (found_crust) then
vpv=vpc
vph=vpc
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_receivers.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_receivers.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_receivers.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -84,6 +84,7 @@
integer i,j,k,ispec,iglob,imin,imax,jmin,jmax,kmin,kmax
#ifdef USE_MPI
+ double precision, dimension(nrec,8) :: array_to_broadcast
integer :: ier
#endif
@@ -210,10 +211,17 @@
call MPI_BCAST(station_name,nrec*MAX_LENGTH_STATION_NAME,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(network_name,nrec*MAX_LENGTH_NETWORK_NAME,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stlat,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stlon,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stele,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stbur,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+ array_to_broadcast(:,1) = stlat(:)
+ array_to_broadcast(:,2) = stlon(:)
+ array_to_broadcast(:,3) = stele(:)
+ array_to_broadcast(:,4) = stbur(:)
+
+ call MPI_BCAST(array_to_broadcast,4*nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+ stlat(:) = array_to_broadcast(:,1)
+ stlon(:) = array_to_broadcast(:,2)
+ stele(:) = array_to_broadcast(:,3)
+ stbur(:) = array_to_broadcast(:,4)
#endif
! loop on all the stations to locate them in the mesh
@@ -688,19 +696,37 @@
#ifdef USE_MPI
call MPI_BCAST(nrec,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(islice_selected_rec,nrec,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ispec_selected_rec,nrec,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(xi_receiver,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(eta_receiver,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(gamma_receiver,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
-
call MPI_BCAST(station_name,nrec*MAX_LENGTH_STATION_NAME,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(network_name,nrec*MAX_LENGTH_NETWORK_NAME,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stlat,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stlon,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(stele,nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(nu,nrec*3*3,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+! integer
+ array_to_broadcast(1:nrec,1) = islice_selected_rec(1:nrec)
+ array_to_broadcast(1:nrec,2) = ispec_selected_rec(1:nrec)
+
+! double precision
+ array_to_broadcast(1:nrec,3) = xi_receiver(1:nrec)
+ array_to_broadcast(1:nrec,4) = eta_receiver(1:nrec)
+ array_to_broadcast(1:nrec,5) = gamma_receiver(1:nrec)
+ array_to_broadcast(1:nrec,6) = stlat(1:nrec)
+ array_to_broadcast(1:nrec,7) = stlon(1:nrec)
+ array_to_broadcast(1:nrec,8) = stele(1:nrec)
+
+ call MPI_BCAST(array_to_broadcast,8*nrec,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+! integer
+ islice_selected_rec(1:nrec) = nint(array_to_broadcast(1:nrec,1))
+ ispec_selected_rec(1:nrec) = nint(array_to_broadcast(1:nrec,2))
+
+! double precision
+ xi_receiver(1:nrec) = array_to_broadcast(1:nrec,3)
+ eta_receiver(1:nrec) = array_to_broadcast(1:nrec,4)
+ gamma_receiver(1:nrec) = array_to_broadcast(1:nrec,5)
+ stlat(1:nrec) = array_to_broadcast(1:nrec,6)
+ stlon(1:nrec) = array_to_broadcast(1:nrec,7)
+ stele(1:nrec) = array_to_broadcast(1:nrec,8)
+
#endif
end subroutine locate_receivers
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_sources.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_sources.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/locate_sources.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -87,6 +87,7 @@
integer i,j,k,ispec,iglob
#ifdef USE_MPI
+ double precision, dimension(NSOURCES,16) :: array_to_broadcast
integer :: ier
#endif
@@ -186,20 +187,55 @@
if(myrank == 0) call get_cmt(yr,jda,ho,mi,sec,t_cmt,hdur,lat,long,depth,moment_tensor,DT,NSOURCES)
! broadcast the information read on the master to the nodes
#ifdef USE_MPI
- call MPI_BCAST(yr,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(jda,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ho,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(mi,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(sec,1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+! integer
+ array_to_broadcast(1,1) = yr
+ array_to_broadcast(1,2) = jda
+ array_to_broadcast(1,3) = ho
+ array_to_broadcast(1,4) = mi
- call MPI_BCAST(t_cmt,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(hdur,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(lat,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(long,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(depth,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+! double precision variable
+ array_to_broadcast(1,5) = sec
- call MPI_BCAST(moment_tensor,6*NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+! double precision arrays
+ array_to_broadcast(:, 6) = t_cmt(:)
+ array_to_broadcast(:, 7) = hdur(:)
+ array_to_broadcast(:, 8) = lat(:)
+ array_to_broadcast(:, 9) = long(:)
+ array_to_broadcast(:,10) = depth(:)
+
+ array_to_broadcast(:,11) = moment_tensor(1,:)
+ array_to_broadcast(:,12) = moment_tensor(2,:)
+ array_to_broadcast(:,13) = moment_tensor(3,:)
+ array_to_broadcast(:,14) = moment_tensor(4,:)
+ array_to_broadcast(:,15) = moment_tensor(5,:)
+ array_to_broadcast(:,16) = moment_tensor(6,:)
+
+ call MPI_BCAST(array_to_broadcast,16*NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+! integer
+ yr = nint(array_to_broadcast(1,1))
+ jda = nint(array_to_broadcast(1,2))
+ ho = nint(array_to_broadcast(1,3))
+ mi = nint(array_to_broadcast(1,4))
+
+! double precision variable
+ sec = array_to_broadcast(1,5)
+
+! double precision arrays
+ t_cmt(:) = array_to_broadcast(:, 6)
+ hdur(:) = array_to_broadcast(:, 7)
+ lat(:) = array_to_broadcast(:, 8)
+ long(:) = array_to_broadcast(:, 9)
+ depth(:) = array_to_broadcast(:,10)
+
+ moment_tensor(1,:) = array_to_broadcast(:,11)
+ moment_tensor(2,:) = array_to_broadcast(:,12)
+ moment_tensor(3,:) = array_to_broadcast(:,13)
+ moment_tensor(4,:) = array_to_broadcast(:,14)
+ moment_tensor(5,:) = array_to_broadcast(:,15)
+ moment_tensor(6,:) = array_to_broadcast(:,16)
+
#endif
! define topology of the control element
@@ -680,11 +716,25 @@
! main process broadcasts the results to all the slices
#ifdef USE_MPI
- call MPI_BCAST(islice_selected_source,NSOURCES,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ispec_selected_source,NSOURCES,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(xi_source,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(eta_source,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(gamma_source,NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+! integer
+ array_to_broadcast(:,1) = islice_selected_source(:)
+ array_to_broadcast(:,2) = ispec_selected_source(:)
+
+! double precision
+ array_to_broadcast(:,3) = xi_source(:)
+ array_to_broadcast(:,4) = eta_source(:)
+ array_to_broadcast(:,5) = gamma_source(:)
+
+ call MPI_BCAST(array_to_broadcast,5*NSOURCES,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+
+! integer
+ islice_selected_source(:) = nint(array_to_broadcast(:,1))
+ ispec_selected_source(:) = nint(array_to_broadcast(:,2))
+
+! double precision
+ xi_source(:) = array_to_broadcast(:,3)
+ eta_source(:) = array_to_broadcast(:,4)
+ gamma_source(:) = array_to_broadcast(:,5)
#endif
! elapsed time since beginning of source detection
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/main_program.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/main_program.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/main_program.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -323,11 +323,10 @@
integer, dimension(NGLOB2DMAX_YMIN_YMAX_IC) :: iboolleft_eta_inner_core,iboolright_eta_inner_core
#endif
- integer :: npoin2D_max_all,NDIM_smaller_buffers
+ integer :: npoin2D_max_all
! receiver information
- integer :: nrec,ios
- character(len=150) :: STATIONS,rec_filename,dummystring
+ integer :: nrec
!---- arrays to assemble between chunks
@@ -361,6 +360,32 @@
type (attenuation_model_variables) AM_V
! attenuation_model_variables
+! arrays from BCAST
+ integer, dimension(NVALUES_bcast_integer) :: bcast_integer
+ double precision, dimension(NVALUES_bcast_double_precision) :: bcast_double_precision
+ logical, dimension(NVALUES_bcast_logical) :: bcast_logical
+
+ character(len=150) MODEL
+
+! computed in read_compute_parameters
+ integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: ner,ratio_sampling_array
+ integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: doubling_index
+ double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: r_bottom,r_top
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_layer_has_a_doubling
+ double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: rmins,rmaxs
+
+! this for all the regions
+ integer, dimension(MAX_NUM_REGIONS) :: NSPEC_computed, &
+ NSPEC2D_XI, NSPEC2D_ETA, &
+ NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX, &
+ NSPEC2D_BOTTOM,NSPEC2D_TOP, &
+ NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
+ NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX, &
+ NGLOB_computed
+
+ integer, dimension(NB_SQUARE_CORNERS,NB_CUT_CASE) :: DIFF_NSPEC1D_RADIAL
+ integer, dimension(NB_SQUARE_EDGES_ONEDIR,NB_CUT_CASE) :: DIFF_NSPEC2D_XI,DIFF_NSPEC2D_ETA
+
! use equivalence statements to reduce total memory size
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle,&
@@ -462,7 +487,11 @@
npoin2D_faces_inner_core,npoin2D_xi_inner_core,npoin2D_eta_inner_core, &
#endif
rmass_ocean_load,normal_top_crust_mantle,ibelm_top_crust_mantle,AM_V, &
- locval,ifseg,copy_ibool_ori,mask_ibool,xstore,ystore,zstore)
+ locval,ifseg,copy_ibool_ori,mask_ibool,xstore,ystore,zstore,nrec, &
+ bcast_integer,bcast_double_precision,bcast_logical,MODEL,ner,ratio_sampling_array,doubling_index, &
+ r_bottom,r_top,rmins,rmaxs,this_layer_has_a_doubling,NSPEC_computed,NSPEC2D_XI,NSPEC2D_ETA, &
+ NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX,NSPEC2D_BOTTOM,NSPEC2D_TOP,NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
+ NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,NGLOB_computed,DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_ETA,DIFF_NSPEC2D_XI)
! synchronize all the processes to make sure everybody has finished creating the mesh
#ifdef USE_MPI
@@ -490,30 +519,7 @@
! size of buffers is the sum of two sizes because we handle two regions in the same MPI call
npoin2D_max_all = max(maxval(npoin2D_xi_crust_mantle(:) + npoin2D_xi_inner_core(:)), &
maxval(npoin2D_eta_crust_mantle(:) + npoin2D_eta_inner_core(:)))
- if(FEWER_MESSAGES_LARGER_BUFFERS) then
- NDIM_smaller_buffers = NDIM
- else
- NDIM_smaller_buffers = 1
- endif
-! read the number of receivers
- rec_filename = 'DATA/STATIONS'
- call get_value_string(STATIONS, 'solver.STATIONS', rec_filename)
-! get total number of receivers
- if(myrank == 0) then
- open(unit=IIN,file=STATIONS,iostat=ios,status='old',action='read')
- nrec = 0
- do while(ios == 0)
- read(IIN,"(a)",iostat=ios) dummystring
- if(ios == 0) nrec = nrec + 1
- enddo
- close(IIN)
- endif
-! broadcast the information read on the master to the nodes
-#ifdef USE_MPI
- call MPI_BCAST(nrec,1,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
-#endif
-
!! DK DK for the merged version, solver inserted here
call specfem3D(myrank,sizeprocs,ichunk_slice,iproc_xi_slice,iproc_eta_slice,NSOURCES, &
NTSTEP_BETWEEN_OUTPUT_SEISMOS,ibool_crust_mantle,ibool_outer_core,ibool_inner_core, &
@@ -523,7 +529,7 @@
kappahstore_crust_mantle,muhstore_crust_mantle,eta_anisostore_crust_mantle,kappavstore_inner_core,muvstore_inner_core, &
rmass_crust_mantle,rmass_outer_core,rmass_inner_core,rmass_ocean_load, &
#ifdef USE_MPI
- NDIM_smaller_buffers,npoin2D_max_all,nrec,addressing,ibathy_topo, &
+ npoin2D_max_all,nrec,addressing,ibathy_topo, &
ibelm_xmin_inner_core,ibelm_xmax_inner_core,ibelm_ymin_inner_core,ibelm_ymax_inner_core,ibelm_bottom_inner_core, &
iboolleft_xi_crust_mantle,iboolright_xi_crust_mantle, iboolleft_eta_crust_mantle,iboolright_eta_crust_mantle, &
iboolleft_xi_outer_core,iboolright_xi_outer_core,iboolleft_eta_outer_core,iboolright_eta_outer_core, &
@@ -539,7 +545,11 @@
#endif
AM_V,xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle,&
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle, &
- gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,displ_crust_mantle)
+ gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,displ_crust_mantle, &
+ bcast_integer,bcast_double_precision,bcast_logical,MODEL,ner,ratio_sampling_array,doubling_index, &
+ r_bottom,r_top,rmins,rmaxs,this_layer_has_a_doubling,NSPEC_computed,NSPEC2D_XI,NSPEC2D_ETA, &
+ NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX,NSPEC2D_BOTTOM,NSPEC2D_TOP,NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
+ NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,NGLOB_computed,DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_ETA,DIFF_NSPEC2D_XI)
! synchronize all the processes to make sure everybody has finished
#ifdef USE_MPI
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -55,7 +55,11 @@
npoin2D_faces_inner_core,npoin2D_xi_inner_core,npoin2D_eta_inner_core, &
#endif
rmass_ocean_load,normal_top_crust_mantle,ibelm_top_crust_mantle,AM_V, &
- locval,ifseg,copy_ibool_ori,mask_ibool,xstore,ystore,zstore)
+ locval,ifseg,copy_ibool_ori,mask_ibool,xstore,ystore,zstore,nrec, &
+ bcast_integer,bcast_double_precision,bcast_logical,MODEL,ner,ratio_sampling_array,doubling_index, &
+ r_bottom,r_top,rmins,rmaxs,this_layer_has_a_doubling,NSPEC,NSPEC2D_XI,NSPEC2D_ETA, &
+ NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX,NSPEC2D_BOTTOM,NSPEC2D_TOP,NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
+ NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,NGLOB,DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_ETA,DIFF_NSPEC2D_XI)
use dyn_array
@@ -269,12 +273,10 @@
! crustal_model_variables
type crustal_model_variables
sequence
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: thlr
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocp
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: velocs
- double precision, dimension(NKEYS_CRUST,NLAYERS_CRUST) :: dens
- character(len=2) abbreviation(NCAP_CRUST/2,NCAP_CRUST)
- character(len=2) code(NKEYS_CRUST)
+ real(kind=4) velocp(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) velocs(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) dens(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
+ real(kind=4) thlr(0:2*NLON_CRUST,0:2*NLAT_CRUST,3:7)
end type crustal_model_variables
type (crustal_model_variables) CM_V
@@ -425,9 +427,9 @@
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: rmins,rmaxs
! arrays for BCAST
- integer, dimension(40) :: bcast_integer
- double precision, dimension(30) :: bcast_double_precision
- logical, dimension(26) :: bcast_logical
+ integer, dimension(NVALUES_bcast_integer) :: bcast_integer
+ double precision, dimension(NVALUES_bcast_double_precision) :: bcast_double_precision
+ logical, dimension(NVALUES_bcast_logical) :: bcast_logical
integer, parameter :: maxker=200
integer, parameter :: maxl=72
@@ -530,6 +532,10 @@
double precision, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT3D) :: Qmu_store
double precision, dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT3D) :: tau_e_store
+! receiver information
+ integer :: nrec,ios
+ character(len=150) :: STATIONS,rec_filename,dummystring
+
!!!!! DK DK for merged version, all the arrays below are allocated statically instead
#ifdef USE_MPI
@@ -647,6 +653,12 @@
! communication pattern for corners between chunks
integer, dimension(NCORNERSCHUNKS_VAL) :: iproc_master_corners,iproc_worker1_corners,iproc_worker2_corners
+
+ double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS,7) :: array_to_broadcast_1
+ integer, dimension(MAX_NUM_REGIONS,12) :: array_to_broadcast_2
+ integer, dimension(NB_SQUARE_EDGES_ONEDIR,NB_CUT_CASE,2) :: array_to_broadcast_3
+ double precision, dimension(2) :: array_to_broadcast_4
+ double precision, dimension(0:NK,0:NS,0:NS,4) :: array_to_broadcast_5
#endif
! get the base pathname for output files
@@ -678,7 +690,7 @@
write(IMAIN,*)
endif
- if (myrank==0) then
+ if (myrank == 0) then
! read the parameter file and compute additional parameters
call read_compute_parameters(MIN_ATTENUATION_PERIOD,MAX_ATTENUATION_PERIOD,NER_CRUST, &
NER_80_MOHO,NER_220_80,NER_400_220,NER_600_400,NER_670_600,NER_771_670, &
@@ -718,6 +730,20 @@
! count the total number of sources in the CMTSOLUTION file
call count_number_of_sources(NSOURCES)
+! read the number of receivers
+ rec_filename = 'DATA/STATIONS'
+ call get_value_string(STATIONS, 'solver.STATIONS', rec_filename)
+! get total number of receivers
+ if(myrank == 0) then
+ open(unit=IIN,file=STATIONS,iostat=ios,status='old',action='read')
+ nrec = 0
+ do while(ios == 0)
+ read(IIN,"(a)",iostat=ios) dummystring
+ if(ios == 0) nrec = nrec + 1
+ enddo
+ close(IIN)
+ endif
+
bcast_integer = (/MIN_ATTENUATION_PERIOD,MAX_ATTENUATION_PERIOD,NER_CRUST, &
NER_80_MOHO,NER_220_80,NER_400_220,NER_600_400,NER_670_600,NER_771_670, &
NER_TOPDDOUBLEPRIME_771,NER_CMB_TOPDDOUBLEPRIME,NER_OUTER_CORE, &
@@ -727,15 +753,17 @@
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,&
SIMULATION_TYPE,REFERENCE_1D_MODEL,THREE_D_MODEL,NPROC,NPROCTOT, &
NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_divide_central_cube,&
- MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP,ifirst_layer_aniso,ilast_layer_aniso/)
+ MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP,ifirst_layer_aniso,ilast_layer_aniso,nrec/)
bcast_logical = (/TRANSVERSE_ISOTROPY,ANISOTROPIC_3D_MANTLE,ANISOTROPIC_INNER_CORE, &
CRUSTAL,ELLIPTICITY,GRAVITY,ONE_CRUST,ROTATION,ISOTROPIC_3D_MANTLE, &
- TOPOGRAPHY,OCEANS,MOVIE_SURFACE,MOVIE_VOLUME,ATTENUATION_3D, &
+ TOPOGRAPHY,OCEANS,MOVIE_SURFACE,MOVIE_VOLUME,MOVIE_COARSE,ATTENUATION_3D, &
RECEIVERS_CAN_BE_BURIED,PRINT_SOURCE_TIME_FUNCTION, &
SAVE_MESH_FILES,ATTENUATION, &
- ABSORBING_CONDITIONS,INCLUDE_CENTRAL_CUBE,INFLATE_CENTRAL_CUBE,SAVE_FORWARD,CASE_3D,&
- CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,SAVE_ALL_SEISMOS_IN_ONE_FILE/)
+ ABSORBING_CONDITIONS,INCLUDE_CENTRAL_CUBE,INFLATE_CENTRAL_CUBE,SAVE_FORWARD,CASE_3D, &
+ OUTPUT_SEISMOS_ASCII_TEXT,OUTPUT_SEISMOS_SAC_ALPHANUM,OUTPUT_SEISMOS_SAC_BINARY, &
+ ROTATE_SEISMOGRAMS_RT,CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,&
+ WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE,USE_BINARY_FOR_LARGE_FILE/)
bcast_double_precision = (/DT,ANGULAR_WIDTH_XI_IN_DEGREES,ANGULAR_WIDTH_ETA_IN_DEGREES,CENTER_LONGITUDE_IN_DEGREES, &
CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH,ROCEAN,RMIDDLE_CRUST, &
@@ -747,46 +775,81 @@
! broadcast the information read on the master to the nodes
#ifdef USE_MPI
- call MPI_BCAST(bcast_integer,40,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(bcast_integer,NVALUES_bcast_integer,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(bcast_double_precision,30,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(bcast_double_precision,NVALUES_bcast_double_precision,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(bcast_logical,26,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(bcast_logical,NVALUES_bcast_logical,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(MODEL,150,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ner,MAX_NUMBER_OF_MESH_LAYERS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ratio_sampling_array,MAX_NUMBER_OF_MESH_LAYERS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(doubling_index,MAX_NUMBER_OF_MESH_LAYERS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(this_layer_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(r_bottom,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(r_top,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(rmins,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(rmaxs,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(rmaxs,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+! integer
+ array_to_broadcast_1(:,1) = ner(:)
+ array_to_broadcast_1(:,2) = ratio_sampling_array(:)
+ array_to_broadcast_1(:,3) = doubling_index(:)
- call MPI_BCAST(this_layer_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
+! double precision
+ array_to_broadcast_1(:,4) = r_bottom(:)
+ array_to_broadcast_1(:,5) = r_top(:)
+ array_to_broadcast_1(:,6) = rmins(:)
+ array_to_broadcast_1(:,7) = rmaxs(:)
- call MPI_BCAST(NSPEC,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_XI,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_ETA,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2DMAX_XMIN_XMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2DMAX_YMIN_YMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_BOTTOM,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_TOP,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC1D_RADIAL,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB1D_RADIAL,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB2DMAX_XMIN_XMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB2DMAX_YMIN_YMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(array_to_broadcast_1,7*MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+! integer
+ ner(:) = nint(array_to_broadcast_1(:,1))
+ ratio_sampling_array(:) = nint(array_to_broadcast_1(:,2))
+ doubling_index(:) = nint(array_to_broadcast_1(:,3))
+
+! double precision
+ r_bottom(:) = array_to_broadcast_1(:,4)
+ r_top(:) = array_to_broadcast_1(:,5)
+ rmins(:) = array_to_broadcast_1(:,6)
+ rmaxs(:) = array_to_broadcast_1(:,7)
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ array_to_broadcast_2(:, 1) = NSPEC(:)
+ array_to_broadcast_2(:, 2) = NSPEC2D_XI(:)
+ array_to_broadcast_2(:, 3) = NSPEC2D_ETA(:)
+ array_to_broadcast_2(:, 4) = NSPEC2DMAX_XMIN_XMAX(:)
+ array_to_broadcast_2(:, 5) = NSPEC2DMAX_YMIN_YMAX(:)
+ array_to_broadcast_2(:, 6) = NSPEC2D_BOTTOM(:)
+ array_to_broadcast_2(:, 7) = NSPEC2D_TOP(:)
+ array_to_broadcast_2(:, 8) = NSPEC1D_RADIAL(:)
+ array_to_broadcast_2(:, 9) = NGLOB1D_RADIAL(:)
+ array_to_broadcast_2(:,10) = NGLOB2DMAX_XMIN_XMAX(:)
+ array_to_broadcast_2(:,11) = NGLOB2DMAX_YMIN_YMAX(:)
+ array_to_broadcast_2(:,12) = NGLOB(:)
+
+ call MPI_BCAST(array_to_broadcast_2,12*MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+
+ NSPEC(:) = array_to_broadcast_2(:, 1)
+ NSPEC2D_XI(:) = array_to_broadcast_2(:, 2)
+ NSPEC2D_ETA(:) = array_to_broadcast_2(:, 3)
+ NSPEC2DMAX_XMIN_XMAX(:) = array_to_broadcast_2(:, 4)
+ NSPEC2DMAX_YMIN_YMAX(:) = array_to_broadcast_2(:, 5)
+ NSPEC2D_BOTTOM(:) = array_to_broadcast_2(:, 6)
+ NSPEC2D_TOP(:) = array_to_broadcast_2(:, 7)
+ NSPEC1D_RADIAL(:) = array_to_broadcast_2(:, 8)
+ NGLOB1D_RADIAL(:) = array_to_broadcast_2(:, 9)
+ NGLOB2DMAX_XMIN_XMAX(:) = array_to_broadcast_2(:,10)
+ NGLOB2DMAX_YMIN_YMAX(:) = array_to_broadcast_2(:,11)
+ NGLOB(:) = array_to_broadcast_2(:,12)
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
call MPI_BCAST(DIFF_NSPEC1D_RADIAL,NB_SQUARE_CORNERS*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(DIFF_NSPEC2D_ETA,NB_SQUARE_EDGES_ONEDIR*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(DIFF_NSPEC2D_XI,NB_SQUARE_EDGES_ONEDIR*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ array_to_broadcast_3(:,:,1) = DIFF_NSPEC2D_XI(:,:)
+ array_to_broadcast_3(:,:,2) = DIFF_NSPEC2D_ETA(:,:)
+ call MPI_BCAST(array_to_broadcast_3,2*NB_SQUARE_EDGES_ONEDIR*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ DIFF_NSPEC2D_XI(:,:) = array_to_broadcast_3(:,:,1)
+ DIFF_NSPEC2D_ETA(:,:) = array_to_broadcast_3(:,:,2)
#endif
- if (myrank /=0) then
-
MIN_ATTENUATION_PERIOD = bcast_integer(1)
MAX_ATTENUATION_PERIOD = bcast_integer(2)
NER_CRUST = bcast_integer(3)
@@ -827,6 +890,7 @@
MOVIE_STOP = bcast_integer(38)
ifirst_layer_aniso = bcast_integer(39)
ilast_layer_aniso = bcast_integer(40)
+ nrec = bcast_integer(41)
TRANSVERSE_ISOTROPY = bcast_logical(1)
ANISOTROPIC_3D_MANTLE = bcast_logical(2)
@@ -841,19 +905,26 @@
OCEANS = bcast_logical(11)
MOVIE_SURFACE = bcast_logical(12)
MOVIE_VOLUME = bcast_logical(13)
- ATTENUATION_3D = bcast_logical(14)
- RECEIVERS_CAN_BE_BURIED = bcast_logical(15)
- PRINT_SOURCE_TIME_FUNCTION = bcast_logical(16)
- SAVE_MESH_FILES = bcast_logical(17)
- ATTENUATION = bcast_logical(18)
- ABSORBING_CONDITIONS = bcast_logical(19)
- INCLUDE_CENTRAL_CUBE = bcast_logical(20)
- INFLATE_CENTRAL_CUBE = bcast_logical(21)
- SAVE_FORWARD = bcast_logical(22)
- CASE_3D = bcast_logical(23)
- CUT_SUPERBRICK_XI = bcast_logical(24)
- CUT_SUPERBRICK_ETA = bcast_logical(25)
- SAVE_ALL_SEISMOS_IN_ONE_FILE = bcast_logical(26)
+ MOVIE_COARSE = bcast_logical(14)
+ ATTENUATION_3D = bcast_logical(15)
+ RECEIVERS_CAN_BE_BURIED = bcast_logical(16)
+ PRINT_SOURCE_TIME_FUNCTION = bcast_logical(17)
+ SAVE_MESH_FILES = bcast_logical(18)
+ ATTENUATION = bcast_logical(19)
+ ABSORBING_CONDITIONS = bcast_logical(20)
+ INCLUDE_CENTRAL_CUBE = bcast_logical(21)
+ INFLATE_CENTRAL_CUBE = bcast_logical(22)
+ SAVE_FORWARD = bcast_logical(23)
+ CASE_3D = bcast_logical(24)
+ OUTPUT_SEISMOS_ASCII_TEXT = bcast_logical(25)
+ OUTPUT_SEISMOS_SAC_ALPHANUM = bcast_logical(26)
+ OUTPUT_SEISMOS_SAC_BINARY = bcast_logical(27)
+ ROTATE_SEISMOGRAMS_RT = bcast_logical(28)
+ CUT_SUPERBRICK_XI = bcast_logical(29)
+ CUT_SUPERBRICK_ETA = bcast_logical(30)
+ WRITE_SEISMOGRAMS_BY_MASTER = bcast_logical(31)
+ SAVE_ALL_SEISMOS_IN_ONE_FILE = bcast_logical(32)
+ USE_BINARY_FOR_LARGE_FILE = bcast_logical(33)
DT = bcast_double_precision(1)
ANGULAR_WIDTH_XI_IN_DEGREES = bcast_double_precision(2)
@@ -886,8 +957,6 @@
MOVIE_NORTH = bcast_double_precision(29)
MOVIE_SOUTH = bcast_double_precision(30)
- endif
-
! check that the code is running with the requested number of processes
if(sizeprocs /= NPROCTOT) call exit_MPI(myrank,'wrong number of MPI processes')
@@ -1088,10 +1157,16 @@
if(myrank == 0) call read_mantle_model(D3MM_V)
! broadcast the information read on the master to the nodes
#ifdef USE_MPI
- call MPI_BCAST(D3MM_V%dvs_a,(NK+1)*(NS+1)*(NS+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(D3MM_V%dvs_b,(NK+1)*(NS+1)*(NS+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(D3MM_V%dvp_a,(NK+1)*(NS+1)*(NS+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(D3MM_V%dvp_b,(NK+1)*(NS+1)*(NS+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+ array_to_broadcast_5(:,:,:,1) = D3MM_V%dvs_a
+ array_to_broadcast_5(:,:,:,2) = D3MM_V%dvs_b
+ array_to_broadcast_5(:,:,:,3) = D3MM_V%dvp_a
+ array_to_broadcast_5(:,:,:,4) = D3MM_V%dvp_b
+ call MPI_BCAST(array_to_broadcast_5,4*(NK+1)*(NS+1)*(NS+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
+ D3MM_V%dvs_a = array_to_broadcast_5(:,:,:,1)
+ D3MM_V%dvs_b = array_to_broadcast_5(:,:,:,2)
+ D3MM_V%dvp_a = array_to_broadcast_5(:,:,:,3)
+ D3MM_V%dvp_b = array_to_broadcast_5(:,:,:,4)
+
call MPI_BCAST(D3MM_V%spknt,NK+1,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(D3MM_V%qq0,(NK+1)*(NK+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(D3MM_V%qq,3*(NK+1)*(NK+1),MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
@@ -1298,15 +1373,13 @@
if(CRUSTAL) then
! the variables read are declared and stored in structure CM_V
- if(myrank == 0) call read_crustal_model(CM_V)
+ if(myrank == 0) call read_smoothed_3D_crustal_model(CM_V)
! broadcast the information read on the master to the nodes
#ifdef USE_MPI
- call MPI_BCAST(CM_V%thlr,NKEYS_CRUST*NLAYERS_CRUST,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(CM_V%velocp,NKEYS_CRUST*NLAYERS_CRUST,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(CM_V%velocs,NKEYS_CRUST*NLAYERS_CRUST,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(CM_V%dens,NKEYS_CRUST*NLAYERS_CRUST,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(CM_V%abbreviation,NCAP_CRUST*NCAP_CRUST,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(CM_V%code,2*NKEYS_CRUST,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(CM_V%velocp,(2*NLON_CRUST+1)*(2*NLAT_CRUST+1)*5,MPI_REAL,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(CM_V%velocs,(2*NLON_CRUST+1)*(2*NLAT_CRUST+1)*5,MPI_REAL,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(CM_V%dens,(2*NLON_CRUST+1)*(2*NLAT_CRUST+1)*5,MPI_REAL,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(CM_V%thlr,(2*NLON_CRUST+1)*(2*NLAT_CRUST+1)*5,MPI_REAL,0,MPI_COMM_WORLD,ier)
#endif
endif
@@ -1318,8 +1391,12 @@
if(ATTENUATION .and. ATTENUATION_3D) then
if(myrank == 0) call read_attenuation_model(MIN_ATTENUATION_PERIOD, MAX_ATTENUATION_PERIOD, AM_V)
#ifdef USE_MPI
- call MPI_BCAST(AM_V%min_period, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
- call MPI_BCAST(AM_V%max_period, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
+ array_to_broadcast_4(1) = AM_V%min_period
+ array_to_broadcast_4(2) = AM_V%max_period
+ call MPI_BCAST(array_to_broadcast_4,2,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD, ier)
+ AM_V%min_period = array_to_broadcast_4(1)
+ AM_V%max_period = array_to_broadcast_4(2)
+
call MPI_BCAST(AM_V%QT_c_source, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
call MPI_BCAST(AM_V%Qtau_s(1), 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
call MPI_BCAST(AM_V%Qtau_s(2), 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
@@ -1330,8 +1407,11 @@
if(ATTENUATION .and. .not. ATTENUATION_3D) then
if(myrank == 0) call read_attenuation_model(MIN_ATTENUATION_PERIOD, MAX_ATTENUATION_PERIOD, AM_V)
#ifdef USE_MPI
- call MPI_BCAST(AM_V%min_period, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
- call MPI_BCAST(AM_V%max_period, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ier)
+ array_to_broadcast_4(1) = AM_V%min_period
+ array_to_broadcast_4(2) = AM_V%max_period
+ call MPI_BCAST(array_to_broadcast_4,2,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD, ier)
+ AM_V%min_period = array_to_broadcast_4(1)
+ AM_V%max_period = array_to_broadcast_4(2)
#endif
call attenuation_model_setup(REFERENCE_1D_MODEL, RICB, RCMB, R670, R220, R80,AM_V,M1066a_V,Mak135_V,Mref_V,SEA1DM_V,AM_S,AS_V)
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -1150,7 +1150,7 @@
ROTATE_SEISMOGRAMS_RT = .false.
WRITE_SEISMOGRAMS_BY_MASTER = .true.
SAVE_ALL_SEISMOS_IN_ONE_FILE = .true.
- USE_BINARY_FOR_LARGE_FILE = .false.
+ USE_BINARY_FOR_LARGE_FILE = .true.
RECEIVERS_CAN_BE_BURIED = .false.
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.F90 2008-09-11 07:25:49 UTC (rev 12862)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.F90 2008-09-11 15:03:34 UTC (rev 12863)
@@ -39,7 +39,7 @@
kappahstore_crust_mantle,muhstore_crust_mantle,eta_anisostore_crust_mantle,kappavstore_inner_core,muvstore_inner_core, &
rmass_crust_mantle,rmass_outer_core,rmass_inner_core,rmass_ocean_load, &
#ifdef USE_MPI
- NDIM_smaller_buffers,npoin2D_max_all,nrec,addressing,ibathy_topo, &
+ npoin2D_max_all,nrec,addressing,ibathy_topo, &
ibelm_xmin_inner_core,ibelm_xmax_inner_core,ibelm_ymin_inner_core,ibelm_ymax_inner_core,ibelm_bottom_inner_core, &
iboolleft_xi_crust_mantle,iboolright_xi_crust_mantle, iboolleft_eta_crust_mantle,iboolright_eta_crust_mantle, &
iboolleft_xi_outer_core,iboolright_xi_outer_core,iboolleft_eta_outer_core,iboolright_eta_outer_core, &
@@ -55,7 +55,11 @@
#endif
AM_V,xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle,&
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle, &
- gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,displ_crust_mantle)
+ gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,displ_crust_mantle, &
+ bcast_integer,bcast_double_precision,bcast_logical,MODEL,ner,ratio_sampling_array,doubling_index, &
+ r_bottom,r_top,rmins,rmaxs,this_layer_has_a_doubling,NSPEC_computed,NSPEC2D_XI,NSPEC2D_ETA, &
+ NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX,NSPEC2D_BOTTOM,NSPEC2D_TOP,NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
+ NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,NGLOB_computed,DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_ETA,DIFF_NSPEC2D_XI)
use dyn_array
@@ -197,8 +201,8 @@
! always three times bigger and therefore scalars can use the first part
! of the vector buffer in memory even if it has an additional index here
! allocate these automatic arrays in the memory stack to avoid memory fragmentation with "allocate()"
- integer :: npoin2D_max_all,NDIM_smaller_buffers
- real(kind=CUSTOM_REAL), dimension(NDIM_smaller_buffers,npoin2D_max_all) :: buffer_send_faces,buffer_received_faces
+ integer :: npoin2D_max_all
+ real(kind=CUSTOM_REAL), dimension(NDIM,npoin2D_max_all) :: buffer_send_faces,buffer_received_faces
#endif
@@ -528,10 +532,10 @@
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: rmins,rmaxs
logical :: CASE_3D
-! arrays for BCAST
- integer, dimension(40) :: bcast_integer
- double precision, dimension(30) :: bcast_double_precision
- logical, dimension(33) :: bcast_logical
+! arrays from BCAST
+ integer, dimension(NVALUES_bcast_integer) :: bcast_integer
+ double precision, dimension(NVALUES_bcast_double_precision) :: bcast_double_precision
+ logical, dimension(NVALUES_bcast_logical) :: bcast_logical
logical :: CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA
integer, dimension(NB_SQUARE_CORNERS,NB_CUT_CASE) :: DIFF_NSPEC1D_RADIAL
@@ -630,76 +634,8 @@
if(err_occurred() /= 0) call exit_MPI(myrank,'an error occurred while reading the parameter file')
- bcast_integer = (/MIN_ATTENUATION_PERIOD,MAX_ATTENUATION_PERIOD,NER_CRUST, &
- NER_80_MOHO,NER_220_80,NER_400_220,NER_600_400,NER_670_600,NER_771_670, &
- NER_TOPDDOUBLEPRIME_771,NER_CMB_TOPDDOUBLEPRIME,NER_OUTER_CORE, &
- NER_TOP_CENTRAL_CUBE_ICB,NEX_XI,NEX_ETA,RMOHO_FICTITIOUS_IN_MESHER, &
- NPROC_XI,NPROC_ETA,NTSTEP_BETWEEN_OUTPUT_SEISMOS, &
- NTSTEP_BETWEEN_READ_ADJSRC,NSTEP,NSOURCES,NTSTEP_BETWEEN_FRAMES, &
- NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,&
- SIMULATION_TYPE,REFERENCE_1D_MODEL,THREE_D_MODEL,NPROC,NPROCTOT, &
- NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_divide_central_cube,&
- MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP,ifirst_layer_aniso,ilast_layer_aniso/)
-
- bcast_logical = (/TRANSVERSE_ISOTROPY,ANISOTROPIC_3D_MANTLE,ANISOTROPIC_INNER_CORE, &
- CRUSTAL,ELLIPTICITY,GRAVITY,ONE_CRUST,ROTATION,ISOTROPIC_3D_MANTLE, &
- TOPOGRAPHY,OCEANS,MOVIE_SURFACE,MOVIE_VOLUME,MOVIE_COARSE,ATTENUATION_3D, &
- RECEIVERS_CAN_BE_BURIED,PRINT_SOURCE_TIME_FUNCTION, &
- SAVE_MESH_FILES,ATTENUATION, &
- ABSORBING_CONDITIONS,INCLUDE_CENTRAL_CUBE,INFLATE_CENTRAL_CUBE,SAVE_FORWARD,CASE_3D, &
- OUTPUT_SEISMOS_ASCII_TEXT,OUTPUT_SEISMOS_SAC_ALPHANUM,OUTPUT_SEISMOS_SAC_BINARY, &
- ROTATE_SEISMOGRAMS_RT,CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,&
- WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE,USE_BINARY_FOR_LARGE_FILE/)
-
- bcast_double_precision = (/DT,ANGULAR_WIDTH_XI_IN_DEGREES,ANGULAR_WIDTH_ETA_IN_DEGREES,CENTER_LONGITUDE_IN_DEGREES, &
- CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH,ROCEAN,RMIDDLE_CRUST, &
- RMOHO,R80,R120,R220,R400,R600,R670,R771,RTOPDDOUBLEPRIME,RCMB,RICB, &
- R_CENTRAL_CUBE,RHO_TOP_OC,RHO_BOTTOM_OC,RHO_OCEANS,HDUR_MOVIE, &
- MOVIE_TOP,MOVIE_BOTTOM,MOVIE_WEST,MOVIE_EAST,MOVIE_NORTH,MOVIE_SOUTH/)
-
endif
-! broadcast the information read on the master to the nodes
-#ifdef USE_MPI
- call MPI_BCAST(bcast_integer,40,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(bcast_double_precision,30,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(bcast_logical,33,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(MODEL,150,MPI_CHARACTER,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(ner,MAX_NUMBER_OF_MESH_LAYERS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(ratio_sampling_array,MAX_NUMBER_OF_MESH_LAYERS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(doubling_index,MAX_NUMBER_OF_MESH_LAYERS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(r_bottom,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(r_top,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(rmins,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(rmaxs,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(this_layer_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(NSPEC_computed,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_XI,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_ETA,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2DMAX_XMIN_XMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2DMAX_YMIN_YMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_BOTTOM,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC2D_TOP,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NSPEC1D_RADIAL,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB1D_RADIAL,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB2DMAX_XMIN_XMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB2DMAX_YMIN_YMAX,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(NGLOB_computed,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
-
- call MPI_BCAST(DIFF_NSPEC1D_RADIAL,NB_SQUARE_CORNERS*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(DIFF_NSPEC2D_ETA,NB_SQUARE_EDGES_ONEDIR*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(DIFF_NSPEC2D_XI,NB_SQUARE_EDGES_ONEDIR*NB_CUT_CASE,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
-#endif
-
- if (myrank /=0) then
-
MIN_ATTENUATION_PERIOD = bcast_integer(1)
MAX_ATTENUATION_PERIOD = bcast_integer(2)
NER_CRUST = bcast_integer(3)
@@ -740,6 +676,10 @@
MOVIE_STOP = bcast_integer(38)
ifirst_layer_aniso = bcast_integer(39)
ilast_layer_aniso = bcast_integer(40)
+ ilast_layer_aniso = bcast_integer(40)
+#ifdef USE_MPI
+ nrec = bcast_integer(41)
+#endif
TRANSVERSE_ISOTROPY = bcast_logical(1)
ANISOTROPIC_3D_MANTLE = bcast_logical(2)
@@ -806,8 +746,6 @@
MOVIE_NORTH = bcast_double_precision(29)
MOVIE_SOUTH = bcast_double_precision(30)
- endif
-
! check simulation pararmeters
if (SIMULATION_TYPE /= 1 .and. SIMULATION_TYPE /= 2 .and. SIMULATION_TYPE /= 3) &
call exit_MPI(myrank, 'SIMULATION_TYPE could be only 1, 2, or 3')
@@ -2357,7 +2295,7 @@
buffer_send_chunkcorners_vector,buffer_recv_chunkcorners_vector, &
NUMMSGS_FACES,NUM_MSG_TYPES,NCORNERSCHUNKS, &
NPROC_XI,NPROC_ETA,NGLOB1D_RADIAL(IREGION_CRUST_MANTLE), &
- NGLOB1D_RADIAL(IREGION_INNER_CORE),NCHUNKS,NDIM_smaller_buffers)
+ NGLOB1D_RADIAL(IREGION_INNER_CORE),NCHUNKS)
#endif
!---
@@ -2515,7 +2453,7 @@
!! DK DK added this for Gordon Bell runs, to save seismograms every 10000 time steps
!! DK DK just in case the full simulation does not finish
- if(PATCH_FOR_GORDON_BELL .and. mod(it,10000) == 0) then
+ if(SAVE_PARTIAL_SEISMOGRAMS .and. mod(it,10000) == 0) then
call write_seismograms(myrank,uxdstore,uydstore,uzdstore,number_receiver_global,station_name, &
network_name,stlat,stlon,stele,nrec,nrec_local,DT,t0,it_end, &
yr_SAC,jda_SAC,ho_SAC,mi_SAC,sec_SAC,t_cmt_SAC, &
More information about the cig-commits
mailing list