[cig-commits] r22664 - in seismo/3D/SPECFEM3D_GLOBE: branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D trunk/src/specfem3D
dkomati1 at geodynamics.org
dkomati1 at geodynamics.org
Wed Jul 24 10:05:50 PDT 2013
Author: dkomati1
Date: 2013-07-24 10:05:49 -0700 (Wed, 24 Jul 2013)
New Revision: 22664
Modified:
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube_block.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_scalar_block.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_vector_block.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_boundary_kernel.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_noDev.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_noDev.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_outer_core_noDev.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_kernels.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_seismograms.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_crust_mantle.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_outer_core.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/create_central_cube_buffers.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_backazimuth.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_event_info.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/initialize_simulation.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/iterate_time.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_receivers.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_regular_points.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_sources.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/noise_tomography.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_buffers_solver.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_solver.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_mesh_databases.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_kernels.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_regular_kernels.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/setup_sources_receivers.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/specfem3D.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_surface.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_volume.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_seismograms.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.F90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/initialize_simulation.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_sources.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_volume.f90
Log:
done with all the easy merges in src/specfem3D
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -510,15 +510,15 @@
!daniel: att - debug - check if mimic is still needed
if( ATTENUATION ) then
- write(IOUT,*) 'logical, parameter :: USE_ATTENUATION_MIMIC = .true.'
+ write(IOUT,*) 'logical, parameter :: PARTIAL_PHYS_DISPERSION_ONLY = .true.'
else
- write(IOUT,*) 'logical, parameter :: USE_ATTENUATION_MIMIC = .false.'
+ write(IOUT,*) 'logical, parameter :: PARTIAL_PHYS_DISPERSION_ONLY = .false.'
endif
else
! calculates full attenuation (phase & amplitude effects) if used
- write(IOUT,*) 'logical, parameter :: USE_ATTENUATION_MIMIC = .false.'
+ write(IOUT,*) 'logical, parameter :: PARTIAL_PHYS_DISPERSION_ONLY = .false.'
endif
! attenuation and/or adjoint simulations
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -30,10 +30,10 @@
receiver_cube_from_slices,ibool_inner_core,idoubling_inner_core, &
ibelm_bottom_inner_core,NSPEC2D_BOTTOM_INNER_CORE,vector_assemble,ndim_assemble,iphase_CC)
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
include 'constants.h'
! include values created by the mesher
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube_block.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube_block.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_central_cube_block.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -32,10 +32,10 @@
! this version of the routine is based on blocking MPI calls
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
include 'constants.h'
! for matching with central cube in inner core
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_scalar_block.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_scalar_block.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_scalar_block.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -45,11 +45,10 @@
! this version of the routine is based on blocking MPI calls
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
include "precision.h"
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_vector_block.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_vector_block.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/assemble_MPI_vector_block.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -51,11 +51,10 @@
! this version of the routine is based on blocking MPI calls
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
include "precision.h"
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_boundary_kernel.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_boundary_kernel.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_boundary_kernel.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -360,7 +360,6 @@
real(kind=CUSTOM_REAL), dimension(NDIM,*) :: displ,accel,b_displ
integer nspec, iregion_code
integer, dimension(NGLLX,NGLLY,NGLLZ,nspec) :: ibool
-! integer, dimension(*) :: idoubling
logical, dimension(*) :: ispec_is_tiso
real(kind=CUSTOM_REAL), dimension(*) :: ystore,zstore
@@ -447,7 +446,7 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- ystore,zstore,ibool,ispec_is_tiso) !idoubling)
+ ystore,zstore,ibool,ispec_is_tiso)
! ----- forward strain -------
temp1(:) = matmul(b_displl(:,:,j,k), hprime_xx(i,:))
@@ -472,7 +471,7 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- ystore,zstore,ibool,ispec_is_tiso) !-- idoubling)
+ ystore,zstore,ibool,ispec_is_tiso)
! ---- precompute K_d for F-S boundaries ----
if (fluid_solid_boundary) then
@@ -554,9 +553,8 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- ystore,zstore,ibool,ispec_is_tiso) !--idoubling)
+ ystore,zstore,ibool,ispec_is_tiso)
-
implicit none
include 'constants.h'
@@ -572,7 +570,6 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store
real(kind=CUSTOM_REAL), dimension(*) :: ystore,zstore
integer, dimension(NGLLX,NGLLY,NGLLZ,*) :: ibool
-! integer, dimension(*) :: idoubling
logical, dimension(*) :: ispec_is_tiso
! --- local variables ---
@@ -651,7 +648,6 @@
c45*duzdxl_plus_duxdzl + c44*duzdyl_plus_duydzl + c34*duzdzl
else if( .not. ispec_is_tiso(ispec) ) then
-!else if(.not. (TRANSVERSE_ISOTROPY_VAL .and. (idoubling(ispec) == IFLAG_80_MOHO .or. idoubling(ispec) == IFLAG_220_80))) then
! isotropic elements
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -42,8 +42,6 @@
tempz1_att,tempz2_att,tempz3_att, &
epsilondev_loc,rho_s_H,is_backward_field)
-! this routine is optimized for NGLLX = NGLLY = NGLLZ = 5 using the Deville et al. (2002) inlined matrix-matrix products
-
implicit none
include "constants.h"
@@ -129,7 +127,7 @@
integer :: ispec_strain
integer :: i,j,k
integer :: int_radius
- integer :: iglob1
+ integer :: iglob
! isotropic element
@@ -255,7 +253,7 @@
sigma_yz = mul*duzdyl_plus_duydzl
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
!daniel: att - debug update
! call compute_element_att_mem_up_cm(ispec,i,j,k, &
@@ -286,13 +284,12 @@
! compute non-symmetric terms for gravity
if(GRAVITY_VAL) then
-
! use mesh coordinates to get theta and phi
! x y and z contain r theta and phi
- iglob1 = ibool(i,j,k,ispec)
+ iglob = ibool(i,j,k,ispec)
- dtheta = dble(ystore(iglob1))
- dphi = dble(zstore(iglob1))
+ dtheta = dble(ystore(iglob))
+ dphi = dble(zstore(iglob))
cos_theta = dcos(dtheta)
sin_theta = dsin(dtheta)
@@ -307,7 +304,7 @@
! get g, rho and dg/dr=dg
! spherical components of the gravitational acceleration
! for efficiency replace with lookup table every 100 m in radial direction
- radius = dble(xstore(iglob1))
+ radius = dble(xstore(iglob))
int_radius = nint(10.d0 * radius * R_EARTH_KM )
minus_g = minus_gravity_table(int_radius)
@@ -335,9 +332,9 @@
if(CUSTOM_REAL == SIZE_REAL) then
! get displacement and multiply by density to compute G tensor
- sx_l = rho * dble(dummyx_loc(i,j,k)) ! dble(displ_crust_mantle(1,iglob1))
- sy_l = rho * dble(dummyy_loc(i,j,k)) ! dble(displ_crust_mantle(2,iglob1))
- sz_l = rho * dble(dummyz_loc(i,j,k)) ! dble(displ_crust_mantle(3,iglob1))
+ sx_l = rho * dble(dummyx_loc(i,j,k)) ! dble(displ_crust_mantle(1,iglob))
+ sy_l = rho * dble(dummyy_loc(i,j,k)) ! dble(displ_crust_mantle(2,iglob))
+ sz_l = rho * dble(dummyz_loc(i,j,k)) ! dble(displ_crust_mantle(3,iglob))
! compute G tensor from s . g and add to sigma (not symmetric)
sigma_xx = sigma_xx + sngl(sy_l*gyl + sz_l*gzl)
@@ -362,9 +359,9 @@
else
! get displacement and multiply by density to compute G tensor
- sx_l = rho * dummyx_loc(i,j,k) ! displ_crust_mantle(1,iglob1)
- sy_l = rho * dummyy_loc(i,j,k) ! displ_crust_mantle(2,iglob1)
- sz_l = rho * dummyz_loc(i,j,k) ! displ_crust_mantle(3,iglob1)
+ sx_l = rho * dummyx_loc(i,j,k) ! displ_crust_mantle(1,iglob)
+ sy_l = rho * dummyy_loc(i,j,k) ! displ_crust_mantle(2,iglob)
+ sz_l = rho * dummyz_loc(i,j,k) ! displ_crust_mantle(3,iglob)
! compute G tensor from s . g and add to sigma (not symmetric)
sigma_xx = sigma_xx + sy_l*gyl + sz_l*gzl
@@ -533,7 +530,7 @@
integer :: ispec_strain
integer :: i,j,k
integer :: int_radius
- integer :: iglob1
+ integer :: iglob
! transverse isotropic element
@@ -631,17 +628,11 @@
! compute either isotropic or anisotropic elements
!
-! note : mesh is built such that anisotropic elements are created first in anisotropic layers,
+! note: the mesh is built such that anisotropic elements are created first in anisotropic layers,
! thus they are listed first ( see in create_regions_mesh.f90: perm_layer() ordering )
! this is therefore still in bounds of 1:NSPECMAX_TISO_MANTLE even if NSPECMAX_TISO is less than NSPEC
- ! uncomment to debug
- !if ( ispec > NSPECMAX_TISO_MANTLE ) then
- ! print*,'error tiso: ispec = ',ispec,'max = ',NSPECMAX_TISO_MANTLE
- ! call exit_mpi(0,'error tiso ispec bounds')
- !endif
-
- ! use Kappa and mu from transversely isotropic model
+ ! use kappa and mu from transversely isotropic model
kappavl = kappavstore(i,j,k,ispec)
muvl = muvstore(i,j,k,ispec)
@@ -671,10 +662,10 @@
! use mesh coordinates to get theta and phi
! ystore and zstore contain theta and phi
- iglob1 = ibool(i,j,k,ispec)
+ iglob = ibool(i,j,k,ispec)
- theta = ystore(iglob1)
- phi = zstore(iglob1)
+ theta = ystore(iglob)
+ phi = zstore(iglob)
! precompute some products to reduce the CPU time
@@ -719,8 +710,7 @@
four_rhovsvsq = 4.0_CUSTOM_REAL*rhovsvsq
four_rhovshsq = 4.0_CUSTOM_REAL*rhovshsq
-
- ! way 2: pre-compute temporary values
+ ! pre-compute temporary values
templ1 = four_rhovsvsq - rhovpvsq + twoetaminone*rhovphsq - four_eta_aniso*rhovsvsq
templ1_cos = rhovphsq - rhovpvsq + costwotheta*templ1
templ2 = four_rhovsvsq - rhovpvsq - rhovphsq + two_eta_aniso*rhovphsq - four_eta_aniso*rhovsvsq
@@ -729,7 +719,7 @@
templ3_two = templ3 - two_rhovshsq - two_rhovsvsq
templ3_cos = templ3_two + costwotheta*templ2
- ! way 2: reordering operations to facilitate compilation, avoiding divisions, using locality for temporary values
+ ! reordering operations to facilitate compilation, avoiding divisions, using locality for temporary values
c11 = rhovphsq*sinphifour &
+ 2.0_CUSTOM_REAL*cosphisq*sinphisq* &
( rhovphsq*costhetasq + sinthetasq*(eta_aniso*rhovphsq + two_rhovsvsq - two_eta_aniso*rhovsvsq) ) &
@@ -851,7 +841,7 @@
c45*duzdxl_plus_duxdzl + c44*duzdyl_plus_duydzl + c34*duzdzl
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
!daniel: att - debug update
! call compute_element_att_mem_up_cm(ispec,i,j,k, &
@@ -882,21 +872,19 @@
! compute non-symmetric terms for gravity
if(GRAVITY_VAL) then
-
- ! use mesh coordinates to get theta and phi
+ ! use mesh coordinates to get theta and phi
! x y and z contain r theta and phi
- iglob1 = ibool(i,j,k,ispec)
+ iglob = ibool(i,j,k,ispec)
- dtheta = dble(ystore(iglob1))
- dphi = dble(zstore(iglob1))
- radius = dble(xstore(iglob1))
+ dtheta = dble(ystore(iglob))
+ dphi = dble(zstore(iglob))
+ radius = dble(xstore(iglob))
cos_theta = dcos(dtheta)
sin_theta = dsin(dtheta)
cos_phi = dcos(dphi)
sin_phi = dsin(dphi)
- ! way 2
cos_theta_sq = cos_theta*cos_theta
sin_theta_sq = sin_theta*sin_theta
cos_phi_sq = cos_phi*cos_phi
@@ -1029,8 +1017,6 @@
tempz1_att,tempz2_att,tempz3_att, &
epsilondev_loc,rho_s_H,is_backward_field)
-! this routine is optimized for NGLLX = NGLLY = NGLLZ = 5 using the Deville et al. (2002) inlined matrix-matrix products
-
implicit none
include "constants.h"
@@ -1121,7 +1107,7 @@
integer :: ispec_strain
integer :: i,j,k
integer :: int_radius
- integer :: iglob1
+ integer :: iglob
! anisotropic elements
@@ -1280,21 +1266,8 @@
c45*duzdxl_plus_duxdzl + c44*duzdyl_plus_duydzl + c34*duzdzl
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
-!daniel: att - debug update
-! call compute_element_att_mem_up_cm(ispec,i,j,k, &
-! R_xx(:,i,j,k,ispec), &
-! R_yy(:,i,j,k,ispec), &
-! R_xy(:,i,j,k,ispec), &
-! R_xz(:,i,j,k,ispec), &
-! R_yz(:,i,j,k,ispec), &
-! epsilondev_loc(:,i,j,k),c44store(i,j,k,ispec),is_backward_field)
-! dummy to avoid compiler warning
- if( is_backward_field ) then
- endif
-
-
! note: fortran passes pointers to array location, thus R_memory(1,1,...) should be fine
call compute_element_att_stress(R_xx(1,i,j,k,ispec), &
R_yy(1,i,j,k,ispec), &
@@ -1312,14 +1285,13 @@
! compute non-symmetric terms for gravity
if(GRAVITY_VAL) then
-
- ! use mesh coordinates to get theta and phi
+ ! use mesh coordinates to get theta and phi
! x y and z contain r theta and phi
- iglob1 = ibool(i,j,k,ispec)
+ iglob = ibool(i,j,k,ispec)
- dtheta = dble(ystore(iglob1))
- dphi = dble(zstore(iglob1))
- radius = dble(xstore(iglob1))
+ dtheta = dble(ystore(iglob))
+ dphi = dble(zstore(iglob))
+ radius = dble(xstore(iglob))
cos_theta = dcos(dtheta)
sin_theta = dsin(dtheta)
@@ -1360,9 +1332,9 @@
if(CUSTOM_REAL == SIZE_REAL) then
! get displacement and multiply by density to compute G tensor
- sx_l = rho * dble(dummyx_loc(i,j,k)) ! dble(displ_crust_mantle(1,iglob1))
- sy_l = rho * dble(dummyy_loc(i,j,k)) ! dble(displ_crust_mantle(2,iglob1))
- sz_l = rho * dble(dummyz_loc(i,j,k)) ! dble(displ_crust_mantle(3,iglob1))
+ sx_l = rho * dble(dummyx_loc(i,j,k)) ! dble(displ_crust_mantle(1,iglob))
+ sy_l = rho * dble(dummyy_loc(i,j,k)) ! dble(displ_crust_mantle(2,iglob))
+ sz_l = rho * dble(dummyz_loc(i,j,k)) ! dble(displ_crust_mantle(3,iglob))
! compute G tensor from s . g and add to sigma (not symmetric)
sigma_xx = sigma_xx + sngl(sy_l*gyl + sz_l*gzl)
@@ -1387,9 +1359,9 @@
else
! get displacement and multiply by density to compute G tensor
- sx_l = rho * dummyx_loc(i,j,k) ! displ_crust_mantle(1,iglob1)
- sy_l = rho * dummyy_loc(i,j,k) ! displ_crust_mantle(2,iglob1)
- sz_l = rho * dummyz_loc(i,j,k) ! displ_crust_mantle(3,iglob1)
+ sx_l = rho * dummyx_loc(i,j,k) ! displ_crust_mantle(1,iglob)
+ sy_l = rho * dummyy_loc(i,j,k) ! displ_crust_mantle(2,iglob)
+ sz_l = rho * dummyz_loc(i,j,k) ! displ_crust_mantle(3,iglob)
! compute G tensor from s . g and add to sigma (not symmetric)
sigma_xx = sigma_xx + sy_l*gyl + sz_l*gzl
@@ -1507,9 +1479,6 @@
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
-!daniel: att - debug predictor
-! use specfem_par,only: tau_sigma_dble,deltat,b_deltat
-
implicit none
include "constants.h"
@@ -1553,8 +1522,8 @@
! IMPROVE we use mu_v here even if there is some anisotropy
! IMPROVE we should probably use an average value instead
-!daniel: att - debug original
do i_SLS = 1,N_SLS
+
! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
if( USE_3D_ATTENUATION_ARRAYS ) then
if(ANISOTROPIC_3D_MANTLE_VAL) then
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -668,7 +668,7 @@
endif ! end of test whether isotropic or anisotropic element
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
do i_SLS = 1,N_SLS
R_xx_val = R_xx(i_SLS,i,j,k,ispec)
R_yy_val = R_yy(i_SLS,i,j,k,ispec)
@@ -887,7 +887,7 @@
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. )) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. )) then
! use Runge-Kutta scheme to march in time
do i_SLS = 1,N_SLS
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -395,7 +395,6 @@
!
if(ANISOTROPIC_3D_MANTLE_VAL) then
! anisotropic element
-
call compute_element_aniso(ispec, &
minus_gravity_table,density_table,minus_deriv_gravity_table, &
xstore,ystore,zstore, &
@@ -415,10 +414,8 @@
tempz1_att,tempz2_att,tempz3_att, &
epsilondev_loc,rho_s_H,is_backward_field)
else
-
- if( .not. ispec_is_tiso(ispec) ) then
+ if(.not. ispec_is_tiso(ispec)) then
! isotropic element
-
call compute_element_iso(ispec, &
minus_gravity_table,density_table,minus_deriv_gravity_table, &
xstore,ystore,zstore, &
@@ -437,7 +434,6 @@
epsilondev_loc,rho_s_H,is_backward_field)
else
! transverse isotropic element
-
call compute_element_tiso(ispec, &
minus_gravity_table,density_table,minus_deriv_gravity_table, &
xstore,ystore,zstore, &
@@ -572,7 +568,7 @@
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
!daniel: att - debug update R_memory variable only if not last time step which will be saved..
! if( .not. it == NSTEP ) then
@@ -597,7 +593,7 @@
epsilondev_yz(:,:,:,ispec) = epsilondev_loc(5,:,:,:)
endif
- enddo ! spectral element loop NSPEC_CRUST_MANTLE
+ enddo ! of spectral element loop NSPEC_CRUST_MANTLE
end subroutine compute_forces_crust_mantle_Dev
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_noDev.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_noDev.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_noDev.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -64,28 +64,7 @@
! done for performance only using static allocation to allow for loop unrolling
include "OUTPUT_FILES/values_from_mesher.h"
-! model_attenuation_variables
-! type model_attenuation_variables
-! sequence
-! double precision min_period, max_period
-! double precision :: QT_c_source ! Source Frequency
-! double precision, dimension(:), pointer :: Qtau_s ! tau_sigma
-! double precision, dimension(:), pointer :: QrDisc ! Discontinutitues Defined
-! double precision, dimension(:), pointer :: Qr ! Radius
-! double precision, dimension(:), pointer :: Qmu ! Shear Attenuation
-! double precision, dimension(:,:), pointer :: Qtau_e ! tau_epsilon
-! double precision, dimension(:), pointer :: Qomsb, Qomsb2 ! one_minus_sum_beta
-! double precision, dimension(:,:), pointer :: Qfc, Qfc2 ! factor_common
-! double precision, dimension(:), pointer :: Qsf, Qsf2 ! scale_factor
-! integer, dimension(:), pointer :: Qrmin ! Max and Mins of idoubling
-! integer, dimension(:), pointer :: Qrmax ! Max and Mins of idoubling
-! integer, dimension(:), pointer :: interval_Q ! Steps
-! integer :: Qn ! Number of points
-! integer dummy_pad ! padding 4 bytes to align the structure
-! end type model_attenuation_variables
-
! array with the local to global mapping per slice
-! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling
logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso
! displacement and acceleration
@@ -469,7 +448,6 @@
else
! do not use transverse isotropy except if element is between d220 and Moho
-! if(.not. (TRANSVERSE_ISOTROPY_VAL .and. (idoubling(ispec)==IFLAG_220_80 .or. idoubling(ispec)==IFLAG_80_MOHO))) then
if( .not. ispec_is_tiso(ispec) ) then
! layer with no transverse isotropy, use kappav and muv
kappal = kappavstore(i,j,k,ispec)
@@ -701,7 +679,7 @@
endif ! end of test whether isotropic or anisotropic element
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. .not. PARTIAL_PHYS_DISPERSION_ONLY_VAL) then
do i_SLS = 1,N_SLS
R_xx_val = R_memory(1,i_SLS,i,j,k,ispec)
R_yy_val = R_memory(2,i_SLS,i,j,k,ispec)
@@ -922,7 +900,7 @@
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
- if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. )) then
+ if(ATTENUATION_VAL .and. .not. PARTIAL_PHYS_DISPERSION_ONLY_VAL) then
call compute_element_strain_att_noDev(ispec,NGLOB_CRUST_MANTLE,NSPEC_CRUST_MANTLE,displ_crust_mantle,veloc_crust_mantle,&
deltat,hprime_xx,hprime_yy,hprime_zz,ibool,&
@@ -958,8 +936,8 @@
BETA_LDDRK(istage) * R_memory_lddrk(i_memory,i_SLS,:,:,:,ispec)
else
R_memory(i_memory,i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_memory(i_memory,i_SLS,:,:,:,ispec) &
- + factor_common_c44_muv(:,:,:) &
- * (betaval(i_SLS) * epsilondev_loc(i_memory,:,:,:) + gammaval(i_SLS) * epsilondev_loc_nplus1(i_memory,:,:,:))
+ + factor_common_c44_muv(:,:,:) &
+ * (betaval(i_SLS) * epsilondev_loc(i_memory,:,:,:) + gammaval(i_SLS) * epsilondev_loc_nplus1(i_memory,:,:,:))
endif
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -424,7 +424,7 @@
endif
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
do i_SLS = 1,N_SLS
R_xx_val = R_xx(i_SLS,i,j,k,ispec)
R_yy_val = R_yy(i_SLS,i,j,k,ispec)
@@ -649,7 +649,7 @@
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. )) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. )) then
do i_SLS = 1,N_SLS
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -191,7 +191,7 @@
integer :: int_radius
integer :: ispec,ispec_strain
integer :: i,j,k
- integer :: iglob1
+ integer :: iglob
! integer :: computed_elements
integer :: num_elements,ispec_p
integer :: iphase
@@ -226,10 +226,10 @@
do k=1,NGLLZ
do j=1,NGLLY
do i=1,NGLLX
- iglob1 = ibool(i,j,k,ispec)
- dummyx_loc(i,j,k) = displ_inner_core(1,iglob1)
- dummyy_loc(i,j,k) = displ_inner_core(2,iglob1)
- dummyz_loc(i,j,k) = displ_inner_core(3,iglob1)
+ iglob = ibool(i,j,k,ispec)
+ dummyx_loc(i,j,k) = displ_inner_core(1,iglob)
+ dummyy_loc(i,j,k) = displ_inner_core(2,iglob)
+ dummyz_loc(i,j,k) = displ_inner_core(3,iglob)
enddo
enddo
enddo
@@ -243,10 +243,10 @@
do k=1,NGLLZ
do j=1,NGLLY
do i=1,NGLLX
- iglob1 = ibool(i,j,k,ispec)
- dummyx_loc_att(i,j,k) = deltat*veloc_inner_core(1,iglob1)
- dummyy_loc_att(i,j,k) = deltat*veloc_inner_core(2,iglob1)
- dummyz_loc_att(i,j,k) = deltat*veloc_inner_core(3,iglob1)
+ iglob = ibool(i,j,k,ispec)
+ dummyx_loc_att(i,j,k) = deltat*veloc_inner_core(1,iglob)
+ dummyy_loc_att(i,j,k) = deltat*veloc_inner_core(2,iglob)
+ dummyz_loc_att(i,j,k) = deltat*veloc_inner_core(3,iglob)
enddo
enddo
enddo
@@ -566,7 +566,7 @@
sigma_xx = c11l*duxdxl + c12l*duydyl + c13l*duzdzl
sigma_yy = c12l*duxdxl + c11l*duydyl + c13l*duzdzl
sigma_zz = c13l*duxdxl + c13l*duydyl + c33l*duzdzl
- sigma_xy = 0.5*(c11l-c12l)*duxdyl_plus_duydxl
+ sigma_xy = 0.5_CUSTOM_REAL*(c11l-c12l)*duxdyl_plus_duydxl
sigma_xz = c44l*duzdxl_plus_duxdzl
sigma_yz = c44l*duzdyl_plus_duydzl
else
@@ -586,7 +586,7 @@
endif
lambdalplus2mul = kappal + FOUR_THIRDS * mul
- lambdal = lambdalplus2mul - 2.*mul
+ lambdal = lambdalplus2mul - 2._CUSTOM_REAL*mul
! compute stress sigma
sigma_xx = lambdalplus2mul*duxdxl + lambdal*duydyl_plus_duzdzl
@@ -600,7 +600,7 @@
endif
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
! note: fortran passes pointers to array location, thus R_memory(1,1,...) should be fine
call compute_element_att_stress(R_xx(1,i,j,k,ispec), &
@@ -622,10 +622,10 @@
! use mesh coordinates to get theta and phi
! x y and z contain r theta and phi
- iglob1 = ibool(i,j,k,ispec)
- radius = dble(xstore(iglob1))
- theta = dble(ystore(iglob1))
- phi = dble(zstore(iglob1))
+ iglob = ibool(i,j,k,ispec)
+ radius = dble(xstore(iglob))
+ theta = dble(ystore(iglob))
+ phi = dble(zstore(iglob))
! make sure radius is never zero even for points at center of cube
! because we later divide by radius
@@ -668,14 +668,14 @@
Hyzl = cos_theta*minus_dg_plus_g_over_radius*sin_phi*sin_theta
! for locality principle, we set iglob again, in order to have it in the cache again
- iglob1 = ibool(i,j,k,ispec)
+ iglob = ibool(i,j,k,ispec)
! distinguish between single and double precision for reals
if(CUSTOM_REAL == SIZE_REAL) then
! get displacement and multiply by density to compute G tensor
- sx_l = rho * dble(displ_inner_core(1,iglob1))
- sy_l = rho * dble(displ_inner_core(2,iglob1))
- sz_l = rho * dble(displ_inner_core(3,iglob1))
+ sx_l = rho * dble(displ_inner_core(1,iglob))
+ sy_l = rho * dble(displ_inner_core(2,iglob))
+ sz_l = rho * dble(displ_inner_core(3,iglob))
! compute G tensor from s . g and add to sigma (not symmetric)
sigma_xx = sigma_xx + sngl(sy_l*gyl + sz_l*gzl)
@@ -700,9 +700,9 @@
else
! get displacement and multiply by density to compute G tensor
- sx_l = rho * displ_inner_core(1,iglob1)
- sy_l = rho * displ_inner_core(2,iglob1)
- sz_l = rho * displ_inner_core(3,iglob1)
+ sx_l = rho * displ_inner_core(1,iglob)
+ sy_l = rho * displ_inner_core(2,iglob)
+ sz_l = rho * displ_inner_core(3,iglob)
! compute G tensor from s . g and add to sigma (not symmetric)
sigma_xx = sigma_xx + sy_l*gyl + sz_l*gzl
@@ -839,8 +839,8 @@
do k=1,NGLLZ
do j=1,NGLLY
do i=1,NGLLX
- iglob1 = ibool(i,j,k,ispec)
- accel_inner_core(:,iglob1) = accel_inner_core(:,iglob1) + sum_terms(:,i,j,k)
+ iglob = ibool(i,j,k,ispec)
+ accel_inner_core(:,iglob) = accel_inner_core(:,iglob) + sum_terms(:,i,j,k)
enddo
enddo
enddo
@@ -859,7 +859,7 @@
! equation (9.59) page 350): Q_\alpha = Q_\mu * 3 * (V_p/V_s)^2 / 4
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
- if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
! updates R_memory
call compute_element_att_memory_ic(ispec,R_xx,R_yy,R_xy,R_xz,R_yz, &
@@ -881,9 +881,9 @@
epsilondev_yz(:,:,:,ispec) = epsilondev_loc(5,:,:,:)
endif
- endif ! end test to exclude fictitious elements in central cube
+ endif ! end of test to exclude fictitious elements in central cube
- enddo ! spectral element loop
+ enddo ! of spectral element loop
end subroutine compute_forces_inner_core_Dev
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_noDev.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_noDev.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_noDev.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -431,7 +431,7 @@
endif
! subtract memory variables if attenuation
- if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. ) ) then
+ if(ATTENUATION_VAL .and. .not. PARTIAL_PHYS_DISPERSION_ONLY_VAL) then
do i_SLS = 1,N_SLS
R_xx_val = R_memory(1,i_SLS,i,j,k,ispec)
R_yy_val = R_memory(2,i_SLS,i,j,k,ispec)
@@ -655,7 +655,7 @@
! therefore Q_\alpha is not zero; for instance for V_p / V_s = sqrt(3)
! we get Q_\alpha = (9 / 4) * Q_\mu = 2.25 * Q_\mu
- if(ATTENUATION_VAL .and. ( PARTIAL_PHYS_DISPERSION_ONLY .eqv. .false. )) then
+ if(ATTENUATION_VAL .and. .not. PARTIAL_PHYS_DISPERSION_ONLY_VAL) then
call compute_element_strain_att_noDev(ispec,NGLOB_INNER_CORE,NSPEC_INNER_CORE,displ_inner_core,&
veloc_inner_core,deltat,hprime_xx,hprime_yy,hprime_zz,ibool,&
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_outer_core_noDev.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_outer_core_noDev.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_outer_core_noDev.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -75,7 +75,7 @@
real(kind=CUSTOM_REAL), dimension(NGLLY,NGLLZ) :: wgllwgll_yz
double precision, dimension(NGLLX,NGLLY,NGLLZ) :: wgll_cube
- logical MOVIE_VOLUME
+ logical :: MOVIE_VOLUME
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: tempx1,tempx2,tempx3
@@ -351,7 +351,7 @@
! note: these calculations are only considered for SIMULATION_TYPE == 1 .and. SAVE_FORWARD
! and one has set MOVIE_VOLUME_TYPE == 4 when MOVIE_VOLUME is .true.;
! in case of SIMULATION_TYPE == 3, it gets overwritten by compute_kernels_outer_core()
- if (NSPEC_OUTER_CORE_ADJOINT /= 1 .and. MOVIE_VOLUME ) then
+ if (NSPEC_OUTER_CORE_ADJOINT /= 1 .and. MOVIE_VOLUME) then
div_displfluid(i,j,k,ispec) = &
minus_rho_g_over_kappa_fluid(int_radius) * (dpotentialdx_with_rot * gxl + &
dpotentialdy_with_rot * gyl + dpotentialdzl * gzl)
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_kernels.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_kernels.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_kernels.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -184,7 +184,6 @@
end subroutine compute_kernels_crust_mantle
-
!
!-------------------------------------------------------------------------------------------------
!
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_seismograms.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_seismograms.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_seismograms.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -218,8 +218,6 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_STRAIN_ONLY) :: &
eps_trace_over_3_crust_mantle
-! real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_STR_OR_ATT) :: &
-! epsilondev_crust_mantle
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_STR_OR_ATT) :: &
epsilondev_xx_crust_mantle,epsilondev_yy_crust_mantle,epsilondev_xy_crust_mantle, &
epsilondev_xz_crust_mantle,epsilondev_yz_crust_mantle
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_crust_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_crust_mantle.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_crust_mantle.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -64,7 +64,7 @@
!integer :: reclen1,reclen2
! note: we use c functions for I/O as they still have a better performance than
- ! fortran, unformatted file I/O. however, using -assume byterecl together with fortran functions
+ ! Fortran, unformatted file I/O. however, using -assume byterecl together with Fortran functions
! comes very close (only ~ 4 % slower ).
!
! tests with intermediate storages (every 8 step) and/or asynchronious
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_outer_core.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_outer_core.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_stacey_outer_core.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -62,7 +62,7 @@
integer :: i,j,k,ispec2D,ispec,iglob
! note: we use c functions for I/O as they still have a better performance than
- ! fortran, unformatted file I/O. however, using -assume byterecl together with fortran functions
+ ! Fortran, unformatted file I/O. however, using -assume byterecl together with Fortran functions
! comes very close (only ~ 4 % slower ).
!
! tests with intermediate storages (every 8 step) and/or asynchronious
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/create_central_cube_buffers.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/create_central_cube_buffers.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/create_central_cube_buffers.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -40,11 +40,10 @@
receiver_cube_from_slices,sender_from_slices_to_cube,ibool_central_cube, &
buffer_slices,buffer_slices2,buffer_all_cube_from_slices)
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
integer, intent(in) :: myrank,iproc_xi,iproc_eta,ichunk, &
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -45,7 +45,9 @@
!! DK DK to Daniel, Jul 2013
!! DK DK to Daniel, Jul 2013
!! DK DK to Daniel, Jul 2013: BEWARE, declared real(kind=CUSTOM_REAL) in trunk and
-!! DK DK to Daniel, Jul 2013: double precision in branch, let us check which one is right
+!! DK DK to Daniel, Jul 2013: double precision in branch.
+!! DK DK to Daniel, Jul 2013 real custom is better, it works fine in the trunk and these arrays are really huge
+!! DK DK to Daniel, Jul 2013 in the crust_mantle region, thus let us not double their size
!! DK DK to Daniel, Jul 2013
!! DK DK to Daniel, Jul 2013
!! DK DK to Daniel, Jul 2013
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_backazimuth.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_backazimuth.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_backazimuth.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -30,6 +30,8 @@
implicit none
+ include "constants.h"
+
double precision the, phe
double precision ths, phs
double precision az,baz,xdeg
@@ -40,24 +42,22 @@
double precision phsrad, sc, sd, ss
double precision temp, therad, thg, thsrad
- double precision, parameter :: rad = 6378.160
- double precision, parameter :: fl = 0.00335293
- double precision, parameter :: twopideg = 360.
- double precision, parameter :: c00 = 1.
- double precision, parameter :: c01 = 0.25
- double precision, parameter :: c02 = -4.6875e-02
- double precision, parameter :: c03 = 1.953125e-02
- double precision, parameter :: c21 = -0.125
- double precision, parameter :: c22 = 3.125e-02
- double precision, parameter :: c23 = -1.46484375e-02
- double precision, parameter :: c42 = -3.90625e-03
- double precision, parameter :: c43 = 2.9296875e-03
- double precision, parameter :: degtokm = 111.3199
- double precision, parameter :: pi = 3.141592654
- double precision, parameter :: TORAD = pi/180.
- double precision, parameter :: TODEG = 1./TORAD
+ double precision, parameter :: rad = 6378.160d0
+ double precision, parameter :: fl = 0.00335293d0
+ double precision, parameter :: twopideg = 360.d0
+ double precision, parameter :: c00 = 1.d0
+ double precision, parameter :: c01 = 0.25d0
+ double precision, parameter :: c02 = -4.6875d-02
+ double precision, parameter :: c03 = 1.953125d-02
+ double precision, parameter :: c21 = -0.125d0
+ double precision, parameter :: c22 = 3.125d-02
+ double precision, parameter :: c23 = -1.46484375d-02
+ double precision, parameter :: c42 = -3.90625d-03
+ double precision, parameter :: c43 = 2.9296875d-03
+ double precision, parameter :: degtokm = 111.3199d0
+ double precision, parameter :: TORAD = DEGREES_TO_RADIANS
+ double precision, parameter :: TODEG = RADIANS_TO_DEGREES
-
!=====================================================================
! PURPOSE: To compute the distance and azimuth between locations.
!=====================================================================
@@ -104,7 +104,7 @@
! (Equations are unstable for latidudes of exactly 0 degrees.)
temp = the
- if( temp == 0. ) temp = 1.0e-08
+ if( temp == 0. ) temp = 1.0d-08
therad = TORAD*temp
pherad = TORAD*phe
@@ -130,7 +130,7 @@
! -- Convert to radians.
temp = Ths
- if( temp == 0. ) temp = 1.0e-08
+ if( temp == 0. ) temp = 1.0d-08
thsrad = TORAD*temp
phsrad = TORAD*Phs
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_event_info.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_event_info.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_event_info.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -38,11 +38,10 @@
elat,elon,depth,mb,cmt_lat, &
cmt_lon,cmt_depth,cmt_hdur,NSOURCES)
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
!--- input or output arguments of the subroutine below
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/initialize_simulation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/initialize_simulation.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/initialize_simulation.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -41,6 +41,7 @@
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
+
integer :: ratio_divide_central_cube
integer :: sizeprocs
integer :: ios
@@ -269,12 +270,12 @@
else
write(IMAIN,*) ' no transverse isotropy'
endif
- if(ANISOTROPIC_INNER_CORE) then
+ if(ANISOTROPIC_INNER_CORE_VAL) then
write(IMAIN,*) ' incorporating anisotropic inner core'
else
write(IMAIN,*) ' no inner-core anisotropy'
endif
- if(ANISOTROPIC_3D_MANTLE) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
write(IMAIN,*) ' incorporating anisotropic mantle'
else
write(IMAIN,*) ' no general mantle anisotropy'
@@ -445,16 +446,13 @@
call exit_MPI(myrank, 'SIMULATION_TYPE can only be 1, 2, or 3')
if (SIMULATION_TYPE /= 1 .and. NSOURCES > 999999) &
- call exit_MPI(myrank, &
- 'for adjoint simulations, NSOURCES <= 999999, if you need more change i6.6 in write_seismograms.f90')
+ call exit_MPI(myrank,'for adjoint simulations, NSOURCES <= 999999, if you need more change i6.6 in write_seismograms.f90')
if((SIMULATION_TYPE == 1 .and. SAVE_FORWARD) .or. SIMULATION_TYPE == 3) then
if ( ATTENUATION_VAL) then
! checks mimic flag:
- ! attenuation for adjoint simulations must have USE_ATTENUATION_MIMIC set by xcreate_header_file
-
-!daniel: att - debug todo check attenuation mimick
- if( USE_ATTENUATION_MIMIC .eqv. .false. ) &
+ ! attenuation for adjoint simulations must have PARTIAL_PHYS_DISPERSION_ONLY set by xcreate_header_file
+ if(.not. PARTIAL_PHYS_DISPERSION_ONLY) &
call exit_MPI(myrank,'error in compiled attenuation parameters, please recompile solver 17b')
! user output
@@ -496,13 +494,13 @@
call exit_MPI(myrank, 'anisotropic model is not implemented for kernel simulations yet')
! checks model for transverse isotropic kernel computation
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
if( ANISOTROPIC_3D_MANTLE_VAL ) then
- call exit_mpi(myrank,'error SAVE_TRANSVERSE_KL: Earth model not supported yet')
+ call exit_mpi(myrank,'error SAVE_TRANSVERSE_KL_ONLY: Earth model not supported yet')
endif
if( SIMULATION_TYPE == 3 ) then
if( .not. ANISOTROPIC_KL ) then
- call exit_mpi(myrank,'error SAVE_TRANSVERSE_KL: needs anisotropic kernel calculations')
+ call exit_mpi(myrank,'error SAVE_TRANSVERSE_KL_ONLY: needs anisotropic kernel calculations')
endif
endif
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/iterate_time.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/iterate_time.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/iterate_time.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -371,9 +371,9 @@
! note: for kernel simulations (SIMULATION_TYPE == 3), attenuation is
! only mimicking effects on phase shifts, but not on amplitudes.
- ! flag USE_ATTENUATION_MIMIC will have to be set to true in this case.
+ ! flag PARTIAL_PHYS_DISPERSION_ONLY will have to be set to true in this case.
!
- ! arrays b_R_xx, ... are not used when USE_ATTENUATION_MIMIC is set,
+ ! arrays b_R_xx, ... are not used when PARTIAL_PHYS_DISPERSION_ONLY is set,
! therefore no need to transfer arrays from GPU to CPU
!if (ATTENUATION) then
!endif
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_receivers.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_receivers.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_receivers.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -34,6 +34,7 @@
yr,jda,ho,mi,sec, &
theta_source,phi_source,NCHUNKS,ELLIPTICITY)
+ use mpi
use constants_solver
use specfem_par,only: &
myrank,DT,NSTEP, &
@@ -46,8 +47,6 @@
implicit none
- ! standard include of the MPI library
- include 'mpif.h'
include "precision.h"
integer nspec,nglob
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_regular_points.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_regular_points.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_regular_points.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -153,17 +153,19 @@
!==============================================================
! how about using single precision for the iterations?
-subroutine locate_reg_points(npoints_slice,points_slice,GRID, &
+subroutine locate_regular_points(npoints_slice,points_slice,GRID, &
NEX_XI,nspec,xstore,ystore,zstore,ibool, &
xigll,yigll,zigll,ispec_reg, &
hxir_reg,hetar_reg,hgammar_reg)
implicit none
+
include 'constants.h'
+ include "OUTPUT_FILES/values_from_mesher.h"
! declarations of regular grid model
integer, intent(in) :: npoints_slice
- integer, dimension(NM_KL_REG_PTS), intent(in) :: points_slice
+ integer, dimension(NM_KL_REG_PTS_VAL), intent(in) :: points_slice
type kl_reg_grid_variables
sequence
@@ -190,10 +192,10 @@
double precision, dimension(NGLLZ), intent(in) :: zigll
! output
- integer, dimension(NM_KL_REG_PTS), intent(out) :: ispec_reg
- real(kind=CUSTOM_REAL), dimension(NGLLX,NM_KL_REG_PTS), intent(out) :: hxir_reg
- real(kind=CUSTOM_REAL), dimension(NGLLY,NM_KL_REG_PTS), intent(out) :: hetar_reg
- real(kind=CUSTOM_REAL), dimension(NGLLZ,NM_KL_REG_PTS), intent(out) :: hgammar_reg
+ integer, dimension(NM_KL_REG_PTS_VAL), intent(out) :: ispec_reg
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NM_KL_REG_PTS_VAL), intent(out) :: hxir_reg
+ real(kind=CUSTOM_REAL), dimension(NGLLY,NM_KL_REG_PTS_VAL), intent(out) :: hetar_reg
+ real(kind=CUSTOM_REAL), dimension(NGLLZ,NM_KL_REG_PTS_VAL), intent(out) :: hgammar_reg
! GLL number of anchors
integer, dimension(NGNOD) :: iaddx, iaddy, iaddr
@@ -354,7 +356,7 @@
! DEBUG
! print *, 'Maximum distance discrepancy ', maxval(dist_final(1:npoints_slice))
-end subroutine locate_reg_points
+end subroutine locate_regular_points
!==============================================================
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_sources.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_sources.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/locate_sources.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -33,6 +33,7 @@
xstore,ystore,zstore, &
ELLIPTICITY,min_tshift_cmt_original)
+ use mpi
use constants_solver
use specfem_par,only: &
NSOURCES,myrank, &
@@ -48,9 +49,6 @@
implicit none
- ! standard include of the MPI library
- include 'mpif.h'
-
include "precision.h"
integer nspec,nglob
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/noise_tomography.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/noise_tomography.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/noise_tomography.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -75,12 +75,12 @@
subroutine read_parameters_noise()
+ use mpi
use specfem_par
use specfem_par_crustmantle
use specfem_par_movie
implicit none
- include 'mpif.h'
include "precision.h"
! local parameters
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -183,7 +183,7 @@
if(ATTENUATION_VAL) then
write(IMAIN,*) 'incorporating attenuation using ',N_SLS,' standard linear solids'
if(ATTENUATION_3D_VAL) write(IMAIN,*) 'using 3D attenuation model'
- if(USE_ATTENUATION_MIMIC ) write(IMAIN,*) 'mimicking effects on velocity only'
+ if(PARTIAL_PHYS_DISPERSION_ONLY ) write(IMAIN,*) 'mimicking effects on velocity only'
else
write(IMAIN,*) 'no attenuation'
endif
@@ -619,7 +619,6 @@
! and that we can neglect the 3D model and use PREM every 100 m in all cases
! this is probably a rather reasonable assumption
if(GRAVITY_VAL) then
-
call make_gravity(nspl_gravity,rspl_gravity,gspl,gspl2,ONE_CRUST)
do int_radius = 1,NRAD_GRAVITY
radius = dble(int_radius) / (R_EARTH_KM * 10.d0)
@@ -1490,7 +1489,7 @@
GRAVITY_VAL, &
ROTATION_VAL, &
ATTENUATION_VAL,ATTENUATION_NEW_VAL, &
- USE_ATTENUATION_MIMIC,USE_3D_ATTENUATION_ARRAYS, &
+ PARTIAL_PHYS_DISPERSION_ONLY,USE_3D_ATTENUATION_ARRAYS, &
COMPUTE_AND_STORE_STRAIN_VAL, &
ANISOTROPIC_3D_MANTLE_VAL,ANISOTROPIC_INNER_CORE_VAL, &
SAVE_BOUNDARY_MESH, &
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_buffers_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_buffers_solver.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_buffers_solver.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -36,11 +36,10 @@
NUMMSGS_FACES,NCORNERSCHUNKS,NPROCTOT,NPROC_XI,NPROC_ETA, &
LOCAL_PATH,NCHUNKS)
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
integer iregion_code,myrank,NCHUNKS
@@ -158,11 +157,10 @@
npoin2D_xi,npoin2D_eta, &
NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,LOCAL_PATH)
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
integer iregion_code,myrank
@@ -313,11 +311,10 @@
NGLOB2DMAX_XY,NGLOB1D_RADIAL, &
NUMMSGS_FACES,NCORNERSCHUNKS,NPROC_XI,NPROC_ETA,LOCAL_PATH)
+ use mpi
+
implicit none
-! standard include of the MPI library
- include 'mpif.h'
-
include "constants.h"
integer iregion_code,myrank
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_solver.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_arrays_solver.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -197,7 +197,7 @@
! mass matrices
!
- ! in the case of stacey boundary conditions, add C*deltat/2 contribution to the mass matrix
+ ! in the case of Stacey boundary conditions, add C*deltat/2 contribution to the mass matrix
! on the Stacey edges for the crust_mantle and outer_core regions but not for the inner_core region
! thus the mass matrix must be replaced by three mass matrices including the "C" damping matrix
!
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_mesh_databases.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_mesh_databases.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_mesh_databases.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -439,6 +439,7 @@
! to couple mantle with outer core
+ use mpi
use specfem_par
use specfem_par_crustmantle
use specfem_par_innercore
@@ -446,8 +447,6 @@
implicit none
- include 'mpif.h'
-
! local parameters
integer :: njunk1,njunk2,njunk3
integer :: ier
@@ -605,6 +604,7 @@
subroutine read_mesh_databases_addressing()
+ use mpi
use specfem_par
use specfem_par_crustmantle
use specfem_par_innercore
@@ -612,8 +612,6 @@
implicit none
- include 'mpif.h'
-
! local parameters
integer, dimension(NCHUNKS_VAL,0:NPROC_XI_VAL-1,0:NPROC_ETA_VAL-1) :: addressing
integer, dimension(0:NPROCTOT_VAL-1) :: ichunk_slice,iproc_xi_slice,iproc_eta_slice
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_kernels.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_kernels.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_kernels.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -64,7 +64,7 @@
scale_kl_rho = scale_t / scale_displ / RHOAV * 1.d9
! allocates temporary arrays
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
! transverse isotropic kernel arrays for file output
allocate(alphav_kl_crust_mantle(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT), &
alphah_kl_crust_mantle(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT), &
@@ -106,7 +106,7 @@
rho_kl_crust_mantle(i,j,k,ispec) = rho_kl_crust_mantle(i,j,k,ispec) * scale_kl_rho
! transverse isotropic kernel calculations
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
! note: transverse isotropic kernels are calculated for all elements
!
! however, the factors A,C,L,N,F are based only on transverse elements
@@ -231,7 +231,7 @@
+ N*an_kl(3) + L*an_kl(4) + F*an_kl(5) &
+ rhonotprime_kl_crust_mantle(i,j,k,ispec)
- ! write the kernel in physical units (01/05/2006)
+ ! write the kernel in physical units
rhonotprime_kl_crust_mantle(i,j,k,ispec) = - rhonotprime_kl_crust_mantle(i,j,k,ispec) * scale_kl
alphav_kl_crust_mantle(i,j,k,ispec) = - alphav_kl_crust_mantle(i,j,k,ispec) * scale_kl
@@ -253,11 +253,11 @@
bulk_sq / alphav_sq * alphav_kl_crust_mantle(i,j,k,ispec) + &
bulk_sq / alphah_sq * alphah_kl_crust_mantle(i,j,k,ispec)
- bulk_betah_kl_crust_mantle(i,j,k,ispec ) = &
+ bulk_betah_kl_crust_mantle(i,j,k,ispec) = &
betah_kl_crust_mantle(i,j,k,ispec) + &
FOUR_THIRDS * betah_sq / alphah_sq * alphah_kl_crust_mantle(i,j,k,ispec)
- bulk_betav_kl_crust_mantle(i,j,k,ispec ) = &
+ bulk_betav_kl_crust_mantle(i,j,k,ispec) = &
betav_kl_crust_mantle(i,j,k,ispec) + &
FOUR_THIRDS * betav_sq / alphav_sq * alphav_kl_crust_mantle(i,j,k,ispec)
! the rest, K_eta and K_rho are the same as above
@@ -270,10 +270,10 @@
!rho_kl_crust_mantle(i,j,k,ispec) = rhonotprime_kl_crust_mantle(i,j,k,ispec) &
! + alpha_kl_crust_mantle(i,j,k,ispec) &
! + beta_kl_crust_mantle(i,j,k,ispec)
- bulk_beta_kl_crust_mantle(i,j,k,ispec) = bulk_betah_kl_crust_mantle(i,j,k,ispec ) &
- + bulk_betav_kl_crust_mantle(i,j,k,ispec )
+ bulk_beta_kl_crust_mantle(i,j,k,ispec) = bulk_betah_kl_crust_mantle(i,j,k,ispec) &
+ + bulk_betav_kl_crust_mantle(i,j,k,ispec)
- endif ! SAVE_TRANSVERSE_KL
+ endif ! SAVE_TRANSVERSE_KL_ONLY
else
@@ -320,7 +320,7 @@
if (ANISOTROPIC_KL) then
! outputs transverse isotropic kernels only
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
! transverse isotropic kernels
! (alpha_v, alpha_h, beta_v, beta_h, eta, rho ) parameterization
open(unit=27,file=trim(prname)//'alphav_kernel.bin',status='unknown',form='unformatted',action='write')
@@ -418,7 +418,7 @@
endif
! cleans up temporary kernel arrays
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
deallocate(alphav_kl_crust_mantle,alphah_kl_crust_mantle, &
betav_kl_crust_mantle,betah_kl_crust_mantle, &
eta_kl_crust_mantle)
@@ -462,7 +462,6 @@
rho_kl_outer_core(i,j,k,ispec) = (rho_kl + alpha_kl) * scale_kl
alpha_kl_outer_core(i,j,k,ispec) = 2 * alpha_kl * scale_kl
-
!deviatoric kernel check
if( deviatoric_outercore ) then
beta_kl = - 2 * beta_kl_outer_core(i,j,k,ispec) ! not using mul, since it's zero for the fluid
@@ -674,7 +673,7 @@
! scaling factors
scale_kl = scale_t/scale_displ * 1.d9
- ! scales approximate hessian
+ ! scales approximate Hessian
hess_kl_crust_mantle(:,:,:,:) = 2._CUSTOM_REAL * hess_kl_crust_mantle(:,:,:,:) * scale_kl
! stores into file
@@ -715,9 +714,11 @@
implicit none
include "constants.h"
- real(kind=CUSTOM_REAL) :: theta_in,phi_in
- real(kind=CUSTOM_REAL),dimension(21) :: cij_kll,cij_kl
+ real(kind=CUSTOM_REAL), intent(in) :: theta_in,phi_in
+ real(kind=CUSTOM_REAL), dimension(21), intent(in) :: cij_kl
+ real(kind=CUSTOM_REAL), dimension(21), intent(out) :: cij_kll
+
double precision :: theta,phi
double precision :: costheta,sintheta,cosphi,sinphi
double precision :: costhetasq,sinthetasq,cosphisq,sinphisq
@@ -767,8 +768,7 @@
sintwothetasq = sintwotheta * sintwotheta
sintwophisq = sintwophi * sintwophi
-
- cij_kll(1) = 1.d0/16.d0* (cij_kl(16) - cij_kl(16)* costwophi + &
+ cij_kll(1) = ONE_SIXTEENTH* (cij_kl(16) - cij_kl(16)* costwophi + &
16.d0* cosphi*cosphisq* costhetafour* (cij_kl(1)* cosphi + cij_kl(6)* sinphi) + &
2.d0* (cij_kl(15) + cij_kl(17))* sintwophi* sintwothetasq - &
2.d0* (cij_kl(16)* cosfourtheta* sinphisq + &
@@ -786,7 +786,7 @@
cij_kl(9)* sinphisq)* sintwotheta + &
sinphi* (-cij_kl(13) + cij_kl(9)* sinphisq)* sinfourtheta))
- cij_kll(2) = 1.d0/4.d0* (costhetasq* (cij_kl(1) + 3.d0* cij_kl(2) + cij_kl(7) - &
+ cij_kll(2) = ONE_FOURTH* (costhetasq* (cij_kl(1) + 3.d0* cij_kl(2) + cij_kl(7) - &
cij_kl(21) + (-cij_kl(1) + cij_kl(2) - cij_kl(7) + &
cij_kl(21))* cosfourphi + (-cij_kl(6) + cij_kl(11))* sinfourphi) + &
4.d0* (cij_kl(8)* cosphisq - cij_kl(15)* cosphi* sinphi + &
@@ -796,7 +796,7 @@
(cij_kl(5) - cij_kl(18))* cosphi* sinphisq + &
cij_kl(4)* sinphisq*sinphi)* sintwotheta)
- cij_kll(3) = 1.d0/8.d0* (sintwophi* (3.d0* cij_kl(15) - cij_kl(17) + &
+ cij_kll(3) = ONE_EIGHTH* (sintwophi* (3.d0* cij_kl(15) - cij_kl(17) + &
4.d0* (cij_kl(2) + cij_kl(21))* costhetasq* sintwophi* sinthetasq) + &
4.d0* cij_kl(12)* sintwothetasq + 4.d0* cij_kl(1)* cosphifour* sintwothetasq + &
2.d0* cosphi*cosphisq* (8.d0* cij_kl(6)* costhetasq* sinphi* sinthetasq + &
@@ -811,17 +811,17 @@
8.d0* cij_kl(11)* costhetasq* sinphi*sinphisq* sinthetasq + &
(-cij_kl(14) + (cij_kl(10) + cij_kl(18))* sinphisq)*sinfourtheta))
- cij_kll(4) = 1.d0/8.d0* (cosphi* costheta *(5.d0* cij_kl(4) - &
+ cij_kll(4) = ONE_EIGHTH* (cosphi* costheta *(5.d0* cij_kl(4) - &
cij_kl(9) + 4.d0* cij_kl(13) - &
3.d0* cij_kl(20) + (cij_kl(4) + 3.d0* cij_kl(9) - &
4.d0* cij_kl(13) + cij_kl(20))* costwotheta) + &
- 1.d0/2.d0* (cij_kl(4) - cij_kl(9) + &
+ ONE_HALF* (cij_kl(4) - cij_kl(9) + &
cij_kl(20))* costhreephi * (costheta + 3.d0* costhreetheta) - &
costheta* (-cij_kl(5) + 5.d0* cij_kl(10) + &
4.d0* cij_kl(14) - 3.d0* cij_kl(18) + &
(3.d0* cij_kl(5) + cij_kl(10) - &
4.d0* cij_kl(14) + cij_kl(18))* costwotheta)* sinphi - &
- 1.d0/2.d0* (cij_kl(5) - cij_kl(10) - cij_kl(18))* (costheta + &
+ ONE_HALF* (cij_kl(5) - cij_kl(10) - cij_kl(18))* (costheta + &
3.d0* costhreetheta)* sinthreephi + &
4.d0* (cij_kl(6) - cij_kl(11))* cosfourphi* costhetasq* sintheta - &
4.d0* (cij_kl(1) + cij_kl(3) - cij_kl(7) - cij_kl(8) + cij_kl(16) - cij_kl(19) + &
@@ -833,7 +833,7 @@
2.d0* cij_kl(17))* sintheta + &
(cij_kl(6) + cij_kl(11) - 2.d0* (cij_kl(15) + cij_kl(17)))* sinthreetheta))
- cij_kll(5) = 1.d0/4.d0* (2.d0* (cij_kl(4) + &
+ cij_kll(5) = ONE_FOURTH* (2.d0* (cij_kl(4) + &
cij_kl(20))* cosphisq* (costwotheta + cosfourtheta)* sinphi + &
2.d0* cij_kl(9)* (costwotheta + cosfourtheta)* sinphi*sinphisq + &
16.d0* cij_kl(1)* cosphifour* costheta*costhetasq* sintheta + &
@@ -851,10 +851,10 @@
(cij_kl(3) - cij_kl(16) + cij_kl(19))* costwophi + &
(cij_kl(15) + cij_kl(17))* sintwophi)* sinfourtheta)
- cij_kll(6) = 1.d0/2.d0* costheta*costhetasq* ((cij_kl(6) + cij_kl(11))* costwophi + &
+ cij_kll(6) = ONE_HALF* costheta*costhetasq* ((cij_kl(6) + cij_kl(11))* costwophi + &
(cij_kl(6) - cij_kl(11))* cosfourphi + 2.d0* (-cij_kl(1) + cij_kl(7))* sintwophi + &
(-cij_kl(1) + cij_kl(2) - cij_kl(7) + cij_kl(21))* sinfourphi) + &
- 1.d0/4.d0* costhetasq* (-(cij_kl(4) + 3* cij_kl(9) + cij_kl(20))* cosphi - &
+ ONE_FOURTH* costhetasq* (-(cij_kl(4) + 3* cij_kl(9) + cij_kl(20))* cosphi - &
3.d0* (cij_kl(4) - cij_kl(9) + cij_kl(20))* costhreephi + &
(3.d0* cij_kl(5) + cij_kl(10) + cij_kl(18))* sinphi + &
3.d0* (cij_kl(5) - cij_kl(10) - cij_kl(18))* sinthreephi)* sintheta + &
@@ -862,12 +862,12 @@
(-cij_kl(3) + cij_kl(8) + cij_kl(16) - cij_kl(19))* sintwophi)* sinthetasq + &
(-cij_kl(13)* cosphi + cij_kl(14)* sinphi)* sintheta*sinthetasq
- cij_kll(7) = cij_kl(7)* cosphifour - cij_kl(11)* cosphi*cosphisq* sinphi + &
+ cij_kll(7) = cij_kl(7)* cosphifour - cij_kl(11)* cosphi*cosphisq* sinphi + &
(cij_kl(2) + cij_kl(21))* cosphisq* sinphisq - &
cij_kl(6)* cosphi* sinphi*sinphisq + &
cij_kl(1)* sinphifour
- cij_kll(8) = 1.d0/2.d0* (2.d0* costhetasq* sinphi* (-cij_kl(15)* cosphi + &
+ cij_kll(8) = ONE_HALF* (2.d0* costhetasq* sinphi* (-cij_kl(15)* cosphi + &
cij_kl(3)* sinphi) + 2.d0* cij_kl(2)* cosphifour* sinthetasq + &
(2.d0* cij_kl(2)* sinphifour + &
(cij_kl(1) + cij_kl(7) - cij_kl(21))* sintwophisq)* sinthetasq + &
@@ -879,7 +879,7 @@
cosphisq* (2.d0* cij_kl(8)* costhetasq + &
(cij_kl(9) - cij_kl(20))* sinphi* sintwotheta))
- cij_kll(9) = cij_kl(11)* cosphifour* sintheta - sinphi*sinphisq* (cij_kl(5)* costheta + &
+ cij_kll(9) = cij_kl(11)* cosphifour* sintheta - sinphi*sinphisq* (cij_kl(5)* costheta + &
cij_kl(6)* sinphi* sintheta) + cosphisq* sinphi* (-(cij_kl(10) + &
cij_kl(18))* costheta + &
3.d0* (cij_kl(6) - cij_kl(11))* sinphi* sintheta) + &
@@ -888,7 +888,7 @@
cosphi*cosphisq* (cij_kl(9)* costheta - 2.d0* (cij_kl(2) - 2.d0* cij_kl(7) + &
cij_kl(21))* sinphi* sintheta)
- cij_kll(10) = 1.d0/4.d0* (4.d0* costwotheta* (cij_kl(10)* cosphi*cosphisq + &
+ cij_kll(10) = ONE_FOURTH* (4.d0* costwotheta* (cij_kl(10)* cosphi*cosphisq + &
(cij_kl(9) - cij_kl(20))* cosphisq* sinphi + &
(cij_kl(5) - cij_kl(18))* cosphi* sinphisq + &
cij_kl(4)* sinphi*sinphisq) + (cij_kl(1) + 3.d0* cij_kl(2) - &
@@ -898,7 +898,7 @@
2.d0* cij_kl(15)* sintwophi + &
(-cij_kl(6) + cij_kl(11))* sinfourphi)* sintwotheta)
- cij_kll(11) = 1.d0/4.d0* (2.d0* costheta* ((cij_kl(6) + cij_kl(11))* costwophi + &
+ cij_kll(11) = ONE_FOURTH* (2.d0* costheta* ((cij_kl(6) + cij_kl(11))* costwophi + &
(-cij_kl(6) + cij_kl(11))* cosfourphi + &
2.d0* (-cij_kl(1) + cij_kl(7))* sintwophi + &
(cij_kl(1) - cij_kl(2) + cij_kl(7) - cij_kl(21))* sinfourphi) + &
@@ -907,7 +907,7 @@
(3.d0* cij_kl(5) + cij_kl(10) + cij_kl(18))* sinphi + &
(-cij_kl(5) + cij_kl(10) + cij_kl(18))* sinthreephi)* sintheta)
- cij_kll(12) = 1.d0/16.d0* (cij_kl(16) - 2.d0* cij_kl(16)* cosfourtheta* sinphisq + &
+ cij_kll(12) = ONE_SIXTEENTH* (cij_kl(16) - 2.d0* cij_kl(16)* cosfourtheta* sinphisq + &
costwophi* (-cij_kl(16) + 8.d0* costheta* sinthetasq* ((cij_kl(3) - &
cij_kl(8) + cij_kl(19))* costheta + &
(cij_kl(5) - cij_kl(10) - cij_kl(18))* cosphi* sintheta)) + &
@@ -928,7 +928,7 @@
cij_kl(19)* sintwothetasq + cij_kl(13)* sinphi* sinfourtheta - &
cij_kl(9)* sinphi*sinphisq* sinfourtheta))
- cij_kll(13) = 1.d0/8.d0* (cosphi* costheta* (cij_kl(4) + 3.d0* cij_kl(9) + &
+ cij_kll(13) = ONE_EIGHTH* (cosphi* costheta* (cij_kl(4) + 3.d0* cij_kl(9) + &
4.d0* cij_kl(13) + cij_kl(20) - (cij_kl(4) + 3.d0* cij_kl(9) - &
4.d0* cij_kl(13) + cij_kl(20))* costwotheta) + 4.d0* (-cij_kl(1) - &
cij_kl(3) + cij_kl(7) + cij_kl(8) + cij_kl(16) - cij_kl(19) + &
@@ -946,7 +946,7 @@
(cij_kl(6) + cij_kl(11) - 2.d0* (cij_kl(15) + &
cij_kl(17)))* sinthreetheta))
- cij_kll(14) = 1.d0/4.d0* (2.d0* cij_kl(13)* (costwotheta + cosfourtheta)* sinphi + &
+ cij_kll(14) = ONE_FOURTH* (2.d0* cij_kl(13)* (costwotheta + cosfourtheta)* sinphi + &
8.d0* costheta*costhetasq* (-2.d0* cij_kl(12) + cij_kl(8)* sinphisq)* sintheta + &
4.d0* (cij_kl(4) + cij_kl(20))* cosphisq* (1.d0 + &
2.d0* costwotheta)* sinphi* sinthetasq + &
@@ -962,8 +962,8 @@
(cij_kl(3) + cij_kl(16) + cij_kl(19) + (cij_kl(3) - cij_kl(16) + &
cij_kl(19))* costwophi + (cij_kl(15) + cij_kl(17))* sintwophi)* sinfourtheta)
- cij_kll(15) = costwophi* costheta* (-cij_kl(17) + (cij_kl(15) + cij_kl(17))* costhetasq) + &
- 1.d0/16.d0* (-((11.d0* cij_kl(4) + cij_kl(9) + 4.d0* cij_kl(13) - &
+ cij_kll(15) = costwophi* costheta* (-cij_kl(17) + (cij_kl(15) + cij_kl(17))* costhetasq) + &
+ ONE_SIXTEENTH* (-((11.d0* cij_kl(4) + cij_kl(9) + 4.d0* cij_kl(13) - &
5.d0* cij_kl(20))* cosphi + (cij_kl(4) - cij_kl(9) + cij_kl(20))* costhreephi - &
(cij_kl(5) + 11.d0* cij_kl(10) + 4.d0* cij_kl(14) - &
5.d0* cij_kl(18))* sinphi + (-cij_kl(5) + cij_kl(10) + &
@@ -979,7 +979,7 @@
(3.d0* cij_kl(5) + cij_kl(10) - 4.d0* cij_kl(14) + cij_kl(18))* sinphi + &
3.d0* (-cij_kl(5) + cij_kl(10) + cij_kl(18))* sinthreephi)* sinthreetheta)
- cij_kll(16) = 1.d0/4.d0*(cij_kl(1) - cij_kl(2) + cij_kl(7) + cij_kl(16) + &
+ cij_kll(16) = ONE_FOURTH*(cij_kl(1) - cij_kl(2) + cij_kl(7) + cij_kl(16) + &
cij_kl(19) + cij_kl(21) + 2.d0*(cij_kl(16) - cij_kl(19))*costwophi* costhetasq + &
(-cij_kl(1) + cij_kl(2) - cij_kl(7) + cij_kl(16) + &
cij_kl(19) - cij_kl(21))*costwotheta - 2.d0* cij_kl(17)* costhetasq* sintwophi + &
@@ -989,7 +989,7 @@
(-cij_kl(4) + cij_kl(9) + cij_kl(20))* sinphi - &
(cij_kl(4) - cij_kl(9) + cij_kl(20))* sinthreephi)* sintwotheta)
- cij_kll(17) = 1.d0/8.d0* (4.d0* costwophi* costheta* (cij_kl(6) + cij_kl(11) - &
+ cij_kll(17) = ONE_EIGHTH* (4.d0* costwophi* costheta* (cij_kl(6) + cij_kl(11) - &
2.d0* cij_kl(15) - (cij_kl(6) + cij_kl(11) - 2.d0* (cij_kl(15) + &
cij_kl(17)))* costwotheta) - (2.d0* cosphi* (-3.d0* cij_kl(4) +&
cij_kl(9) + 2.d0* cij_kl(13) + cij_kl(20) + (cij_kl(4) - cij_kl(9) + &
@@ -1005,7 +1005,7 @@
(3.d0* cij_kl(5) + cij_kl(10) - 4.d0* cij_kl(14) + cij_kl(18))* sinphi + &
3.d0* (-cij_kl(5) + cij_kl(10) + cij_kl(18))* sinthreephi)* sinthreetheta)
- cij_kll(18) = 1.d0/2.d0* ((cij_kl(5) - cij_kl(10) + cij_kl(18))* cosphi* costwotheta - &
+ cij_kll(18) = ONE_HALF* ((cij_kl(5) - cij_kl(10) + cij_kl(18))* cosphi* costwotheta - &
(cij_kl(5) - cij_kl(10) - cij_kl(18))* costhreephi* costwotheta - &
2.d0* (cij_kl(4) - cij_kl(9) + &
(cij_kl(4) - cij_kl(9) + cij_kl(20))* costwophi)* costwotheta* sinphi + &
@@ -1015,7 +1015,7 @@
cij_kl(17)* sintwophi + &
(-cij_kl(6) + cij_kl(11))* sinfourphi)* sintwotheta)
- cij_kll(19) = 1.d0/4.d0* (cij_kl(16) - cij_kl(16)* costwophi + &
+ cij_kll(19) = ONE_FOURTH* (cij_kl(16) - cij_kl(16)* costwophi + &
(-cij_kl(15) + cij_kl(17))* sintwophi + &
4.d0* cij_kl(12)* sintwothetasq + &
2.d0* (2.d0* cij_kl(1)* cosphifour* sintwothetasq + &
@@ -1029,7 +1029,7 @@
cosphi* (8.d0* cij_kl(11)* costhetasq* sinphi*sinphisq* sinthetasq + &
(-cij_kl(14) + (cij_kl(10) + cij_kl(18))* sinphisq)* sinfourtheta)))
- cij_kll(20) = 1.d0/8.d0* (2.d0* cosphi* costheta* (-3.d0* cij_kl(4) - cij_kl(9) + &
+ cij_kll(20) = ONE_EIGHTH* (2.d0* cosphi* costheta* (-3.d0* cij_kl(4) - cij_kl(9) + &
4.d0* cij_kl(13) + cij_kl(20) + (cij_kl(4) + 3.d0* cij_kl(9) - &
4.d0* cij_kl(13) + cij_kl(20))* costwotheta) + &
(cij_kl(4) - cij_kl(9) + cij_kl(20))* costhreephi* (costheta + &
@@ -1049,7 +1049,7 @@
2.d0* cij_kl(17))* sintheta + &
(cij_kl(6) + cij_kl(11) - 2.d0* (cij_kl(15) + cij_kl(17)))* sinthreetheta))
- cij_kll(21) = 1.d0/4.d0* (cij_kl(1) - cij_kl(2) + cij_kl(7) + cij_kl(16) + &
+ cij_kll(21) = ONE_FOURTH* (cij_kl(1) - cij_kl(2) + cij_kl(7) + cij_kl(16) + &
cij_kl(19) + cij_kl(21) - 2.d0* (cij_kl(1) - cij_kl(2) + cij_kl(7) - &
cij_kl(21))* cosfourphi* costhetasq + &
(cij_kl(1) - cij_kl(2) + cij_kl(7) - cij_kl(16) - cij_kl(19) + &
@@ -1062,3 +1062,4 @@
cij_kl(20))* sinthreephi)* sintwotheta)
end subroutine rotate_kernels_dble
+
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_regular_kernels.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_regular_kernels.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_regular_kernels.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -35,7 +35,6 @@
kappavstore_crust_mantle,ibool_crust_mantle, &
kappahstore_crust_mantle,muhstore_crust_mantle, &
eta_anisostore_crust_mantle,ispec_is_tiso_crust_mantle, &
- ! --idoubling_crust_mantle, &
LOCAL_PATH)
implicit none
@@ -46,14 +45,14 @@
integer myrank
integer, intent(in) :: npoints_slice
- real, dimension(NGLLX, NM_KL_REG_PTS), intent(in) :: hxir_reg
- real, dimension(NGLLY, NM_KL_REG_PTS), intent(in) :: hetar_reg
- real, dimension(NGLLZ, NM_KL_REG_PTS), intent(in) :: hgammar_reg
- integer, dimension(NM_KL_REG_PTS), intent(in) :: ispec_reg
+ real, dimension(NGLLX, NM_KL_REG_PTS_VAL), intent(in) :: hxir_reg
+ real, dimension(NGLLY, NM_KL_REG_PTS_VAL), intent(in) :: hetar_reg
+ real, dimension(NGLLZ, NM_KL_REG_PTS_VAL), intent(in) :: hgammar_reg
+ integer, dimension(NM_KL_REG_PTS_VAL), intent(in) :: ispec_reg
double precision :: scale_t,scale_displ
- real(kind=CUSTOM_REAL), dimension(21,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT) :: &
+ real(kind=CUSTOM_REAL), dimension(21,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT_ANISO_KL) :: &
cijkl_kl_crust_mantle
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ADJOINT) :: &
@@ -68,7 +67,6 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPECMAX_TISO_MANTLE) :: &
kappahstore_crust_mantle,muhstore_crust_mantle,eta_anisostore_crust_mantle
-! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso_crust_mantle
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
@@ -120,7 +118,7 @@
beta_kl_crust_mantle_reg(npoints_slice), &
alpha_kl_crust_mantle_reg(npoints_slice))
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
! transverse isotropic kernel arrays for file output
allocate(alphav_kl_crust_mantle(npoints_slice), &
alphah_kl_crust_mantle(npoints_slice), &
@@ -145,7 +143,7 @@
do ipoint = 1, npoints_slice
ispec = ispec_reg(ipoint)
if (ANISOTROPIC_KL) then
- if ( SAVE_TRANSVERSE_KL ) then
+ if ( SAVE_TRANSVERSE_KL_ONLY ) then
alphav_kl_crust_mantle(ipoint) = 0.0
alphah_kl_crust_mantle(ipoint) = 0.0
betav_kl_crust_mantle(ipoint) = 0.0
@@ -194,7 +192,7 @@
rho_kl = rho_kl_crust_mantle(i,j,k,ispec) * scale_kl_rho * hlagrange
! transverse isotropic kernel calculations
- if( SAVE_TRANSVERSE_KL ) then
+ if( SAVE_TRANSVERSE_KL_ONLY ) then
! note: transverse isotropic kernels are calculated for all elements
!
! however, the factors A,C,L,N,F are based only on transverse elements
@@ -207,9 +205,6 @@
! Get A,C,F,L,N,eta from kappa,mu
! element can have transverse isotropy if between d220 and Moho
- !if( .not. (TRANSVERSE_ISOTROPY_VAL .and. &
- ! (idoubling_crust_mantle(ispec) == IFLAG_80_MOHO .or. &
- ! idoubling_crust_mantle(ispec) == IFLAG_220_80))) then
if( .not. ispec_is_tiso_crust_mantle(ispec) ) then
! layer with no transverse isotropy
@@ -354,7 +349,7 @@
rho_kl_crust_mantle_reg(ipoint) = rho_kl_crust_mantle_reg(ipoint) + rho_kl
- endif ! SAVE_TRANSVERSE_KL
+ endif ! SAVE_TRANSVERSE_KL_ONLY
else
@@ -400,8 +395,8 @@
! do some transforms that are independent of GLL points
if (ANISOTROPIC_KL) then
- if (SAVE_TRANSVERSE_KL) then
- ! write the kernel in physical units (01/05/2006)
+ if (SAVE_TRANSVERSE_KL_ONLY) then
+ ! write the kernel in physical units
rhonotprime_kl_crust_mantle(ipoint) = - rhonotprime_kl_crust_mantle(ipoint) * scale_kl
alphav_kl_crust_mantle(ipoint) = - alphav_kl_crust_mantle(ipoint) * scale_kl
@@ -432,7 +427,7 @@
if (ANISOTROPIC_KL) then
! outputs transverse isotropic kernels only
- if (SAVE_TRANSVERSE_KL) then
+ if (SAVE_TRANSVERSE_KL_ONLY) then
! transverse isotropic kernels
! (alpha_v, alpha_h, beta_v, beta_h, eta, rho ) parameterization
open(unit=27,file=trim(prname)//'alphav_kernel.bin',status='unknown',form='unformatted',action='write')
@@ -529,7 +524,7 @@
endif
! cleans up temporary kernel arrays
- if (SAVE_TRANSVERSE_KL) then
+ if (SAVE_TRANSVERSE_KL_ONLY) then
deallocate(alphav_kl_crust_mantle,alphah_kl_crust_mantle, &
betav_kl_crust_mantle,betah_kl_crust_mantle, &
eta_kl_crust_mantle)
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/setup_sources_receivers.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/setup_sources_receivers.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/setup_sources_receivers.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -443,7 +443,6 @@
end subroutine setup_receivers
-
!
!-------------------------------------------------------------------------------------------------
!
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/specfem3D.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/specfem3D.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/specfem3D.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -29,6 +29,7 @@
program xspecfem3D
+ use mpi
use specfem_par
use specfem_par_crustmantle
use specfem_par_innercore
@@ -37,9 +38,6 @@
implicit none
- ! standard include of the MPI library
- include 'mpif.h'
-
! local parameters
integer :: ier
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_surface.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_surface.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_surface.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -27,13 +27,13 @@
subroutine write_movie_surface()
+ use mpi
use specfem_par
use specfem_par_crustmantle
use specfem_par_movie
implicit none
- include 'mpif.h'
include "precision.h"
! local parameters
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_volume.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_volume.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_movie_volume.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -54,7 +54,6 @@
! create name of database
write(prname,'(a,i6.6,a)') trim(LOCAL_TMP_PATH)//'/'//'proc',myrank,'_'
-
open(unit=IOUT,file=trim(prname)//'movie3D_info.txt', &
status='unknown',iostat=ier)
if( ier /= 0 ) call exit_mpi(myrank,'error opening file movie3D_info.txt')
@@ -474,13 +473,13 @@
if(NDIM /= 3) call exit_MPI(myrank,'write_movie_volume requires NDIM = 3')
! allocates arrays
-!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppress it
+!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppressing it
!! DK DK (it does not appear in the trunk version of the same routine)
-!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppress it
+!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppressing it
!! DK DK (it does not appear in the trunk version of the same routine)
-!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppress it
+!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppressing it
!! DK DK (it does not appear in the trunk version of the same routine)
-!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppress it
+!! DK DK to Daniel: "vector_scaled" is a big array, we could consider suppressing it
!! DK DK (it does not appear in the trunk version of the same routine)
allocate(store_val3d_N(npoints_3dmovie), &
store_val3d_E(npoints_3dmovie), &
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_seismograms.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_seismograms.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/write_seismograms.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -27,12 +27,9 @@
subroutine write_seismograms()
-! computes the maximum of the norm of the displacement
-! in all the slices using an MPI reduction
-! and output timestamp file to check that simulation is running fine
-
use specfem_par
use specfem_par_crustmantle
+
implicit none
! update position in seismograms
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.F90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.F90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -117,7 +117,8 @@
! gravity
double precision, dimension(NRAD_GRAVITY) :: minus_gravity_table,density_table,minus_deriv_gravity_table
-! local parameters
+ ! local parameters
+
! Deville
! manually inline the calls to the Deville et al. (2002) routines
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: &
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/initialize_simulation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/initialize_simulation.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/initialize_simulation.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -421,7 +421,7 @@
call exit_MPI(myrank, 'SIMULATION_TYPE can only be 1, 2, or 3')
if (SIMULATION_TYPE /= 1 .and. NSOURCES > 999999) &
- call exit_MPI(myrank, 'for adjoint simulations, NSOURCES <= 999999, if you need more change i6.6 in write_seismograms.f90')
+ call exit_MPI(myrank,'for adjoint simulations, NSOURCES <= 999999, if you need more change i6.6 in write_seismograms.f90')
if(UNDO_ATTENUATION .and. PARTIAL_PHYS_DISPERSION_ONLY) &
call exit_MPI(myrank,'cannot have both UNDO_ATTENUATION and PARTIAL_PHYS_DISPERSION_ONLY')
@@ -464,10 +464,10 @@
if (ATTENUATION_VAL .or. SIMULATION_TYPE /= 1 .or. SAVE_FORWARD &
.or. (MOVIE_VOLUME .and. SIMULATION_TYPE /= 3)) then
if( COMPUTE_AND_STORE_STRAIN_VAL .neqv. .true. ) &
- call exit_MPI(myrank, 'error in compiled compute_and_store_strain parameter, please recompile solver 19')
+ call exit_MPI(myrank, 'error in compiled COMPUTE_AND_STORE_STRAIN_VAL parameter, please recompile solver 19')
else
if( COMPUTE_AND_STORE_STRAIN_VAL .neqv. .false. ) &
- call exit_MPI(myrank, 'error in compiled compute_and_store_strain parameter, please recompile solver 20')
+ call exit_MPI(myrank, 'error in compiled COMPUTE_AND_STORE_STRAIN_VAL parameter, please recompile solver 20')
endif
if (SIMULATION_TYPE == 3 .and. (ANISOTROPIC_3D_MANTLE_VAL .or. ANISOTROPIC_INNER_CORE_VAL)) &
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_sources.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_sources.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_sources.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -29,7 +29,7 @@
!---- locate_sources finds the correct position of the sources
!----
- subroutine locate_sources(NSOURCES,myrank,nspec,nglob,ibool,&
+ subroutine locate_sources(NSOURCES,myrank,nspec,nglob,ibool, &
xstore,ystore,zstore,xigll,yigll,zigll, &
NPROCTOT,ELLIPTICITY,TOPOGRAPHY, &
sec,tshift_cmt,min_tshift_cmt_original,yr,jda,ho,mi,theta_source,phi_source, &
@@ -700,7 +700,7 @@
write(IMAIN,*) 'position of the source that will be used:'
write(IMAIN,*)
write(IMAIN,*) ' latitude: ',(PI_OVER_TWO-colat_source)*RADIANS_TO_DEGREES
- write(IMAIN,*) ' longitude: ',phi_source(isource)*180.0d0/PI
+ write(IMAIN,*) ' longitude: ',phi_source(isource)*RADIANS_TO_DEGREES
write(IMAIN,*) ' depth: ',(r0-r_found_source)*R_EARTH/1000.0d0,' km'
write(IMAIN,*)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -79,6 +79,7 @@
NIT, ibool_crust_mantle, ibelm_top_crust_mantle, &
xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle, &
irec_master_noise,normal_x_noise,normal_y_noise,normal_z_noise,mask_noise)
+
use mpi
implicit none
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -687,7 +687,8 @@
integer SIMULATION_TYPE
! local parameters
- integer njunk1,njunk2,njunk3,ier
+ integer :: njunk1,njunk2,njunk3
+ integer :: ier
character(len=150) prname
! user output
@@ -702,7 +703,9 @@
! Stacey put back
open(unit=27,file=prname(1:len_trim(prname))//'boundary.bin', &
- status='old',form='unformatted',action='read')
+ status='old',form='unformatted',action='read',iostat=ier)
+ if( ier /= 0 ) call exit_mpi(myrank,'error opening crust_mantle boundary.bin file')
+
read(27) nspec2D_xmin_crust_mantle
read(27) nspec2D_xmax_crust_mantle
read(27) nspec2D_ymin_crust_mantle
@@ -748,7 +751,9 @@
! Stacey put back
open(unit=27,file=prname(1:len_trim(prname))//'boundary.bin', &
- status='old',form='unformatted',action='read')
+ status='old',form='unformatted',action='read',iostat=ier)
+ if( ier /= 0 ) call exit_mpi(myrank,'error opening outer_core boundary.bin file')
+
read(27) nspec2D_xmin_outer_core
read(27) nspec2D_xmax_outer_core
read(27) nspec2D_ymin_outer_core
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_volume.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_volume.f90 2013-07-24 16:01:33 UTC (rev 22663)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_volume.f90 2013-07-24 17:05:49 UTC (rev 22664)
@@ -68,7 +68,9 @@
mask_3dmovie(:,:,:,:)=.false.
! create name of database
- open(unit=IOUT,file=trim(prname)//'movie3D_info.txt',status='unknown')
+ open(unit=IOUT,file=trim(prname)//'movie3D_info.txt', &
+ status='unknown',iostat=ier)
+ if( ier /= 0 ) call exit_mpi(myrank,'error opening file movie3D_info.txt')
!find and count points within given region for storing movie
do ispec = 1,NSPEC_CRUST_MANTLE
More information about the CIG-COMMITS
mailing list