[cig-commits] r22227 - in seismo/3D/SPECFEM3D_GLOBE/branches/undo_att: . src/specfem3D
xie.zhinan at geodynamics.org
xie.zhinan at geodynamics.org
Wed Jun 12 01:31:01 PDT 2013
Author: xie.zhinan
Date: 2013-06-12 01:31:00 -0700 (Wed, 12 Jun 2013)
New Revision: 22227
Added:
seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/part4_save_kernel.F90
Modified:
seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/declaration_part_for_backward_wavefield_simulation.f90
seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/compute_kernels.f90
seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/specfem3D.F90
Log:
restructuring the code
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/declaration_part_for_backward_wavefield_simulation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/declaration_part_for_backward_wavefield_simulation.f90 2013-06-11 23:37:32 UTC (rev 22226)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/declaration_part_for_backward_wavefield_simulation.f90 2013-06-12 08:31:00 UTC (rev 22227)
@@ -25,8 +25,6 @@
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_OUTER_CORE_ADJOINT) :: b_vector_displ_outer_core
- real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_OUTER_CORE_ADJOINT) :: b_vector_displ_outer_core !ZN
-
real(kind=CUSTOM_REAL), dimension(:,:,:), allocatable :: b_buffer_send_faces,b_buffer_received_faces
integer :: b_iphase,b_iphase_CC,b_icall
Added: seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/part4_save_kernel.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/part4_save_kernel.F90 (rev 0)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/part4_save_kernel.F90 2013-06-12 08:31:00 UTC (rev 22227)
@@ -0,0 +1,60 @@
+
+ ! synchronize all processes
+ call MPI_BARRIER(MPI_COMM_WORLD,ier)
+ if( ier /= 0 ) call exit_mpi(myrank,'error synchronize saving forward')
+
+ ! dump kernel arrays
+ if (SIMULATION_TYPE == 3) then
+
+ ! crust mantle
+ if (SAVE_REGULAR_KL) then
+ call save_regular_kernels_crust_mantle(myrank, &
+ npoints_slice, hxir_reg, hetar_reg, hgammar_reg, &
+ scale_t,scale_displ, &
+ cijkl_kl_crust_mantle,rho_kl_crust_mantle, &
+ alpha_kl_crust_mantle,beta_kl_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle, &
+ rhostore_crust_mantle,muvstore_crust_mantle, &
+ kappavstore_crust_mantle,ibool_crust_mantle, &
+ kappahstore_crust_mantle,muhstore_crust_mantle, &
+ eta_anisostore_crust_mantle,ispec_is_tiso_crust_mantle,LOCAL_PATH)
+ else
+ call save_kernels_crust_mantle(myrank,scale_t,scale_displ, &
+ cijkl_kl_crust_mantle,rho_kl_crust_mantle, &
+ alpha_kl_crust_mantle,beta_kl_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle, &
+ rhostore_crust_mantle,muvstore_crust_mantle, &
+ kappavstore_crust_mantle,ibool_crust_mantle, &
+ kappahstore_crust_mantle,muhstore_crust_mantle, &
+ eta_anisostore_crust_mantle,ispec_is_tiso_crust_mantle,LOCAL_PATH)
+ endif
+
+ ! noise strength kernel
+ if (NOISE_TOMOGRAPHY == 3) then
+ call save_kernels_strength_noise(myrank,LOCAL_PATH,Sigma_kl_crust_mantle)
+ endif
+
+ ! outer core
+ call save_kernels_outer_core(myrank,scale_t,scale_displ, &
+ rho_kl_outer_core,alpha_kl_outer_core, &
+ rhostore_outer_core,kappavstore_outer_core, &
+ deviatoric_outercore,nspec_beta_kl_outer_core,beta_kl_outer_core,LOCAL_PATH)
+
+ ! inner core
+ call save_kernels_inner_core(myrank,scale_t,scale_displ, &
+ rho_kl_inner_core,beta_kl_inner_core,alpha_kl_inner_core, &
+ rhostore_inner_core,muvstore_inner_core,kappavstore_inner_core,LOCAL_PATH)
+
+ ! boundary kernel
+ if (SAVE_BOUNDARY_MESH) then
+ call save_kernels_boundary_kl(myrank,scale_t,scale_displ, &
+ moho_kl,d400_kl,d670_kl,cmb_kl,icb_kl,LOCAL_PATH,HONOR_1D_SPHERICAL_MOHO)
+ endif
+
+ ! approximate hessian
+ if( APPROXIMATE_HESS_KL) then
+ stop "APPROXIMATE_HESS_KL is not support for undoing att" !ZN
+ call save_kernels_hessian(myrank,scale_t,scale_displ,hess_kl_crust_mantle,LOCAL_PATH)
+ endif
+ endif
+
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/compute_kernels.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/compute_kernels.f90 2013-06-11 23:37:32 UTC (rev 22226)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/compute_kernels.f90 2013-06-12 08:31:00 UTC (rev 22227)
@@ -1003,5 +1003,326 @@
end subroutine compute_kernels_hessian
+!=====================================================================
+ subroutine compute_stain_crust_mantle(displ_crust_mantle,hprime_xx,hprime_yy,hprime_zz,ibool,&
+ xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz, &
+ epsilondev,epsilon_trace_over_3)
+ implicit none
+
+ include "constants.h"
+
+! include values created by the mesher
+! done for performance only using static allocation to allow for loop unrolling
+ include "OUTPUT_FILES/values_from_mesher.h"
+
+! array with derivatives of Lagrange polynomials and precalculated products
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLX) :: hprime_xx
+ real(kind=CUSTOM_REAL), dimension(NGLLY,NGLLY) :: hprime_yy
+ real(kind=CUSTOM_REAL), dimension(NGLLZ,NGLLZ) :: hprime_zz
+
+
+! arrays with mesh parameters per slice
+ integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: &
+ xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz
+
+ real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilon_trace_over_3
+
+ real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilondev
+
+ real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) :: displ_crust_mantle
+
+!local parameters
+ integer ispec,iglob,ispec_strain
+ integer i,j,k,l
+
+ real(kind=CUSTOM_REAL) tempx1l,tempx2l,tempx3l
+ real(kind=CUSTOM_REAL) tempy1l,tempy2l,tempy3l
+ real(kind=CUSTOM_REAL) tempz1l,tempz2l,tempz3l
+
+ real(kind=CUSTOM_REAL) hp1,hp2,hp3
+
+ real(kind=CUSTOM_REAL) xixl,xiyl,xizl,etaxl,etayl,etazl,gammaxl,gammayl,gammazl,jacobianl
+
+ real(kind=CUSTOM_REAL) duxdxl,duxdyl,duxdzl,duydxl,duydyl,duydzl,duzdxl,duzdyl,duzdzl
+
+ real(kind=CUSTOM_REAL) duxdxl_plus_duydyl,duxdxl_plus_duzdzl,duydyl_plus_duzdzl
+ real(kind=CUSTOM_REAL) duxdyl_plus_duydxl,duzdxl_plus_duxdzl,duzdyl_plus_duydzl
+
+ real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ) :: epsilondev_loc
+
+
+ do ispec = 1,NSPEC_CRUST_MANTLE
+
+ do k=1,NGLLZ
+ do j=1,NGLLY
+ do i=1,NGLLX
+
+ tempx1l = 0._CUSTOM_REAL
+ tempx2l = 0._CUSTOM_REAL
+ tempx3l = 0._CUSTOM_REAL
+
+ tempy1l = 0._CUSTOM_REAL
+ tempy2l = 0._CUSTOM_REAL
+ tempy3l = 0._CUSTOM_REAL
+
+ tempz1l = 0._CUSTOM_REAL
+ tempz2l = 0._CUSTOM_REAL
+ tempz3l = 0._CUSTOM_REAL
+
+ do l=1,NGLLX
+ hp1 = hprime_xx(i,l)
+ iglob = ibool(l,j,k,ispec)
+ tempx1l = tempx1l + displ_crust_mantle(1,iglob)*hp1
+ tempy1l = tempy1l + displ_crust_mantle(2,iglob)*hp1
+ tempz1l = tempz1l + displ_crust_mantle(3,iglob)*hp1
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ enddo
+
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ do l=1,NGLLY
+ hp2 = hprime_yy(j,l)
+ iglob = ibool(i,l,k,ispec)
+ tempx2l = tempx2l + displ_crust_mantle(1,iglob)*hp2
+ tempy2l = tempy2l + displ_crust_mantle(2,iglob)*hp2
+ tempz2l = tempz2l + displ_crust_mantle(3,iglob)*hp2
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ enddo
+
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ do l=1,NGLLZ
+ hp3 = hprime_zz(k,l)
+ iglob = ibool(i,j,l,ispec)
+ tempx3l = tempx3l + displ_crust_mantle(1,iglob)*hp3
+ tempy3l = tempy3l + displ_crust_mantle(2,iglob)*hp3
+ tempz3l = tempz3l + displ_crust_mantle(3,iglob)*hp3
+ enddo
+
+! get derivatives of ux, uy and uz with respect to x, y and z
+
+ xixl = xix(i,j,k,ispec)
+ xiyl = xiy(i,j,k,ispec)
+ xizl = xiz(i,j,k,ispec)
+ etaxl = etax(i,j,k,ispec)
+ etayl = etay(i,j,k,ispec)
+ etazl = etaz(i,j,k,ispec)
+ gammaxl = gammax(i,j,k,ispec)
+ gammayl = gammay(i,j,k,ispec)
+ gammazl = gammaz(i,j,k,ispec)
+
+! compute the jacobian
+ jacobianl = 1._CUSTOM_REAL / (xixl*(etayl*gammazl-etazl*gammayl) &
+ - xiyl*(etaxl*gammazl-etazl*gammaxl) &
+ + xizl*(etaxl*gammayl-etayl*gammaxl))
+
+ duxdxl = xixl*tempx1l + etaxl*tempx2l + gammaxl*tempx3l
+ duxdyl = xiyl*tempx1l + etayl*tempx2l + gammayl*tempx3l
+ duxdzl = xizl*tempx1l + etazl*tempx2l + gammazl*tempx3l
+
+ duydxl = xixl*tempy1l + etaxl*tempy2l + gammaxl*tempy3l
+ duydyl = xiyl*tempy1l + etayl*tempy2l + gammayl*tempy3l
+ duydzl = xizl*tempy1l + etazl*tempy2l + gammazl*tempy3l
+
+ duzdxl = xixl*tempz1l + etaxl*tempz2l + gammaxl*tempz3l
+ duzdyl = xiyl*tempz1l + etayl*tempz2l + gammayl*tempz3l
+ duzdzl = xizl*tempz1l + etazl*tempz2l + gammazl*tempz3l
+
+! precompute some sums to save CPU time
+ duxdxl_plus_duydyl = duxdxl + duydyl
+ duxdxl_plus_duzdzl = duxdxl + duzdzl
+ duydyl_plus_duzdzl = duydyl + duzdzl
+ duxdyl_plus_duydxl = duxdyl + duydxl
+ duzdxl_plus_duxdzl = duzdxl + duxdzl
+ duzdyl_plus_duydzl = duzdyl + duydzl
+
+! compute deviatoric strain
+ if(NSPEC_CRUST_MANTLE_STRAIN_ONLY == 1) then
+ ispec_strain = 1
+ else
+ ispec_strain = ispec
+ endif
+ epsilon_trace_over_3(i,j,k,ispec_strain) = ONE_THIRD * (duxdxl + duydyl + duzdzl)
+ epsilondev_loc(1,i,j,k) = duxdxl - epsilon_trace_over_3(i,j,k,ispec_strain)
+ epsilondev_loc(2,i,j,k) = duydyl - epsilon_trace_over_3(i,j,k,ispec_strain)
+ epsilondev_loc(3,i,j,k) = 0.5 * duxdyl_plus_duydxl
+ epsilondev_loc(4,i,j,k) = 0.5 * duzdxl_plus_duxdzl
+ epsilondev_loc(5,i,j,k) = 0.5 * duzdyl_plus_duydzl
+
+ enddo ! NGLLX
+ enddo ! NGLLY
+ enddo ! NGLLZ
+
+ !epsilondev(:,:,:,:,ispec) = epsilondev_loc(:,:,:,:)
+ do k=1,NGLLZ
+ do j=1,NGLLY
+ do i=1,NGLLX
+ epsilondev(:,i,j,k,ispec) = epsilondev_loc(:,i,j,k)
+ enddo
+ enddo
+ enddo
+
+ enddo ! spectral element loop NSPEC_CRUST_MANTLE
+
+ end subroutine compute_stain_crust_mantle
+
+
+!=====================================================================
+
+ subroutine compute_stain_inner_core(displ_inner_core,hprime_xx,hprime_yy,hprime_zz,ibool,&
+ xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz, &
+ epsilondev,epsilon_trace_over_3)
+
+ implicit none
+
+ include "constants.h"
+
+! include values created by the mesher
+! done for performance only using static allocation to allow for loop unrolling
+ include "OUTPUT_FILES/values_from_mesher.h"
+
+! array with derivatives of Lagrange polynomials and precalculated products
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLX) :: hprime_xx
+ real(kind=CUSTOM_REAL), dimension(NGLLY,NGLLY) :: hprime_yy
+ real(kind=CUSTOM_REAL), dimension(NGLLZ,NGLLZ) :: hprime_zz
+
+
+! arrays with mesh parameters per slice
+ integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: ibool
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: &
+ xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz
+
+ real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: epsilon_trace_over_3
+
+ real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: epsilondev
+
+ real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_INNER_CORE) :: displ_inner_core
+
+!local parameters
+ integer ispec,iglob,ispec_strain
+ integer i,j,k,l
+
+ real(kind=CUSTOM_REAL) tempx1l,tempx2l,tempx3l
+ real(kind=CUSTOM_REAL) tempy1l,tempy2l,tempy3l
+ real(kind=CUSTOM_REAL) tempz1l,tempz2l,tempz3l
+
+ real(kind=CUSTOM_REAL) hp1,hp2,hp3
+
+ real(kind=CUSTOM_REAL) xixl,xiyl,xizl,etaxl,etayl,etazl,gammaxl,gammayl,gammazl,jacobianl
+
+ real(kind=CUSTOM_REAL) duxdxl,duxdyl,duxdzl,duydxl,duydyl,duydzl,duzdxl,duzdyl,duzdzl
+
+ real(kind=CUSTOM_REAL) duxdxl_plus_duydyl,duxdxl_plus_duzdzl,duydyl_plus_duzdzl
+ real(kind=CUSTOM_REAL) duxdyl_plus_duydxl,duzdxl_plus_duxdzl,duzdyl_plus_duydzl
+
+ real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ) :: epsilondev_loc
+
+
+ do ispec = 1,NSPEC_INNER_CORE
+
+ do k=1,NGLLZ
+ do j=1,NGLLY
+ do i=1,NGLLX
+
+ tempx1l = 0._CUSTOM_REAL
+ tempx2l = 0._CUSTOM_REAL
+ tempx3l = 0._CUSTOM_REAL
+
+ tempy1l = 0._CUSTOM_REAL
+ tempy2l = 0._CUSTOM_REAL
+ tempy3l = 0._CUSTOM_REAL
+
+ tempz1l = 0._CUSTOM_REAL
+ tempz2l = 0._CUSTOM_REAL
+ tempz3l = 0._CUSTOM_REAL
+
+ do l=1,NGLLX
+ hp1 = hprime_xx(i,l)
+ iglob = ibool(l,j,k,ispec)
+ tempx1l = tempx1l + displ_inner_core(1,iglob)*hp1
+ tempy1l = tempy1l + displ_inner_core(2,iglob)*hp1
+ tempz1l = tempz1l + displ_inner_core(3,iglob)*hp1
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ enddo
+
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ do l=1,NGLLY
+ hp2 = hprime_yy(j,l)
+ iglob = ibool(i,l,k,ispec)
+ tempx2l = tempx2l + displ_inner_core(1,iglob)*hp2
+ tempy2l = tempy2l + displ_inner_core(2,iglob)*hp2
+ tempz2l = tempz2l + displ_inner_core(3,iglob)*hp2
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ enddo
+
+!!! can merge these loops because NGLLX = NGLLY = NGLLZ do l=1,NGLLZ
+ hp3 = hprime_zz(k,l)
+ iglob = ibool(i,j,l,ispec)
+ tempx3l = tempx3l + displ_inner_core(1,iglob)*hp3
+ tempy3l = tempy3l + displ_inner_core(2,iglob)*hp3
+ tempz3l = tempz3l + displ_inner_core(3,iglob)*hp3
+ enddo
+
+! get derivatives of ux, uy and uz with respect to x, y and z
+
+ xixl = xix(i,j,k,ispec)
+ xiyl = xiy(i,j,k,ispec)
+ xizl = xiz(i,j,k,ispec)
+ etaxl = etax(i,j,k,ispec)
+ etayl = etay(i,j,k,ispec)
+ etazl = etaz(i,j,k,ispec)
+ gammaxl = gammax(i,j,k,ispec)
+ gammayl = gammay(i,j,k,ispec)
+ gammazl = gammaz(i,j,k,ispec)
+
+! compute the jacobian
+ jacobianl = 1._CUSTOM_REAL / (xixl*(etayl*gammazl-etazl*gammayl) &
+ - xiyl*(etaxl*gammazl-etazl*gammaxl) &
+ + xizl*(etaxl*gammayl-etayl*gammaxl))
+
+ duxdxl = xixl*tempx1l + etaxl*tempx2l + gammaxl*tempx3l
+ duxdyl = xiyl*tempx1l + etayl*tempx2l + gammayl*tempx3l
+ duxdzl = xizl*tempx1l + etazl*tempx2l + gammazl*tempx3l
+
+ duydxl = xixl*tempy1l + etaxl*tempy2l + gammaxl*tempy3l
+ duydyl = xiyl*tempy1l + etayl*tempy2l + gammayl*tempy3l
+ duydzl = xizl*tempy1l + etazl*tempy2l + gammazl*tempy3l
+
+ duzdxl = xixl*tempz1l + etaxl*tempz2l + gammaxl*tempz3l
+ duzdyl = xiyl*tempz1l + etayl*tempz2l + gammayl*tempz3l
+ duzdzl = xizl*tempz1l + etazl*tempz2l + gammazl*tempz3l
+
+! precompute some sums to save CPU time
+ duxdxl_plus_duydyl = duxdxl + duydyl
+ duxdxl_plus_duzdzl = duxdxl + duzdzl
+ duydyl_plus_duzdzl = duydyl + duzdzl
+ duxdyl_plus_duydxl = duxdyl + duydxl
+ duzdxl_plus_duxdzl = duzdxl + duxdzl
+ duzdyl_plus_duydzl = duzdyl + duydzl
+
+! compute deviatoric strain
+ if(NSPEC_INNER_CORE_STRAIN_ONLY == 1) then
+ ispec_strain = 1
+ else
+ ispec_strain = ispec
+ endif
+ epsilon_trace_over_3(i,j,k,ispec_strain) = ONE_THIRD * (duxdxl + duydyl + duzdzl)
+ epsilondev_loc(1,i,j,k) = duxdxl - epsilon_trace_over_3(i,j,k,ispec_strain)
+ epsilondev_loc(2,i,j,k) = duydyl - epsilon_trace_over_3(i,j,k,ispec_strain)
+ epsilondev_loc(3,i,j,k) = 0.5 * duxdyl_plus_duydxl
+ epsilondev_loc(4,i,j,k) = 0.5 * duzdxl_plus_duxdzl
+ epsilondev_loc(5,i,j,k) = 0.5 * duzdyl_plus_duydzl
+
+ enddo ! NGLLX
+ enddo ! NGLLY
+ enddo ! NGLLZ
+
+ !epsilondev(:,:,:,:,ispec) = epsilondev_loc(:,:,:,:)
+ do k=1,NGLLZ
+ do j=1,NGLLY
+ do i=1,NGLLX
+ epsilondev(:,i,j,k,ispec) = epsilondev_loc(:,i,j,k)
+ enddo
+ enddo
+ enddo
+
+ enddo ! spectral element loop NSPEC_CRUST_MANTLE
+
+ end subroutine compute_stain_inner_core
+
+
+
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/specfem3D.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/specfem3D.F90 2013-06-11 23:37:32 UTC (rev 22226)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/undo_att/src/specfem3D/specfem3D.F90 2013-06-12 08:31:00 UTC (rev 22227)
@@ -2215,6 +2215,16 @@
endif
+ if(SIMULATION_TYPE == 2)then
+ !!add this part
+!ZN !ZN we need to be careful to arrange this part
+!ZN ! save source derivatives for adjoint simulations
+!ZN if (SIMULATION_TYPE == 2 .and. nrec_local > 0) then
+!ZN call save_kernels_source_derivatives(nrec_local,NSOURCES,scale_displ,scale_t, &
+!ZN nu_source,moment_der,sloc_der,stshift_der,shdur_der,number_receiver_global)
+!ZN endif
+ endif
+
if(SIMULATION_TYPE == 3)then
it = 0
@@ -2264,8 +2274,34 @@
include "part1_undo_att.F90"
+ call compute_stain_crust_mantle(displ_crust_mantle,hprime_xx,hprime_yy,hprime_zz,ibool_crust_mantle,&
+ 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, &
+ epsilondev_crust_mantle,eps_trace_over_3_crust_mantle)
+
+ call compute_stain_crust_mantle(b_displ_crust_mantle,hprime_xx,hprime_yy,hprime_zz,ibool_crust_mantle,&
+ 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, &
+ b_epsilondev_crust_mantle,b_eps_trace_over_3_crust_mantle)
+
+ call compute_stain_inner_core(displ_inner_core,hprime_xx,hprime_yy,hprime_zz,ibool_inner_core,&
+ xix_inner_core,xiy_inner_core,xiz_inner_core,&
+ etax_inner_core,etay_inner_core,etaz_inner_core,&
+ gammax_inner_core,gammay_inner_core,gammaz_inner_core, &
+ epsilondev_inner_core,eps_trace_over_3_inner_core)
+
+ call compute_stain_inner_core(b_displ_inner_core,hprime_xx,hprime_yy,hprime_zz,ibool_inner_core,&
+ xix_inner_core,xiy_inner_core,xiz_inner_core,&
+ etax_inner_core,etay_inner_core,etaz_inner_core,&
+ gammax_inner_core,gammay_inner_core,gammaz_inner_core, &
+ b_epsilondev_inner_core,b_eps_trace_over_3_inner_core)
+
include "part3_kernel_computation.F90"
+ include "part4_save_kernel.F90"
+
enddo
enddo ! end of main time loop
@@ -2279,7 +2315,7 @@
!-------------------------------------------------------------------------------------------------
!-------------------------------------------------------------------------------------------------
!-------------------------------------------------------------------------------------------------
-
+ !ZN need to be remove for undoing att
! synchronize all processes, waits until all processes have written their seismograms
call MPI_BARRIER(MPI_COMM_WORLD,ier)
if( ier /= 0 ) call exit_mpi(myrank,'error synchronize after time loop')
@@ -2373,6 +2409,9 @@
call MPI_BARRIER(MPI_COMM_WORLD,ier)
if( ier /= 0 ) call exit_mpi(myrank,'error synchronize saving forward')
+#ifdef UNDO_ATT
+ !ZN we move this part of code(save kernels)inside the time loop above
+#else
! dump kernel arrays
if (SIMULATION_TYPE == 3) then
@@ -2427,11 +2466,13 @@
endif
endif
+ !ZN we need to be careful to arrange this part
! save source derivatives for adjoint simulations
if (SIMULATION_TYPE == 2 .and. nrec_local > 0) then
call save_kernels_source_derivatives(nrec_local,NSOURCES,scale_displ,scale_t, &
nu_source,moment_der,sloc_der,stshift_der,shdur_der,number_receiver_global)
endif
+#endif
! frees dynamically allocated memory
! mpi buffers
More information about the CIG-COMMITS
mailing list