[cig-commits] r12601 - in seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta: DATA OUTPUT_FILES setup src
dkomati1 at geodynamics.org
dkomati1 at geodynamics.org
Sat Aug 9 16:20:10 PDT 2008
Author: dkomati1
Date: 2008-08-09 16:20:10 -0700 (Sat, 09 Aug 2008)
New Revision: 12601
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_header_file.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/memory_eval.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/model_prem.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/ok_until_here.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90
Log:
fixed many problems to make anisotropic models (both 1D and 3D) work.
Several parameters still do not work, for instance TOPOGRAPHY, GRAVITY, ROTATION
and OCEANS.
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file 2008-08-09 23:20:10 UTC (rev 12601)
@@ -4,7 +4,7 @@
SAVE_FORWARD = .false. # save last frame of forward simulation or not
# number of chunks (1,2,3 or 6)
-NCHUNKS = 1
+NCHUNKS = 6
# angular width of the first chunk (not used if full sphere with six chunks)
ANGULAR_WIDTH_XI_IN_DEGREES = 90.d0 # angular size of a chunk
@@ -31,8 +31,7 @@
# fully 3D models:
# transversely_isotropic_prem_plus_3D_crust_2.0, 3D_anisotropic, 3D_attenuation,
# s20rts, s362ani, s362iso, s362wmani, s362ani_prem, s29ea, s29ea,sea99_jp3d1994,sea99,jp3d1994
-#MODEL = 1D_transversely_isotropic_prem_onecrust
-MODEL = 1D_isotropic_prem_onecrust
+MODEL = 1D_transversely_isotropic_prem_onecrust
# parameters describing the Earth model
OCEANS = .false.
@@ -80,7 +79,7 @@
NUMBER_OF_THIS_RUN = 1
# interval at which we output time step info and max of norm of displacement
-NTSTEP_BETWEEN_OUTPUT_INFO = 1000
+NTSTEP_BETWEEN_OUTPUT_INFO = 100
# interval in time steps for temporary writing of seismograms
NTSTEP_BETWEEN_OUTPUT_SEISMOS = 5000000
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h 2008-08-09 23:20:10 UTC (rev 12601)
@@ -75,7 +75,7 @@
! approximate static memory needed by the solver:
! ----------------------------------------------
!
- ! size of static arrays per slice = 7.239703834056854E-002 GB
+ ! size of static arrays per slice = 7.561429217457771E-002 GB
!
! (should be below and typically equal to 80% of 1.5 GB = 1.2 GB on pangu
! at Caltech, and below and typically equal to 85% of 2 GB = 1.7 GB
@@ -83,8 +83,8 @@
! (if significantly more, the job will not run by lack of memory)
! (if significantly less, you waste a significant amount of memory)
!
- ! size of static arrays for all slices = 0.289588153362274 GB
- ! = 2.828009310178459E-004 TB
+ ! size of static arrays for all slices = 0.302457168698311 GB
+ ! = 2.953683288069442E-004 TB
!
integer, parameter :: NEX_XI_VAL = 64
@@ -101,7 +101,7 @@
integer, parameter :: NSPECMAX_ANISO_IC = 1
integer, parameter :: NSPECMAX_ISO_MANTLE = 7616
- integer, parameter :: NSPECMAX_TISO_MANTLE = 1
+ integer, parameter :: NSPECMAX_TISO_MANTLE = 2304
integer, parameter :: NSPECMAX_ANISO_MANTLE = 1
integer, parameter :: NSPEC_CRUST_MANTLE_ATTENUAT = 1
@@ -129,7 +129,7 @@
integer, parameter :: NGLOB_CRUST_MANTLE_OCEANS = 1
- logical, parameter :: TRANSVERSE_ISOTROPY_VAL = .false.
+ logical, parameter :: TRANSVERSE_ISOTROPY_VAL = .true.
logical, parameter :: ANISOTROPIC_3D_MANTLE_VAL = .false.
@@ -139,7 +139,7 @@
logical, parameter :: ATTENUATION_3D_VAL = .false.
- logical, parameter :: ELLIPTICITY_VAL = .false.
+ logical, parameter :: ELLIPTICITY_VAL = .true.
logical, parameter :: GRAVITY_VAL = .false.
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/setup/constants.h 2008-08-09 23:20:10 UTC (rev 12601)
@@ -56,9 +56,9 @@
! local file unit for output of buffers
integer, parameter :: IOUT_BUFFERS = 35
! uncomment this to write messages to a text file
-! integer, parameter :: IMAIN = 42
+ integer, parameter :: IMAIN = 42
! uncomment this to write messages to the screen (slows down the code)
- integer, parameter :: IMAIN = ISTANDARD_OUTPUT
+! integer, parameter :: IMAIN = ISTANDARD_OUTPUT
! R_EARTH is the radius of the bottom of the oceans (radius of Earth in m)
double precision, parameter :: R_EARTH = 6371000.d0
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_header_file.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_header_file.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_header_file.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -42,7 +42,8 @@
NPROC_XI,NPROC_ETA,NTSTEP_BETWEEN_OUTPUT_SEISMOS, &
NTSTEP_BETWEEN_READ_ADJSRC,NSTEP,NSOURCES,NTSTEP_BETWEEN_FRAMES, &
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,SIMULATION_TYPE, &
- REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP
+ REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP, &
+ ifirst_layer_aniso,ilast_layer_aniso
double precision DT,ANGULAR_WIDTH_XI_IN_DEGREES,ANGULAR_WIDTH_ETA_IN_DEGREES,CENTER_LONGITUDE_IN_DEGREES, &
CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH,ROCEAN,RMIDDLE_CRUST, &
@@ -66,7 +67,7 @@
integer NPROC,NPROCTOT,NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_divide_central_cube
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: ner,ratio_sampling_array
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: doubling_index
- logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_region_has_a_doubling
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_layer_has_a_doubling
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: r_bottom,r_top
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: rmins,rmaxs
@@ -136,11 +137,12 @@
NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX,NSPEC2D_BOTTOM,NSPEC2D_TOP, &
NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,NGLOB, &
- ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_region_has_a_doubling,rmins,rmaxs,CASE_3D, &
+ ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_layer_has_a_doubling,rmins,rmaxs,CASE_3D, &
OUTPUT_SEISMOS_ASCII_TEXT,OUTPUT_SEISMOS_SAC_ALPHANUM,OUTPUT_SEISMOS_SAC_BINARY, &
ROTATE_SEISMOGRAMS_RT,ratio_divide_central_cube,HONOR_1D_SPHERICAL_MOHO,CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,&
DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_XI,DIFF_NSPEC2D_ETA,&
- WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE,USE_BINARY_FOR_LARGE_FILE,.false.)
+ WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE, &
+ USE_BINARY_FOR_LARGE_FILE,ifirst_layer_aniso,ilast_layer_aniso,.false.)
! count the total number of sources in the CMTSOLUTION file
call count_number_of_sources(NSOURCES)
@@ -157,7 +159,7 @@
! evaluate the amount of static memory needed by the solver
call memory_eval(OCEANS,ABSORBING_CONDITIONS,ATTENUATION,ANISOTROPIC_3D_MANTLE,&
TRANSVERSE_ISOTROPY,ANISOTROPIC_INNER_CORE,ROTATION,&
- ONE_CRUST,doubling_index,this_region_has_a_doubling,&
+ ONE_CRUST,doubling_index,this_layer_has_a_doubling,&
ner,NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_sampling_array,&
NSPEC,nglob,SIMULATION_TYPE,MOVIE_VOLUME,SAVE_FORWARD, &
NSPECMAX_ANISO_IC,NSPECMAX_ISO_MANTLE,NSPECMAX_TISO_MANTLE, &
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -38,7 +38,7 @@
ATTENUATION,ATTENUATION_3D, &
NCHUNKS,INCLUDE_CENTRAL_CUBE,ABSORBING_CONDITIONS,REFERENCE_1D_MODEL,THREE_D_MODEL, &
R_CENTRAL_CUBE,RICB,RHO_OCEANS,RCMB,R670,RMOHO,RTOPDDOUBLEPRIME,R600,R220,R771,R400,R120,R80,RMIDDLE_CRUST,ROCEAN, &
- ner,ratio_sampling_array,doubling_index,r_bottom,r_top,this_region_has_a_doubling,CASE_3D, &
+ ner,ratio_sampling_array,doubling_index,r_bottom,r_top,this_layer_has_a_doubling,CASE_3D, &
AMM_V,AM_V,M1066a_V,Mak135_V,Mref_V,SEA1DM_V,D3MM_V,JP3DM_V,SEA99M_V,CM_V, AM_S, AS_V, &
numker,numhpa,numcof,ihpa,lmax,nylm, &
lmxhpa,itypehpa,ihpakern,numcoe,ivarkern, &
@@ -66,7 +66,7 @@
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
- rho_vp,rho_vs, Qmu_store, tau_e_store)
+ rho_vp,rho_vs, Qmu_store, tau_e_store,ifirst_layer_aniso,ilast_layer_aniso)
! create the different regions of the mesh
@@ -123,7 +123,7 @@
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: ner,ratio_sampling_array
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: r_bottom,r_top
- logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_region_has_a_doubling
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_layer_has_a_doubling
integer :: ignod,ner_without_doubling,ispec_superbrick,ilayer,ilayer_loop,ix_elem,iy_elem,iz_elem, &
ifirst_layer,ilast_layer,ratio_divide_central_cube
@@ -537,11 +537,9 @@
logical :: USE_ONE_LAYER_SB,CASE_3D
integer :: nspec_sb
- integer NUMBER_OF_MESH_LAYERS,layer_shift,cpt,first_layer_aniso,last_layer_aniso,FIRST_ELT_NON_ANISO
+ integer NUMBER_OF_MESH_LAYERS,layer_shift,ifirst_layer_aniso,ilast_layer_aniso
double precision, dimension(:,:), allocatable :: stretch_tab
- integer :: nb_layer_above_aniso,FIRST_ELT_ABOVE_ANISO
-
integer, parameter :: maxker=200
integer, parameter :: maxl=72
integer, parameter :: maxcoe=2000
@@ -690,27 +688,31 @@
endif
-! to consider anisotropic elements first and to build the mesh from the bottom to the top of the region
- if (ONE_CRUST) then
- first_layer_aniso=2
- last_layer_aniso=3
- nb_layer_above_aniso = 1
- else
- first_layer_aniso=3
- last_layer_aniso=4
- nb_layer_above_aniso = 2
- endif
- permutation_layer(ifirst_layer:ilast_layer) = (/ (i, i=ilast_layer,ifirst_layer,-1) /)
if(iregion_code == IREGION_CRUST_MANTLE) then
- cpt=3
- permutation_layer(1)=first_layer_aniso
- permutation_layer(2)=last_layer_aniso
- do i = ilast_layer,ifirst_layer,-1
- if (i/=first_layer_aniso .and. i/=last_layer_aniso) then
- permutation_layer(cpt) = i
- cpt=cpt+1
+
+! create anisotropic (transversely isotropic) layers first to save memory when
+! storing the anisotropic arrays
+ ilayer = 0
+ do ilayer_loop = ifirst_layer_aniso,ilast_layer_aniso
+ ilayer = ilayer + 1
+ permutation_layer(ilayer) = ilayer_loop
+ enddo
+
+! and then create all the isotropic layers
+ do ilayer_loop = ifirst_layer,ilast_layer
+ if(ilayer_loop < ifirst_layer_aniso .or. ilayer_loop > ilast_layer_aniso) then
+ ilayer = ilayer + 1
+ permutation_layer(ilayer) = ilayer_loop
endif
enddo
+
+ else
+
+! use identity permutation for regions that do not have transversely isotropic layer
+ do ilayer_loop = ifirst_layer,ilast_layer
+ permutation_layer(ilayer_loop) = ilayer_loop
+ enddo
+
endif
! initialize mesh arrays
@@ -754,19 +756,12 @@
rmin = rmins(ilayer)
rmax = rmaxs(ilayer)
- if(iregion_code == IREGION_CRUST_MANTLE .and. ilayer_loop==3) then
- FIRST_ELT_NON_ANISO = ispec+1
- endif
- if(iregion_code == IREGION_CRUST_MANTLE .and. ilayer_loop==(ilast_layer-nb_layer_above_aniso+1)) then
- FIRST_ELT_ABOVE_ANISO = ispec+1
- endif
-
ner_without_doubling = ner(ilayer)
! if there is a doubling at the top of this region, we implement it in the last two layers of elements
! and therefore we suppress two layers of regular elements here
USE_ONE_LAYER_SB = .false.
- if(this_region_has_a_doubling(ilayer)) then
+ if(this_layer_has_a_doubling(ilayer)) then
if (ner(ilayer) == 1) then
ner_without_doubling = ner_without_doubling - 1
USE_ONE_LAYER_SB = .true.
@@ -890,7 +885,7 @@
! We have imposed that NEX be a multiple of 16 therefore we know that we can always create
! these 2 x 2 blocks because NEX_PER_PROC_XI / ratio_sampling_array(ilayer) and
! NEX_PER_PROC_ETA / ratio_sampling_array(ilayer) are always divisible by 2.
- if(this_region_has_a_doubling(ilayer)) then
+ if(this_layer_has_a_doubling(ilayer)) then
if (USE_ONE_LAYER_SB) then
call define_superbrick_one_layer(x_superbrick,y_superbrick,z_superbrick,ibool_superbrick,iboun_sb)
nspec_sb = NSPEC_SUPERBRICK_1L
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_model.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -860,16 +860,25 @@
!! DK DK changed this for merged version
rhostore_local(i,j,k) = sngl(rho)
-!! DK DK added this for merged version
+! store this for the solid regions: crust_mantle or inner_core
if(iregion_code /= IREGION_OUTER_CORE) then
+
kappavstore(i,j,k,ispec) = sngl(rho*(vpv*vpv - 4.d0*vsv*vsv/3.d0))
- if(iregion_code == IREGION_CRUST_MANTLE .and. NSPECMAX_TISO_MANTLE > 1) &
- kappahstore(i,j,k,ispec) = sngl(rho*(vph*vph - 4.d0*vsh*vsh/3.d0))
muvstore(i,j,k,ispec) = sngl(rho*vsv*vsv)
- if(iregion_code == IREGION_CRUST_MANTLE .and. NSPECMAX_TISO_MANTLE > 1) muhstore(i,j,k,ispec) = sngl(rho*vsh*vsh)
- if(iregion_code == IREGION_CRUST_MANTLE .and. NSPECMAX_TISO_MANTLE > 1) eta_anisostore(i,j,k,ispec) = sngl(eta_aniso)
+
+! store this as well for anisotropic elements in the mantle
+ if(iregion_code == IREGION_CRUST_MANTLE .and. NSPECMAX_TISO_MANTLE > 1) then
+! only store the anisotropic layers, because these arrays have purposely been
+! declared with a reduced size to save memory
+ if(ispec <= NSPECMAX_TISO_MANTLE) then
+ kappahstore(i,j,k,ispec) = sngl(rho*(vph*vph - 4.d0*vsh*vsh/3.d0))
+ muhstore(i,j,k,ispec) = sngl(rho*vsh*vsh)
+ eta_anisostore(i,j,k,ispec) = sngl(eta_aniso)
+ endif
+ endif
+
else
-!! DK DK added this for merged version
+! store only this in the fluid outer core
kappavstore_local(i,j,k) = sngl(rho*(vpv*vpv - 4.d0*vsv*vsv/3.d0))
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/memory_eval.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/memory_eval.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/memory_eval.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -29,7 +29,7 @@
subroutine memory_eval(OCEANS,ABSORBING_CONDITIONS,ATTENUATION,ANISOTROPIC_3D_MANTLE,&
TRANSVERSE_ISOTROPY,ANISOTROPIC_INNER_CORE,ROTATION,&
- ONE_CRUST,doubling_index,this_region_has_a_doubling,&
+ ONE_CRUST,doubling_index,this_layer_has_a_doubling,&
ner,NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_sampling_array,&
NSPEC,nglob,SIMULATION_TYPE,MOVIE_VOLUME,SAVE_FORWARD, &
NSPECMAX_ANISO_IC,NSPECMAX_ISO_MANTLE,NSPECMAX_TISO_MANTLE, &
@@ -55,7 +55,7 @@
integer, dimension(MAX_NUM_REGIONS), intent(in) :: NSPEC, nglob
integer, intent(in) :: NEX_PER_PROC_XI,NEX_PER_PROC_ETA,SIMULATION_TYPE
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS), intent(in) :: doubling_index
- logical, dimension(MAX_NUMBER_OF_MESH_LAYERS), intent(in) :: this_region_has_a_doubling
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS), intent(in) :: this_layer_has_a_doubling
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS), intent(in) :: ner,ratio_sampling_array
! output
@@ -90,7 +90,7 @@
do ilayer = 1, NUMBER_OF_MESH_LAYERS
if(doubling_index(ilayer) == IFLAG_220_80 .or. doubling_index(ilayer) == IFLAG_80_MOHO) then
ner_without_doubling = ner(ilayer)
- if(this_region_has_a_doubling(ilayer)) then
+ if(this_layer_has_a_doubling(ilayer)) then
ner_without_doubling = ner_without_doubling - 2
ispec_aniso = ispec_aniso + &
(NSPEC_DOUBLING_SUPERBRICK*(NEX_PER_PROC_XI/ratio_sampling_array(ilayer)/2)* &
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -288,7 +288,7 @@
! correct number of spectral elements in each block depending on chunk type
- integer nspec_aniso,npointot
+ integer nspec_tiso,npointot
! allocate these automatic arrays in the memory stack to avoid memory fragmentation with "allocate()"
! use the size of the largest region (crust_mantle) and therefore largest possible array
@@ -354,7 +354,8 @@
NPROC_XI,NPROC_ETA,NTSTEP_BETWEEN_OUTPUT_SEISMOS, &
NTSTEP_BETWEEN_READ_ADJSRC,NSTEP,NSOURCES,NTSTEP_BETWEEN_FRAMES, &
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,SIMULATION_TYPE, &
- REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP
+ REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP, &
+ ifirst_layer_aniso,ilast_layer_aniso
double precision DT,ANGULAR_WIDTH_XI_IN_DEGREES,ANGULAR_WIDTH_ETA_IN_DEGREES,CENTER_LONGITUDE_IN_DEGREES, &
CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH,ROCEAN,RMIDDLE_CRUST, &
@@ -394,14 +395,14 @@
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: ner,ratio_sampling_array
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: doubling_index
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: r_bottom,r_top
- logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_region_has_a_doubling
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_layer_has_a_doubling
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: rmins,rmaxs
! memory size of all the static arrays
! double precision :: static_memory_size
! arrays for BCAST
- integer, dimension(38) :: bcast_integer
+ integer, dimension(40) :: bcast_integer
double precision, dimension(30) :: bcast_double_precision
logical, dimension(26) :: bcast_logical
@@ -506,11 +507,12 @@
NSPEC2DMAX_XMIN_XMAX,NSPEC2DMAX_YMIN_YMAX,NSPEC2D_BOTTOM,NSPEC2D_TOP, &
NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX,NGLOB, &
- ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_region_has_a_doubling,rmins,rmaxs,CASE_3D, &
+ ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_layer_has_a_doubling,rmins,rmaxs,CASE_3D, &
OUTPUT_SEISMOS_ASCII_TEXT,OUTPUT_SEISMOS_SAC_ALPHANUM,OUTPUT_SEISMOS_SAC_BINARY, &
ROTATE_SEISMOGRAMS_RT,ratio_divide_central_cube,HONOR_1D_SPHERICAL_MOHO,CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,&
DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_XI,DIFF_NSPEC2D_ETA,&
- WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE,USE_BINARY_FOR_LARGE_FILE,.false.)
+ WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE, &
+ USE_BINARY_FOR_LARGE_FILE,ifirst_layer_aniso,ilast_layer_aniso,.false.)
if(err_occurred() /= 0) call exit_MPI(myrank,'an error occurred while reading the parameter file')
@@ -526,7 +528,7 @@
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,&
SIMULATION_TYPE,REFERENCE_1D_MODEL,THREE_D_MODEL,NPROC,NPROCTOT, &
NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_divide_central_cube,&
- MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP/)
+ MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP,ifirst_layer_aniso,ilast_layer_aniso/)
bcast_logical = (/TRANSVERSE_ISOTROPY,ANISOTROPIC_3D_MANTLE,ANISOTROPIC_INNER_CORE, &
CRUSTAL,ELLIPTICITY,GRAVITY,ONE_CRUST,ROTATION,ISOTROPIC_3D_MANTLE, &
@@ -545,7 +547,7 @@
endif
! broadcast the information read on the master to the nodes
- call MPI_BCAST(bcast_integer,38,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(bcast_integer,40,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(bcast_double_precision,30,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
@@ -563,7 +565,7 @@
call MPI_BCAST(rmaxs,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(rmaxs,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(this_region_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(this_layer_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(NSPEC,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(NSPEC2D_XI,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
@@ -622,6 +624,8 @@
MOVIE_VOLUME_TYPE = bcast_integer(36)
MOVIE_START = bcast_integer(37)
MOVIE_STOP = bcast_integer(38)
+ ifirst_layer_aniso = bcast_integer(39)
+ ilast_layer_aniso = bcast_integer(40)
TRANSVERSE_ISOTROPY = bcast_logical(1)
ANISOTROPIC_3D_MANTLE = bcast_logical(2)
@@ -1181,7 +1185,7 @@
if(iregion_code == IREGION_CRUST_MANTLE) then
! crust_mantle
call create_regions_mesh(iregion_code,ibool_crust_mantle,idoubling_crust_mantle, &
- xstore,ystore,zstore,rmins,rmaxs,iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_aniso, &
+ xstore,ystore,zstore,rmins,rmaxs,iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_tiso, &
volume_local,area_local_bottom,area_local_top,nspl,rspl,espl,espl2,nglob(iregion_code),npointot, &
NEX_XI,NEX_PER_PROC_XI,NEX_PER_PROC_ETA,NSPEC2DMAX_XMIN_XMAX(iregion_code),NSPEC2DMAX_YMIN_YMAX(iregion_code), &
NSPEC2D_BOTTOM(iregion_code),NSPEC2D_TOP(iregion_code),ELLIPTICITY,TOPOGRAPHY,TRANSVERSE_ISOTROPY,ANISOTROPIC_3D_MANTLE, &
@@ -1190,7 +1194,7 @@
myrank,OCEANS,ibathy_topo,rotation_matrix,ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD, &
ATTENUATION,ATTENUATION_3D,NCHUNKS,INCLUDE_CENTRAL_CUBE,ABSORBING_CONDITIONS,REFERENCE_1D_MODEL,THREE_D_MODEL, &
R_CENTRAL_CUBE,RICB,RHO_OCEANS,RCMB,R670,RMOHO,RTOPDDOUBLEPRIME,R600,R220,R771,R400,R120,R80,RMIDDLE_CRUST,ROCEAN, &
- ner,ratio_sampling_array,doubling_index,r_bottom, r_top,this_region_has_a_doubling,CASE_3D, &
+ ner,ratio_sampling_array,doubling_index,r_bottom, r_top,this_layer_has_a_doubling,CASE_3D, &
AMM_V,AM_V,M1066a_V,Mak135_V,Mref_V,SEA1DM_V,D3MM_V,JP3DM_V,SEA99M_V,CM_V, AM_S, AS_V, &
numker,numhpa,numcof,ihpa,lmax,nylm,lmxhpa,itypehpa,ihpakern,numcoe,ivarkern, &
nconpt,iver,iconpt,conpt,xlaspl,xlospl,radspl, &
@@ -1220,12 +1224,12 @@
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
- rho_vp,rho_vs, Qmu_store, tau_e_store)
+ rho_vp,rho_vs, Qmu_store, tau_e_store,ifirst_layer_aniso,ilast_layer_aniso)
else if(iregion_code == IREGION_OUTER_CORE) then
! outer_core
call create_regions_mesh(iregion_code,ibool_outer_core,idoubling_outer_core, &
- xstore,ystore,zstore,rmins,rmaxs,iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_aniso, &
+ xstore,ystore,zstore,rmins,rmaxs,iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_tiso, &
volume_local,area_local_bottom,area_local_top,nspl,rspl,espl,espl2,nglob(iregion_code),npointot, &
NEX_XI,NEX_PER_PROC_XI,NEX_PER_PROC_ETA,NSPEC2DMAX_XMIN_XMAX(iregion_code), &
NSPEC2DMAX_YMIN_YMAX(iregion_code),NSPEC2D_BOTTOM(iregion_code),NSPEC2D_TOP(iregion_code), &
@@ -1235,7 +1239,7 @@
myrank,OCEANS,ibathy_topo,rotation_matrix,ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD, &
ATTENUATION,ATTENUATION_3D,NCHUNKS,INCLUDE_CENTRAL_CUBE,ABSORBING_CONDITIONS,REFERENCE_1D_MODEL,THREE_D_MODEL, &
R_CENTRAL_CUBE,RICB,RHO_OCEANS,RCMB,R670,RMOHO,RTOPDDOUBLEPRIME,R600,R220,R771,R400,R120,R80,RMIDDLE_CRUST,ROCEAN, &
- ner,ratio_sampling_array,doubling_index,r_bottom, r_top,this_region_has_a_doubling,CASE_3D, &
+ ner,ratio_sampling_array,doubling_index,r_bottom, r_top,this_layer_has_a_doubling,CASE_3D, &
AMM_V,AM_V,M1066a_V,Mak135_V,Mref_V,SEA1DM_V,D3MM_V,JP3DM_V,SEA99M_V,CM_V, AM_S, AS_V, &
numker,numhpa,numcof,ihpa,lmax,nylm,lmxhpa,itypehpa,ihpakern,numcoe,ivarkern, &
nconpt,iver,iconpt,conpt,xlaspl,xlospl,radspl, &
@@ -1265,12 +1269,12 @@
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
- rho_vp,rho_vs, Qmu_store, tau_e_store)
+ rho_vp,rho_vs, Qmu_store, tau_e_store,ifirst_layer_aniso,ilast_layer_aniso)
else if(iregion_code == IREGION_INNER_CORE) then
! inner_core
call create_regions_mesh(iregion_code,ibool_inner_core,idoubling_inner_core, &
- xstore,ystore,zstore,rmins,rmaxs,iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_aniso, &
+ xstore,ystore,zstore,rmins,rmaxs,iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_tiso, &
volume_local,area_local_bottom,area_local_top,nspl,rspl,espl,espl2,nglob(iregion_code),npointot, &
NEX_XI,NEX_PER_PROC_XI,NEX_PER_PROC_ETA,NSPEC2DMAX_XMIN_XMAX(iregion_code), &
NSPEC2DMAX_YMIN_YMAX(iregion_code),NSPEC2D_BOTTOM(iregion_code),NSPEC2D_TOP(iregion_code), &
@@ -1280,7 +1284,7 @@
myrank,OCEANS,ibathy_topo,rotation_matrix,ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD, &
ATTENUATION,ATTENUATION_3D,NCHUNKS,INCLUDE_CENTRAL_CUBE,ABSORBING_CONDITIONS,REFERENCE_1D_MODEL,THREE_D_MODEL, &
R_CENTRAL_CUBE,RICB,RHO_OCEANS,RCMB,R670,RMOHO,RTOPDDOUBLEPRIME,R600,R220,R771,R400,R120,R80,RMIDDLE_CRUST,ROCEAN, &
- ner,ratio_sampling_array,doubling_index,r_bottom, r_top,this_region_has_a_doubling,CASE_3D, &
+ ner,ratio_sampling_array,doubling_index,r_bottom, r_top,this_layer_has_a_doubling,CASE_3D, &
AMM_V,AM_V,M1066a_V,Mak135_V,Mref_V,SEA1DM_V,D3MM_V,JP3DM_V,SEA99M_V,CM_V, AM_S, AS_V, &
numker,numhpa,numcof,ihpa,lmax,nylm,lmxhpa,itypehpa,ihpakern,numcoe,ivarkern, &
nconpt,iver,iconpt,conpt,xlaspl,xlospl,radspl, &
@@ -1309,18 +1313,18 @@
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
- rho_vp,rho_vs, Qmu_store, tau_e_store)
+ rho_vp,rho_vs, Qmu_store, tau_e_store,ifirst_layer_aniso,ilast_layer_aniso)
else
stop 'DK DK incorrect region in merged code'
endif
! store number of anisotropic elements found in the mantle
- if(nspec_aniso /= 0 .and. iregion_code /= IREGION_CRUST_MANTLE) &
- call exit_MPI(myrank,'found anisotropic elements outside of the mantle')
+ if(nspec_tiso /= 0 .and. iregion_code /= IREGION_CRUST_MANTLE) &
+ call exit_MPI(myrank,'found transversely isotropic elements outside of the mantle')
- if(iregion_code == IREGION_CRUST_MANTLE .and. nspec_aniso == 0) &
- call exit_MPI(myrank,'found no anisotropic elements in the mantle')
+ if(iregion_code == IREGION_CRUST_MANTLE .and. nspec_tiso == 0) &
+ call exit_MPI(myrank,'found no transversely isotropic elements in the mantle')
! use MPI reduction to compute total area and volume
!! DK DK suppressed for now in the merged version, for simplicity
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/model_prem.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/model_prem.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/model_prem.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -418,7 +418,7 @@
Qkappa=57827.0d0
else if(r > R220 .and. r <= R80) then
-! anisotropy in PREM only above 220 km
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
rho=2.6910d0+0.6924d0*x
vpv=0.8317d0+7.2180d0*x
@@ -446,7 +446,7 @@
! use PREM crust
if(r > R80 .and. r <= RMOHO) then
-! anisotropy in PREM only above 220 km
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
rho=2.6910d0+0.6924d0*x
vpv=0.8317d0+7.2180d0*x
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/ok_until_here.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/ok_until_here.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/ok_until_here.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -1,5 +1,6 @@
print *,'ok until here'
+ call flush(6) ! flush the screen buffer (supported at least in Intel ifort)
call MPI_BARRIER(MPI_COMM_WORLD,ier)
! stop all the MPI processes, and exit
call MPI_ABORT(MPI_COMM_WORLD,30,ier)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/read_compute_parameters.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -52,13 +52,13 @@
NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX, &
NGLOB, &
- ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_region_has_a_doubling,rmins,rmaxs,CASE_3D, &
+ ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_layer_has_a_doubling,rmins,rmaxs,CASE_3D, &
OUTPUT_SEISMOS_ASCII_TEXT,OUTPUT_SEISMOS_SAC_ALPHANUM,OUTPUT_SEISMOS_SAC_BINARY, &
ROTATE_SEISMOGRAMS_RT,ratio_divide_central_cube,HONOR_1D_SPHERICAL_MOHO,CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,&
DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_XI,DIFF_NSPEC2D_ETA,&
- WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE,USE_BINARY_FOR_LARGE_FILE,EMULATE_ONLY)
+ WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE, &
+ USE_BINARY_FOR_LARGE_FILE,ifirst_layer_aniso,ilast_layer_aniso,EMULATE_ONLY)
-
implicit none
include "constants.h"
@@ -72,7 +72,7 @@
NTSTEP_BETWEEN_READ_ADJSRC,NSTEP,NTSTEP_BETWEEN_FRAMES, &
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,SIMULATION_TYPE, &
REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP, &
- NEX_XI_read,NEX_ETA_read,NPROC_XI_read,NPROC_ETA_read
+ NEX_XI_read,NEX_ETA_read,NPROC_XI_read,NPROC_ETA_read,ifirst_layer_aniso,ilast_layer_aniso
double precision DT,ANGULAR_WIDTH_XI_IN_DEGREES,ANGULAR_WIDTH_ETA_IN_DEGREES,CENTER_LONGITUDE_IN_DEGREES, &
CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH,ROCEAN,RMIDDLE_CRUST, &
@@ -115,7 +115,7 @@
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: ner,ratio_sampling_array
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: r_bottom,r_top
- logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_region_has_a_doubling
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_layer_has_a_doubling
integer :: ielem,elem_doubling_mantle,elem_doubling_middle_outer_core,elem_doubling_bottom_outer_core
double precision :: DEPTH_SECOND_DOUBLING_REAL,DEPTH_THIRD_DOUBLING_REAL, &
@@ -1251,6 +1251,10 @@
ner(13) = elem_doubling_middle_outer_core
ner(14) = NER_TOP_CENTRAL_CUBE_ICB
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
+ ifirst_layer_aniso = 1
+ ilast_layer_aniso = 4
+
! value of the doubling ratio in each radial region of the mesh
ratio_sampling_array(1:9) = 1
ratio_sampling_array(10:12) = 2
@@ -1265,9 +1269,9 @@
doubling_index(14) = IFLAG_INNER_CORE_NORMAL
! define the three regions in which we implement a mesh doubling at the top of that region
- this_region_has_a_doubling(:) = .false.
- this_region_has_a_doubling(10) = .true.
- this_region_has_a_doubling(13) = .true.
+ this_layer_has_a_doubling(:) = .false.
+ this_layer_has_a_doubling(10) = .true.
+ this_layer_has_a_doubling(13) = .true.
lastdoubling_layer = 13
! define the top and bottom radii of all the regions of the mesh in the radial direction
@@ -1372,6 +1376,10 @@
ner(12) = elem_doubling_middle_outer_core
ner(13) = NER_TOP_CENTRAL_CUBE_ICB
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
+ ifirst_layer_aniso = 2
+ ilast_layer_aniso = 3
+
! value of the doubling ratio in each radial region of the mesh
ratio_sampling_array(1) = 1
ratio_sampling_array(2:8) = 2
@@ -1388,10 +1396,10 @@
doubling_index(13) = IFLAG_INNER_CORE_NORMAL
! define the three regions in which we implement a mesh doubling at the top of that region
- this_region_has_a_doubling(:) = .false.
- this_region_has_a_doubling(2) = .true.
- this_region_has_a_doubling(9) = .true.
- this_region_has_a_doubling(12) = .true.
+ this_layer_has_a_doubling(:) = .false.
+ this_layer_has_a_doubling(2) = .true.
+ this_layer_has_a_doubling(9) = .true.
+ this_layer_has_a_doubling(12) = .true.
lastdoubling_layer = 12
! define the top and bottom radii of all the regions of the mesh in the radial direction
@@ -1502,6 +1510,10 @@
ner(13) = elem_doubling_middle_outer_core
ner(14) = NER_TOP_CENTRAL_CUBE_ICB
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
+ ifirst_layer_aniso = 3
+ ilast_layer_aniso = 4
+
! value of the doubling ratio in each radial region of the mesh
ratio_sampling_array(1:2) = 1
ratio_sampling_array(3:9) = 2
@@ -1518,11 +1530,11 @@
doubling_index(14) = IFLAG_INNER_CORE_NORMAL
! define the three regions in which we implement a mesh doubling at the top of that region
- this_region_has_a_doubling(:) = .false.
- this_region_has_a_doubling(3) = .true.
- this_region_has_a_doubling(10) = .true.
- this_region_has_a_doubling(13) = .true.
- this_region_has_a_doubling(14) = .false.
+ this_layer_has_a_doubling(:) = .false.
+ this_layer_has_a_doubling(3) = .true.
+ this_layer_has_a_doubling(10) = .true.
+ this_layer_has_a_doubling(13) = .true.
+ this_layer_has_a_doubling(14) = .false.
lastdoubling_layer = 13
! define the top and bottom radii of all the regions of the mesh in the radial direction
@@ -1638,6 +1650,10 @@
ner(14) = elem_doubling_bottom_outer_core
ner(15) = NER_TOP_CENTRAL_CUBE_ICB
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
+ ifirst_layer_aniso = 1
+ ilast_layer_aniso = 4
+
! value of the doubling ratio in each radial region of the mesh
ratio_sampling_array(1:9) = 1
ratio_sampling_array(10:12) = 2
@@ -1653,10 +1669,10 @@
doubling_index(15) = IFLAG_INNER_CORE_NORMAL
! define the three regions in which we implement a mesh doubling at the top of that region
- this_region_has_a_doubling(:) = .false.
- this_region_has_a_doubling(10) = .true.
- this_region_has_a_doubling(13) = .true.
- this_region_has_a_doubling(14) = .true.
+ this_layer_has_a_doubling(:) = .false.
+ this_layer_has_a_doubling(10) = .true.
+ this_layer_has_a_doubling(13) = .true.
+ this_layer_has_a_doubling(14) = .true.
lastdoubling_layer = 14
! define the top and bottom radii of all the regions of the mesh in the radial direction
@@ -1765,6 +1781,10 @@
ner(13) = elem_doubling_bottom_outer_core
ner(14) = NER_TOP_CENTRAL_CUBE_ICB
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
+ ifirst_layer_aniso = 2
+ ilast_layer_aniso = 3
+
! value of the doubling ratio in each radial region of the mesh
ratio_sampling_array(1) = 1
ratio_sampling_array(2:8) = 2
@@ -1782,11 +1802,11 @@
doubling_index(14) = IFLAG_INNER_CORE_NORMAL
! define the three regions in which we implement a mesh doubling at the top of that region
- this_region_has_a_doubling(:) = .false.
- this_region_has_a_doubling(2) = .true.
- this_region_has_a_doubling(9) = .true.
- this_region_has_a_doubling(12) = .true.
- this_region_has_a_doubling(13) = .true.
+ this_layer_has_a_doubling(:) = .false.
+ this_layer_has_a_doubling(2) = .true.
+ this_layer_has_a_doubling(9) = .true.
+ this_layer_has_a_doubling(12) = .true.
+ this_layer_has_a_doubling(13) = .true.
lastdoubling_layer = 13
! define the top and bottom radii of all the regions of the mesh in the radial direction
@@ -1901,6 +1921,10 @@
ner(14) = elem_doubling_bottom_outer_core
ner(15) = NER_TOP_CENTRAL_CUBE_ICB
+! anisotropy in PREM only between 220 km and the Moho (bottom of the crust)
+ ifirst_layer_aniso = 3
+ ilast_layer_aniso = 4
+
! value of the doubling ratio in each radial region of the mesh
ratio_sampling_array(1:2) = 1
ratio_sampling_array(3:9) = 2
@@ -1918,11 +1942,11 @@
doubling_index(15) = IFLAG_INNER_CORE_NORMAL
! define the three regions in which we implement a mesh doubling at the top of that region
- this_region_has_a_doubling(:) = .false.
- this_region_has_a_doubling(3) = .true.
- this_region_has_a_doubling(10) = .true.
- this_region_has_a_doubling(13) = .true.
- this_region_has_a_doubling(14) = .true.
+ this_layer_has_a_doubling(:) = .false.
+ this_layer_has_a_doubling(3) = .true.
+ this_layer_has_a_doubling(10) = .true.
+ this_layer_has_a_doubling(13) = .true.
+ this_layer_has_a_doubling(14) = .true.
lastdoubling_layer = 14
! define the top and bottom radii of all the regions of the mesh in the radial direction
@@ -2107,7 +2131,7 @@
tmp_sum_nglob2D_xi = 0
tmp_sum_nglob2D_eta = 0
do iter_layer = ifirst_region, ilast_region
- if (this_region_has_a_doubling(iter_layer)) then
+ if (this_layer_has_a_doubling(iter_layer)) then
if (iter_region == IREGION_OUTER_CORE .and. iter_layer == lastdoubling_layer) then
! simple brick
divider = 1
@@ -2279,7 +2303,7 @@
endif
tmp_sum = 0;
do iter_layer = ifirst_region, ilast_region
- if (this_region_has_a_doubling(iter_layer)) then
+ if (this_layer_has_a_doubling(iter_layer)) then
if (ner(iter_layer) == 1) then
nb_lay_sb = 1
nspec_sb = NSPEC_SUPERBRICK_1L
@@ -2375,7 +2399,7 @@
nglob_center_edge = 0
nglob_corner_edge = 0
nglob_border_edge = 0
- if (this_region_has_a_doubling(iter_layer)) then
+ if (this_layer_has_a_doubling(iter_layer)) then
if (iter_region == IREGION_OUTER_CORE .and. iter_layer == lastdoubling_layer .and. &
(CUT_SUPERBRICK_XI .or. CUT_SUPERBRICK_ETA)) then
doubling = 1
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90 2008-08-09 00:22:31 UTC (rev 12600)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90 2008-08-09 23:20:10 UTC (rev 12601)
@@ -378,7 +378,8 @@
NPROC_XI,NPROC_ETA,NTSTEP_BETWEEN_OUTPUT_SEISMOS,&
NTSTEP_BETWEEN_READ_ADJSRC,NSTEP,NSOURCES,NTSTEP_BETWEEN_FRAMES, &
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,SIMULATION_TYPE, &
- REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP
+ REFERENCE_1D_MODEL,THREE_D_MODEL,MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP, &
+ ifirst_layer_aniso,ilast_layer_aniso
double precision DT,ANGULAR_WIDTH_XI_IN_DEGREES,ANGULAR_WIDTH_ETA_IN_DEGREES,CENTER_LONGITUDE_IN_DEGREES, &
CENTER_LATITUDE_IN_DEGREES,GAMMA_ROTATION_AZIMUTH,ROCEAN,RMIDDLE_CRUST, &
@@ -434,12 +435,12 @@
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: ner,ratio_sampling_array
integer, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: doubling_index
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: r_bottom,r_top
- logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_region_has_a_doubling
+ logical, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: this_layer_has_a_doubling
double precision, dimension(MAX_NUMBER_OF_MESH_LAYERS) :: rmins,rmaxs
logical :: CASE_3D
! arrays for BCAST
- integer, dimension(38) :: bcast_integer
+ integer, dimension(40) :: bcast_integer
double precision, dimension(30) :: bcast_double_precision
logical, dimension(33) :: bcast_logical
@@ -529,11 +530,12 @@
NSPEC1D_RADIAL,NGLOB1D_RADIAL, &
NGLOB2DMAX_XMIN_XMAX,NGLOB2DMAX_YMIN_YMAX, &
NGLOB_computed, &
- ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_region_has_a_doubling,rmins,rmaxs,CASE_3D, &
+ ratio_sampling_array, ner, doubling_index,r_bottom,r_top,this_layer_has_a_doubling,rmins,rmaxs,CASE_3D, &
OUTPUT_SEISMOS_ASCII_TEXT,OUTPUT_SEISMOS_SAC_ALPHANUM,OUTPUT_SEISMOS_SAC_BINARY, &
ROTATE_SEISMOGRAMS_RT,ratio_divide_central_cube,HONOR_1D_SPHERICAL_MOHO,CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,&
DIFF_NSPEC1D_RADIAL,DIFF_NSPEC2D_XI,DIFF_NSPEC2D_ETA,&
- WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE,USE_BINARY_FOR_LARGE_FILE,.false.)
+ WRITE_SEISMOGRAMS_BY_MASTER,SAVE_ALL_SEISMOS_IN_ONE_FILE, &
+ USE_BINARY_FOR_LARGE_FILE,ifirst_layer_aniso,ilast_layer_aniso,.false.)
if(err_occurred() /= 0) call exit_MPI(myrank,'an error occurred while reading the parameter file')
@@ -546,7 +548,7 @@
NTSTEP_BETWEEN_OUTPUT_INFO,NUMBER_OF_RUNS,NUMBER_OF_THIS_RUN,NCHUNKS,&
SIMULATION_TYPE,REFERENCE_1D_MODEL,THREE_D_MODEL,NPROC,NPROCTOT, &
NEX_PER_PROC_XI,NEX_PER_PROC_ETA,ratio_divide_central_cube,&
- MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP/)
+ MOVIE_VOLUME_TYPE,MOVIE_START,MOVIE_STOP,ifirst_layer_aniso,ilast_layer_aniso/)
bcast_logical = (/TRANSVERSE_ISOTROPY,ANISOTROPIC_3D_MANTLE,ANISOTROPIC_INNER_CORE, &
CRUSTAL,ELLIPTICITY,GRAVITY,ONE_CRUST,ROTATION,ISOTROPIC_3D_MANTLE, &
@@ -567,7 +569,7 @@
endif
! broadcast the information read on the master to the nodes
- call MPI_BCAST(bcast_integer,38,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(bcast_integer,40,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(bcast_double_precision,30,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
@@ -584,7 +586,7 @@
call MPI_BCAST(rmins,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(rmaxs,MAX_NUMBER_OF_MESH_LAYERS,MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ier)
- call MPI_BCAST(this_region_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
+ call MPI_BCAST(this_layer_has_a_doubling,MAX_NUMBER_OF_MESH_LAYERS,MPI_LOGICAL,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(NSPEC_computed,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
call MPI_BCAST(NSPEC2D_XI,MAX_NUM_REGIONS,MPI_INTEGER,0,MPI_COMM_WORLD,ier)
@@ -643,6 +645,8 @@
MOVIE_VOLUME_TYPE = bcast_integer(36)
MOVIE_START = bcast_integer(37)
MOVIE_STOP = bcast_integer(38)
+ ifirst_layer_aniso = bcast_integer(39)
+ ilast_layer_aniso = bcast_integer(40)
TRANSVERSE_ISOTROPY = bcast_logical(1)
ANISOTROPIC_3D_MANTLE = bcast_logical(2)
More information about the cig-commits
mailing list