[cig-commits] r18139 - in seismo/3D/SPECFEM3D_GLOBE/trunk: EXAMPLES/noise_examples src/specfem3D

yangl at geodynamics.org yangl at geodynamics.org
Fri Mar 25 07:28:40 PDT 2011


Author: yangl
Date: 2011-03-25 07:28:40 -0700 (Fri, 25 Mar 2011)
New Revision: 18139

Modified:
   seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/noise_examples/pre-processing
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90
   seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
Log:
NOISE GLOBAL bugs fixed

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/noise_examples/pre-processing
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/noise_examples/pre-processing	2011-03-24 22:25:01 UTC (rev 18138)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/noise_examples/pre-processing	2011-03-25 14:28:40 UTC (rev 18139)
@@ -9,6 +9,7 @@
 
 cd ../../
 
+rm -rf   SEM zzz_job_info NOISE_TOMOGRAPHY ZZZ_*
 mkdir -p SEM zzz_job_info NOISE_TOMOGRAPHY
 
 cp $dir/adj_traveltime_filter.f90        ./

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90	2011-03-24 22:25:01 UTC (rev 18138)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90	2011-03-25 14:28:40 UTC (rev 18139)
@@ -215,14 +215,14 @@
   subroutine check_parameters_noise(myrank,NOISE_TOMOGRAPHY,SIMULATION_TYPE,SAVE_FORWARD, &
                                     NUMBER_OF_RUNS, NUMBER_OF_THIS_RUN,ROTATE_SEISMOGRAMS_RT, &
                                     SAVE_ALL_SEISMOS_IN_ONE_FILE, USE_BINARY_FOR_LARGE_FILE, &
-                                    MOVIE_COARSE,LOCAL_PATH,NSPEC_TOP)
+                                    MOVIE_COARSE,LOCAL_PATH,NSPEC_TOP,NSTEP)
   implicit none
   include 'mpif.h'
   include "precision.h"
   include "constants.h"
   include "OUTPUT_FILES/values_from_mesher.h"
   ! input parameters
-  integer :: myrank,NOISE_TOMOGRAPHY,SIMULATION_TYPE,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN, NSPEC_TOP
+  integer :: myrank,NOISE_TOMOGRAPHY,SIMULATION_TYPE,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NSPEC_TOP,NSTEP
   logical :: SAVE_FORWARD,ROTATE_SEISMOGRAMS_RT,SAVE_ALL_SEISMOS_IN_ONE_FILE, USE_BINARY_FOR_LARGE_FILE,MOVIE_COARSE
   character(len=150) :: LOCAL_PATH
   ! output parameters
@@ -278,11 +278,11 @@
 
   if (NOISE_TOMOGRAPHY/=0) then
      ! save/read the surface movie using the same c routine as we do for absorbing boundaries (file ID is 9)
-     reclen=CUSTOM_REAL*NDIM*NGLLX*NGLLY*NSPEC_TOP
+     reclen=CUSTOM_REAL*NDIM*NGLLX*NGLLY*NSPEC_TOP*NSTEP
      write(outputname,"('/proc',i6.6,'_surface_movie')") myrank
-     if (NOISE_TOMOGRAPHY==1) call open_file_abs_w(9,trim(LOCAL_PATH)//outputname,len_trim(trim(LOCAL_PATH)//outputname),reclen)
-     if (NOISE_TOMOGRAPHY==2) call open_file_abs_r(9,trim(LOCAL_PATH)//outputname,len_trim(trim(LOCAL_PATH)//outputname),reclen)
-     if (NOISE_TOMOGRAPHY==3) call open_file_abs_r(9,trim(LOCAL_PATH)//outputname,len_trim(trim(LOCAL_PATH)//outputname),reclen)
+     if (NOISE_TOMOGRAPHY==1) call open_file_abs_w(9,trim(LOCAL_PATH)//trim(outputname),len_trim(trim(LOCAL_PATH)//trim(outputname)),reclen)
+     if (NOISE_TOMOGRAPHY==2) call open_file_abs_r(9,trim(LOCAL_PATH)//trim(outputname),len_trim(trim(LOCAL_PATH)//trim(outputname)),reclen)
+     if (NOISE_TOMOGRAPHY==3) call open_file_abs_r(9,trim(LOCAL_PATH)//trim(outputname),len_trim(trim(LOCAL_PATH)//trim(outputname)),reclen)
   endif
 
   end subroutine check_parameters_noise
@@ -425,69 +425,80 @@
 ! subroutine for NOISE TOMOGRAPHY
 ! step 1: calculate the "ensemble forward source"
 ! save surface movie (displacement) at every time steps, for step 2 & 3.
-  subroutine noise_save_surface_movie_original(myrank,nmovie_points,displ_crust_mantle, &
-                    xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle, &
-                    store_val_x,store_val_y,store_val_z, &
-                    store_val_ux,store_val_uy,store_val_uz, &
-                    ibelm_top_crust_mantle,ibool_crust_mantle,nspec_top, &
-                    NIT,it,LOCAL_PATH)
-  implicit none
-  include 'mpif.h'
-  include "precision.h"
-  include "constants.h"
-  include "OUTPUT_FILES/values_from_mesher.h"
-  ! input parameters
-  integer :: myrank,nmovie_points,nspec_top,NIT,it
-  integer, dimension(NSPEC2D_TOP_CM) :: ibelm_top_crust_mantle
-  integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
-  real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) ::  displ_crust_mantle
-  real(kind=CUSTOM_REAL), dimension(NGLOB_CRUST_MANTLE) :: &
-        xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle
-  character(len=150) :: LOCAL_PATH
-  ! output parameters
-  ! local parameters
-  integer :: ipoin,ispec2D,ispec,i,j,k,iglob
-  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: &
-      store_val_x,store_val_y,store_val_z, &
-      store_val_ux,store_val_uy,store_val_uz
-  character(len=150) :: outputname
+!
+! there are two subroutines --- noise_save_surface_movie_original & noise_save_surface_movie
+!    noise_save_surface_movie_original is implemented at first, which creates one file at each time step
+!    noise_save_surface_movie is implemented later, which utilizes 'src/shared/write_c_binary.c' for faster I/O,
+!                                                   which creates one file for the all time steps
+!
+! by this modification, the efficiency is greatly improved
+! and now, it should be OK to run NOISE_TOMOGRAPHY on a cluster with global storage
 
+!!!!! original implementation, not used anymore (but kept here for references)
+!  subroutine noise_save_surface_movie_original(myrank,nmovie_points,displ_crust_mantle, &
+!                    xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle, &
+!                    store_val_x,store_val_y,store_val_z, &
+!                    store_val_ux,store_val_uy,store_val_uz, &
+!                    ibelm_top_crust_mantle,ibool_crust_mantle,nspec_top, &
+!                    NIT,it,LOCAL_PATH)
+!  implicit none
+!  include 'mpif.h'
+!  include "precision.h"
+!  include "constants.h"
+!  include "OUTPUT_FILES/values_from_mesher.h"
+!  ! input parameters
+!  integer :: myrank,nmovie_points,nspec_top,NIT,it
+!  integer, dimension(NSPEC2D_TOP_CM) :: ibelm_top_crust_mantle
+!  integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
+!  real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) ::  displ_crust_mantle
+!  real(kind=CUSTOM_REAL), dimension(NGLOB_CRUST_MANTLE) :: &
+!        xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle
+!  character(len=150) :: LOCAL_PATH
+!  ! output parameters
+!  ! local parameters
+!  integer :: ipoin,ispec2D,ispec,i,j,k,iglob
+!  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: &
+!      store_val_x,store_val_y,store_val_z, &
+!      store_val_ux,store_val_uy,store_val_uz
+!  character(len=150) :: outputname
+!
+!
+!  ! get coordinates of surface mesh and surface displacement
+!  ipoin = 0
+!  do ispec2D = 1, nspec_top ! NSPEC2D_TOP(IREGION_CRUST_MANTLE)
+!    ispec = ibelm_top_crust_mantle(ispec2D)
+!
+!    k = NGLLZ
+!
+!    ! loop on all the points inside the element
+!    do j = 1,NGLLY,NIT
+!      do i = 1,NGLLX,NIT
+!        ipoin = ipoin + 1
+!        iglob = ibool_crust_mantle(i,j,k,ispec)
+!        store_val_x(ipoin) = xstore_crust_mantle(iglob)
+!        store_val_y(ipoin) = ystore_crust_mantle(iglob)
+!        store_val_z(ipoin) = zstore_crust_mantle(iglob)
+!        store_val_ux(ipoin) = displ_crust_mantle(1,iglob)
+!        store_val_uy(ipoin) = displ_crust_mantle(2,iglob)
+!        store_val_uz(ipoin) = displ_crust_mantle(3,iglob)
+!      enddo
+!    enddo
+!
+!  enddo
+!
+!  ! save surface motion to disk
+!  ! LOCAL storage is better than GLOBAL, because we have to save the 'movie' at every time step
+!  ! also note that the surface movie does NOT have to be shared with other nodes/CPUs
+!  ! change LOCAL_PATH specified in "DATA/Par_file"
+!    write(outputname,"('/proc',i6.6,'_surface_movie',i6.6)") myrank, it
+!    open(unit=IOUT_NOISE,file=trim(LOCAL_PATH)//outputname,status='unknown',form='unformatted',action='write')
+!    write(IOUT_NOISE) store_val_ux
+!    write(IOUT_NOISE) store_val_uy
+!    write(IOUT_NOISE) store_val_uz
+!    close(IOUT_NOISE)
+!  end subroutine noise_save_surface_movie_original
 
-  ! get coordinates of surface mesh and surface displacement
-  ipoin = 0
-  do ispec2D = 1, nspec_top ! NSPEC2D_TOP(IREGION_CRUST_MANTLE)
-    ispec = ibelm_top_crust_mantle(ispec2D)
-
-    k = NGLLZ
-
-    ! loop on all the points inside the element
-    do j = 1,NGLLY,NIT
-      do i = 1,NGLLX,NIT
-        ipoin = ipoin + 1
-        iglob = ibool_crust_mantle(i,j,k,ispec)
-        store_val_x(ipoin) = xstore_crust_mantle(iglob)
-        store_val_y(ipoin) = ystore_crust_mantle(iglob)
-        store_val_z(ipoin) = zstore_crust_mantle(iglob)
-        store_val_ux(ipoin) = displ_crust_mantle(1,iglob)
-        store_val_uy(ipoin) = displ_crust_mantle(2,iglob)
-        store_val_uz(ipoin) = displ_crust_mantle(3,iglob)
-      enddo
-    enddo
-
-  enddo
-
-  ! save surface motion to disk
-  ! LOCAL storage is better than GLOBAL, because we have to save the 'movie' at every time step
-  ! also note that the surface movie does NOT have to be shared with other nodes/CPUs
-  ! change LOCAL_PATH specified in "DATA/Par_file"
-    write(outputname,"('/proc',i6.6,'_surface_movie',i6.6)") myrank, it
-    open(unit=IOUT_NOISE,file=trim(LOCAL_PATH)//outputname,status='unknown',form='unformatted',action='write')
-    write(IOUT_NOISE) store_val_ux
-    write(IOUT_NOISE) store_val_uy
-    write(IOUT_NOISE) store_val_uz
-    close(IOUT_NOISE)
-  end subroutine noise_save_surface_movie_original
-
+!!!!! the improved version, with some dummy variables (i.e., variables that are NOT used in the this improved version)
   subroutine noise_save_surface_movie(myrank,nmovie_points,displ_crust_mantle, &
                     xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle, &
                     store_val_x,store_val_y,store_val_z, &
@@ -534,7 +545,6 @@
   ! LOCAL storage is better than GLOBAL, because we have to save the 'movie' at every time step
   ! also note that the surface movie does NOT have to be shared with other nodes/CPUs
   ! change LOCAL_PATH specified in "DATA/Par_file"
-  write(outputname,"('/proc',i6.6,'_surface_movie')") myrank
   call write_abs(9,SURFACE_MOVIE,CUSTOM_REAL*NDIM*NGLLX*NGLLY*nspec_top,it)
   deallocate(SURFACE_MOVIE)
   end subroutine noise_save_surface_movie
@@ -548,72 +558,83 @@
 ! read surface movie (displacement) at every time steps, injected as the source of "ensemble forward wavefield"
 ! in step 2, call noise_read_add_surface_movie(..., NSTEP-it+1 ,...)
 ! in step 3, call noise_read_add_surface_movie(..., it ,...)
-  subroutine noise_read_add_surface_movie_original(myrank,nmovie_points,accel_crust_mantle, &
-                  normal_x_noise,normal_y_noise,normal_z_noise,mask_noise, &
-                  store_val_ux,store_val_uy,store_val_uz, &
-                  ibelm_top_crust_mantle,ibool_crust_mantle,nspec_top, &
-                  NIT,it,LOCAL_PATH,jacobian2D_top_crust_mantle,wgllwgll_xy)
-  implicit none
-  include 'mpif.h'
-  include "precision.h"
-  include "constants.h"
-  include "OUTPUT_FILES/values_from_mesher.h"
-  ! input parameters
-  integer :: myrank,nmovie_points,nspec_top,NIT,it
-  integer, dimension(NSPEC2D_TOP_CM) :: ibelm_top_crust_mantle
-  integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
-  real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NSPEC2D_TOP_CM) :: jacobian2D_top_crust_mantle
-  real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) :: accel_crust_mantle ! both input and output
-  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: normal_x_noise,normal_y_noise,normal_z_noise, mask_noise
-  real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY) :: wgllwgll_xy
-  character(len=150) :: LOCAL_PATH
-  ! output parameters
-  ! local parameters
-  integer :: ipoin,ispec2D,ispec,i,j,k,iglob,ios
-  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: store_val_ux,store_val_uy,store_val_uz
-  real(kind=CUSTOM_REAL) :: eta
-  character(len=150) :: outputname
+!
+! there are two subroutines --- noise_read_add_surface_movie_original & noise_read_add_surface_movie
+!    noise_read_add_surface_movie_original is implemented at first, which creates one file at each time step
+!    noise_read_add_surface_movie is implemented later, which utilizes 'src/shared/write_c_binary.c' for faster I/O,
+!                                                   which creates one file for the all time steps
+!
+! by this modification, the efficiency is greatly improved
+! and now, it should be OK to run NOISE_TOMOGRAPHY on a cluster with global storage
 
+!!!!! original implementation, not used anymore (but kept here for references)
+!  subroutine noise_read_add_surface_movie_original(myrank,nmovie_points,accel_crust_mantle, &
+!                  normal_x_noise,normal_y_noise,normal_z_noise,mask_noise, &
+!                  store_val_ux,store_val_uy,store_val_uz, &
+!                  ibelm_top_crust_mantle,ibool_crust_mantle,nspec_top, &
+!                  NIT,it,LOCAL_PATH,jacobian2D_top_crust_mantle,wgllwgll_xy)
+!  implicit none
+!  include 'mpif.h'
+!  include "precision.h"
+!  include "constants.h"
+!  include "OUTPUT_FILES/values_from_mesher.h"
+!  ! input parameters
+!  integer :: myrank,nmovie_points,nspec_top,NIT,it
+!  integer, dimension(NSPEC2D_TOP_CM) :: ibelm_top_crust_mantle
+!  integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
+!  real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NSPEC2D_TOP_CM) :: jacobian2D_top_crust_mantle
+!  real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) :: accel_crust_mantle ! both input and output
+!  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: normal_x_noise,normal_y_noise,normal_z_noise, mask_noise
+!  real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY) :: wgllwgll_xy
+!  character(len=150) :: LOCAL_PATH
+!  ! output parameters
+!  ! local parameters
+!  integer :: ipoin,ispec2D,ispec,i,j,k,iglob,ios
+!  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: store_val_ux,store_val_uy,store_val_uz
+!  real(kind=CUSTOM_REAL) :: eta
+!  character(len=150) :: outputname
+!
+!
+!  ! read surface movie
+!  write(outputname,"('/proc',i6.6,'_surface_movie',i6.6)") myrank, it
+!  open(unit=IIN_NOISE,file=trim(LOCAL_PATH)//outputname,status='old',form='unformatted',action='read',iostat=ios)
+!  if( ios /= 0)  call exit_MPI(myrank,'file '//trim(outputname)//' does NOT exist!')
+!  read(IIN_NOISE) store_val_ux
+!  read(IIN_NOISE) store_val_uy
+!  read(IIN_NOISE) store_val_uz
+!  close(IIN_NOISE)
+!
+!  ! get coordinates of surface mesh and surface displacement
+!  ipoin = 0
+!  do ispec2D = 1, nspec_top ! NSPEC2D_TOP(IREGION_CRUST_MANTLE)
+!    ispec = ibelm_top_crust_mantle(ispec2D)
+!
+!    k = NGLLZ
+!
+!    ! loop on all the points inside the element
+!    do j = 1,NGLLY,NIT
+!      do i = 1,NGLLX,NIT
+!        ipoin = ipoin + 1
+!        iglob = ibool_crust_mantle(i,j,k,ispec)
+!
+!        eta = store_val_ux(ipoin) * normal_x_noise(ipoin) + &
+!              store_val_uy(ipoin) * normal_y_noise(ipoin) + &
+!              store_val_uz(ipoin) * normal_z_noise(ipoin)
+!
+!        accel_crust_mantle(1,iglob) = accel_crust_mantle(1,iglob) + eta * mask_noise(ipoin) * normal_x_noise(ipoin) &
+!                                                      * wgllwgll_xy(i,j) * jacobian2D_top_crust_mantle(i,j,ispec2D)
+!        accel_crust_mantle(2,iglob) = accel_crust_mantle(2,iglob) + eta * mask_noise(ipoin) * normal_y_noise(ipoin) &
+!                                                      * wgllwgll_xy(i,j) * jacobian2D_top_crust_mantle(i,j,ispec2D)
+!        accel_crust_mantle(3,iglob) = accel_crust_mantle(3,iglob) + eta * mask_noise(ipoin) * normal_z_noise(ipoin) &
+!                                                      * wgllwgll_xy(i,j) * jacobian2D_top_crust_mantle(i,j,ispec2D)
+!      enddo
+!    enddo
+!
+!  enddo
+!
+!  end subroutine noise_read_add_surface_movie_original
 
-  ! read surface movie
-  write(outputname,"('/proc',i6.6,'_surface_movie',i6.6)") myrank, it
-  open(unit=IIN_NOISE,file=trim(LOCAL_PATH)//outputname,status='old',form='unformatted',action='read',iostat=ios)
-  if( ios /= 0)  call exit_MPI(myrank,'file '//trim(outputname)//' does NOT exist!')
-  read(IIN_NOISE) store_val_ux
-  read(IIN_NOISE) store_val_uy
-  read(IIN_NOISE) store_val_uz
-  close(IIN_NOISE)
-
-  ! get coordinates of surface mesh and surface displacement
-  ipoin = 0
-  do ispec2D = 1, nspec_top ! NSPEC2D_TOP(IREGION_CRUST_MANTLE)
-    ispec = ibelm_top_crust_mantle(ispec2D)
-
-    k = NGLLZ
-
-    ! loop on all the points inside the element
-    do j = 1,NGLLY,NIT
-      do i = 1,NGLLX,NIT
-        ipoin = ipoin + 1
-        iglob = ibool_crust_mantle(i,j,k,ispec)
-
-        eta = store_val_ux(ipoin) * normal_x_noise(ipoin) + &
-              store_val_uy(ipoin) * normal_y_noise(ipoin) + &
-              store_val_uz(ipoin) * normal_z_noise(ipoin)
-
-        accel_crust_mantle(1,iglob) = accel_crust_mantle(1,iglob) + eta * mask_noise(ipoin) * normal_x_noise(ipoin) &
-                                                      * wgllwgll_xy(i,j) * jacobian2D_top_crust_mantle(i,j,ispec2D)
-        accel_crust_mantle(2,iglob) = accel_crust_mantle(2,iglob) + eta * mask_noise(ipoin) * normal_y_noise(ipoin) &
-                                                      * wgllwgll_xy(i,j) * jacobian2D_top_crust_mantle(i,j,ispec2D)
-        accel_crust_mantle(3,iglob) = accel_crust_mantle(3,iglob) + eta * mask_noise(ipoin) * normal_z_noise(ipoin) &
-                                                      * wgllwgll_xy(i,j) * jacobian2D_top_crust_mantle(i,j,ispec2D)
-      enddo
-    enddo
-
-  enddo
-
-  end subroutine noise_read_add_surface_movie_original
-
+!!!!! the improved version, with some dummy variables (i.e., variables that are NOT used in the this improved version)
   subroutine noise_read_add_surface_movie(myrank,nmovie_points,accel_crust_mantle, &
                   normal_x_noise,normal_y_noise,normal_z_noise,mask_noise, &
                   store_val_ux,store_val_uy,store_val_uz, &
@@ -644,7 +665,6 @@
 
   allocate(SURFACE_MOVIE(NDIM,NGLLX,NGLLY,nspec_top))
   ! read surface movie
-  write(outputname,"('/proc',i6.6,'_surface_movie')") myrank
   call read_abs(9,SURFACE_MOVIE,CUSTOM_REAL*NDIM*NGLLX*NGLLY*nspec_top,it)
   ! get coordinates of surface mesh and surface displacement
   ipoin = 0
@@ -683,71 +703,83 @@
 
 ! subroutine for NOISE TOMOGRAPHY
 ! step 3: constructing noise source strength kernel
-  subroutine compute_kernels_strength_noise_original(myrank,ibool_crust_mantle, &
-                          Sigma_kl_crust_mantle,displ_crust_mantle,deltat,it, &
-                          nmovie_points,normal_x_noise,normal_y_noise,normal_z_noise, &
-                          nspec_top,ibelm_top_crust_mantle,LOCAL_PATH)
-  implicit none
-  include "constants.h"
-  include "OUTPUT_FILES/values_from_mesher.h"
-  ! input parameters
-  integer :: myrank,nmovie_points,it,nspec_top
-  integer, dimension(NSPEC2D_TOP_CM) :: ibelm_top_crust_mantle
-  integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
-  real(kind=CUSTOM_REAL) :: deltat
-  real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE_ADJOINT) :: displ_crust_mantle
-  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: normal_x_noise,normal_y_noise,normal_z_noise
-  character(len=150) :: LOCAL_PATH
-  ! output parameters
-  real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT) :: &
-    Sigma_kl_crust_mantle
-  ! local parameters
-  integer :: i,j,k,ispec,iglob,ipoin,ispec2D,ios
-  real(kind=CUSTOM_REAL) :: eta
-  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: store_val_ux,store_val_uy,store_val_uz
-  character(len=150) :: outputname
+!
+! there are two subroutines --- compute_kernels_strength_noise_original & compute_kernels_strength_noise
+!    compute_kernels_strength_noise_original is implemented at first, which creates one file at each time step
+!    compute_kernels_strength_noise is implemented later, which utilizes 'src/shared/write_c_binary.c' for faster I/O,
+!                                                         which creates only one file for the all time steps
+!
+! by this modification, the efficiency is greatly improved
+! and now, it should be OK to run NOISE_TOMOGRAPHY on a cluster with global storage
 
+!!!!! original implementation, not used anymore (but kept here for references)
+!  subroutine compute_kernels_strength_noise_original(myrank,ibool_crust_mantle, &
+!                          Sigma_kl_crust_mantle,displ_crust_mantle,deltat,it, &
+!                          nmovie_points,normal_x_noise,normal_y_noise,normal_z_noise, &
+!                          nspec_top,ibelm_top_crust_mantle,LOCAL_PATH)
+!  implicit none
+!  include "constants.h"
+!  include "OUTPUT_FILES/values_from_mesher.h"
+!  ! input parameters
+!  integer :: myrank,nmovie_points,it,nspec_top
+!  integer, dimension(NSPEC2D_TOP_CM) :: ibelm_top_crust_mantle
+!  integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
+!  real(kind=CUSTOM_REAL) :: deltat
+!  real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE_ADJOINT) :: displ_crust_mantle
+!  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: normal_x_noise,normal_y_noise,normal_z_noise
+!  character(len=150) :: LOCAL_PATH
+!  ! output parameters
+!  real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT) :: &
+!    Sigma_kl_crust_mantle
+!  ! local parameters
+!  integer :: i,j,k,ispec,iglob,ipoin,ispec2D,ios
+!  real(kind=CUSTOM_REAL) :: eta
+!  real(kind=CUSTOM_REAL), dimension(nmovie_points) :: store_val_ux,store_val_uy,store_val_uz
+!  character(len=150) :: outputname
+!
+!
+!  ! read surface movie, needed for Sigma_kl_crust_mantle
+!  write(outputname,"('/proc',i6.6,'_surface_movie',i6.6)") myrank, it
+!  open(unit=IIN_NOISE,file=trim(LOCAL_PATH)//outputname,status='old',form='unformatted',action='read',iostat=ios)
+!  if( ios /= 0)  call exit_MPI(myrank,'file '//trim(outputname)//' does NOT exist!')
+!  call read_abs(9,SURFACE_MOVIE,CUSTOM_REAL*NDIM*NGLLX*NGLLY*nspec_top,it)
+!
+!  read(IIN_NOISE) store_val_ux
+!  read(IIN_NOISE) store_val_uy
+!  read(IIN_NOISE) store_val_uz
+!  close(IIN_NOISE)
+!
+!  ! noise source strength kernel
+!  ! to keep similar structure to other kernels, the source strength kernel is saved as a volumetric kernel
+!  ! but only updated at the surface, because the noise is generated there
+!  ipoin = 0
+!  do ispec2D = 1, nspec_top
+!    ispec = ibelm_top_crust_mantle(ispec2D)
+!
+!    k = NGLLZ
+!
+!    ! loop on all the points inside the element
+!    do j = 1,NGLLY
+!      do i = 1,NGLLX
+!        ipoin = ipoin + 1
+!        iglob = ibool_crust_mantle(i,j,k,ispec)
+!
+!        eta = store_val_ux(ipoin) * normal_x_noise(ipoin) + &
+!              store_val_uy(ipoin) * normal_y_noise(ipoin) + &
+!              store_val_uz(ipoin) * normal_z_noise(ipoin)
+!
+!        Sigma_kl_crust_mantle(i,j,k,ispec) =  Sigma_kl_crust_mantle(i,j,k,ispec) &
+!           + deltat * eta * ( normal_x_noise(ipoin) * displ_crust_mantle(1,iglob) &
+!                            + normal_y_noise(ipoin) * displ_crust_mantle(2,iglob) &
+!                            + normal_z_noise(ipoin) * displ_crust_mantle(3,iglob) )
+!      enddo
+!    enddo
+!
+!  enddo
+!
+!  end subroutine compute_kernels_strength_noise_original
 
-  ! read surface movie, needed for Sigma_kl_crust_mantle
-  write(outputname,"('/proc',i6.6,'_surface_movie',i6.6)") myrank, it
-  open(unit=IIN_NOISE,file=trim(LOCAL_PATH)//outputname,status='old',form='unformatted',action='read',iostat=ios)
-  if( ios /= 0)  call exit_MPI(myrank,'file '//trim(outputname)//' does NOT exist!')
-
-  read(IIN_NOISE) store_val_ux
-  read(IIN_NOISE) store_val_uy
-  read(IIN_NOISE) store_val_uz
-  close(IIN_NOISE)
-
-  ! noise source strength kernel
-  ! to keep similar structure to other kernels, the source strength kernel is saved as a volumetric kernel
-  ! but only updated at the surface, because the noise is generated there
-  ipoin = 0
-  do ispec2D = 1, nspec_top
-    ispec = ibelm_top_crust_mantle(ispec2D)
-
-    k = NGLLZ
-
-    ! loop on all the points inside the element
-    do j = 1,NGLLY
-      do i = 1,NGLLX
-        ipoin = ipoin + 1
-        iglob = ibool_crust_mantle(i,j,k,ispec)
-
-        eta = store_val_ux(ipoin) * normal_x_noise(ipoin) + &
-              store_val_uy(ipoin) * normal_y_noise(ipoin) + &
-              store_val_uz(ipoin) * normal_z_noise(ipoin)
-
-        Sigma_kl_crust_mantle(i,j,k,ispec) =  Sigma_kl_crust_mantle(i,j,k,ispec) &
-           + deltat * eta * ( normal_x_noise(ipoin) * displ_crust_mantle(1,iglob) &
-                            + normal_y_noise(ipoin) * displ_crust_mantle(2,iglob) &
-                            + normal_z_noise(ipoin) * displ_crust_mantle(3,iglob) )
-      enddo
-    enddo
-
-  enddo
-
-  end subroutine compute_kernels_strength_noise_original
-
+!!!!! the improved version, with some dummy variables (i.e., variables that are NOT used in the this improved version)
   subroutine compute_kernels_strength_noise(myrank,ibool_crust_mantle, &
                           Sigma_kl_crust_mantle,displ_crust_mantle,deltat,it, &
                           nmovie_points,normal_x_noise,normal_y_noise,normal_z_noise, &
@@ -776,7 +808,6 @@
 
   allocate(SURFACE_MOVIE(NDIM,NGLLX,NGLLY,nspec_top))
   ! read surface movie, needed for Sigma_kl_crust_mantle
-  write(outputname,"('/proc',i6.6,'_surface_movie')") myrank
   call read_abs(9,SURFACE_MOVIE,CUSTOM_REAL*NDIM*NGLLX*NGLLY*nspec_top,it)
   ! noise source strength kernel
   ! to keep similar structure to other kernels, the source strength kernel is saved as a volumetric kernel

Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90	2011-03-24 22:25:01 UTC (rev 18138)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90	2011-03-25 14:28:40 UTC (rev 18139)
@@ -1792,11 +1792,10 @@
                                   xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle, &
                                   irec_master_noise,normal_x_noise,normal_y_noise,normal_z_noise,mask_noise)
 
-       if (myrank == 0) &
        call check_parameters_noise(myrank,NOISE_TOMOGRAPHY,SIMULATION_TYPE,SAVE_FORWARD, &
                                   NUMBER_OF_RUNS, NUMBER_OF_THIS_RUN,ROTATE_SEISMOGRAMS_RT, &
                                   SAVE_ALL_SEISMOS_IN_ONE_FILE, USE_BINARY_FOR_LARGE_FILE, &
-                                  MOVIE_COARSE,LOCAL_PATH,NSPEC2D_TOP(IREGION_CRUST_MANTLE))
+                                  MOVIE_COARSE,LOCAL_PATH,NSPEC2D_TOP(IREGION_CRUST_MANTLE),NSTEP)
     endif
 !>YANGL
 



More information about the CIG-COMMITS mailing list