[cig-commits] r18167 - in seismo/3D/SPECFEM3D_GLOBE/trunk: setup src/meshfem3D src/shared src/specfem3D
danielpeter at geodynamics.org
danielpeter at geodynamics.org
Sun Apr 3 21:19:03 PDT 2011
Author: danielpeter
Date: 2011-04-03 21:19:03 -0700 (Sun, 03 Apr 2011)
New Revision: 18167
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/save_kernels.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
Log:
adds constant USE_FULL_TISO_MANTLE to simulate transverse isotropy in all mantle elements (more expensive, computation takes ~45% longer)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in 2011-04-04 04:19:03 UTC (rev 18167)
@@ -256,8 +256,11 @@
logical, parameter :: APPROXIMATE_HESS_KL = .false.
! output kernel mask to zero out source region
- logical,parameter :: SAVE_SOURCE_MASK = .false.
+ logical, parameter :: SAVE_SOURCE_MASK = .false.
+! forces transverse isotropy for all mantle elements
+! (default is to use transverse isotropy only between MOHO and 220 )
+ logical, parameter :: USE_FULL_TISO_MANTLE = .false.
!!-----------------------------------------------------------
!!
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -39,7 +39,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
- xigll,yigll,zigll)
+ xigll,yigll,zigll,ispec_is_tiso)
use meshfem3D_models_par
@@ -107,15 +107,21 @@
double precision:: yigll(NGLLY)
double precision:: zigll(NGLLZ)
+ logical, dimension(nspec) :: ispec_is_tiso
+
! Parameter used to decide whether this element is in the crust or not
logical:: elem_in_crust,elem_in_mantle
-
+ ! flag for transverse isotropic elements
+ logical:: elem_is_tiso
+
! add topography of the Moho *before* adding the 3D crustal velocity model so that the streched
! mesh gets assigned the right model values
elem_in_crust = .false.
elem_in_mantle = .false.
+ elem_is_tiso = .false.
if( iregion_code == IREGION_CRUST_MANTLE ) then
if( CRUSTAL .and. CASE_3D ) then
+ ! 3D crustal models
if( idoubling(ispec) == IFLAG_CRUST &
.or. idoubling(ispec) == IFLAG_220_80 &
.or. idoubling(ispec) == IFLAG_80_MOHO ) then
@@ -130,11 +136,49 @@
call moho_stretching_honor_crust(myrank, &
xelm,yelm,zelm,RMOHO_FICTITIOUS_IN_MESHER,&
R220,RMIDDLE_CRUST,elem_in_crust,elem_in_mantle)
- endif
+ endif
+ else
+ ! element below 220km
+ ! sets element flag for mantle
+ elem_in_mantle = .true.
endif
+ else
+ ! 1D crust, no stretching
+ ! sets element flags
+ if( idoubling(ispec) == IFLAG_CRUST ) then
+ elem_in_crust = .true.
+ else
+ elem_in_mantle = .true.
+ endif
endif
- endif
-
+
+ ! sets transverse isotropic flag for elements in mantle
+ if( TRANSVERSE_ISOTROPY ) then
+ ! modifies tiso to have it for all mantle elements
+ ! preferred for example, when using 1Dref (STW model)
+ ! note: this will increase the computation time by ~ 45 %
+ if( USE_FULL_TISO_MANTLE ) then
+ ! fully transverse isotropic mantle
+ ! preferred for harvard (kustowski's) models using STW 1D reference, i.e.
+ ! THREE_D_MODEL_S362ANI
+ ! THREE_D_MODEL_S362WMANI
+ ! THREE_D_MODEL_S29EA
+ ! THREE_D_MODEL_GLL
+ if( elem_in_mantle ) elem_is_tiso = .true.
+ else if( idoubling(ispec)==IFLAG_220_80 .or. idoubling(ispec)==IFLAG_80_MOHO ) then
+ ! default case:
+ ! models use only transverse isotropy between moho and 220 km depth
+ elem_is_tiso = .true.
+ ! checks mantle flag to be sure
+ if( elem_in_mantle .eqv. .false. ) stop 'error mantle flag confused between moho and 220'
+ endif
+ endif
+
+ endif ! IREGION_CRUST_MANTLE
+
+ ! sets element tiso flag
+ ispec_is_tiso(ispec) = elem_is_tiso
+
! interpolates and stores GLL point locations
call compute_element_GLL_locations(xelm,yelm,zelm,ispec,nspec, &
xstore,ystore,zstore,shape3D)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -39,8 +39,9 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
- rho_vp,rho_vs,ABSORBING_CONDITIONS,ACTUALLY_STORE_ARRAYS,xigll,yigll,zigll)
+ nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ rho_vp,rho_vs,ABSORBING_CONDITIONS,ACTUALLY_STORE_ARRAYS,xigll,yigll,zigll, &
+ ispec_is_tiso)
! creates the inner core cube of the mesh
@@ -120,6 +121,8 @@
logical :: ACTUALLY_STORE_ARRAYS,ABSORBING_CONDITIONS
+ logical, dimension(nspec) :: ispec_is_tiso
+
!local parameters
double precision, dimension(NGNOD) :: xelm,yelm,zelm
! parameters needed to store the radii of the grid points in the spherically symmetric Earth
@@ -267,7 +270,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
- xigll,yigll,zigll)
+ xigll,yigll,zigll,ispec_is_tiso)
enddo
enddo
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -51,7 +51,8 @@
normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,&
ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot, &
- CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,offset_proc_xi,offset_proc_eta)
+ CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,offset_proc_xi,offset_proc_eta, &
+ ispec_is_tiso)
! adds doubling elements to the different regions of the mesh
@@ -154,6 +155,8 @@
integer :: offset_proc_xi,offset_proc_eta
logical :: CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA
+ logical, dimension(nspec) :: ispec_is_tiso
+
! local parameters
double precision, dimension(NGLOB_DOUBLING_SUPERBRICK) :: x_superbrick,y_superbrick,z_superbrick
double precision, dimension(NGNOD) :: offset_x,offset_y,offset_z
@@ -344,7 +347,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
- xigll,yigll,zigll)
+ xigll,yigll,zigll,ispec_is_tiso)
! boundary mesh
if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -226,6 +226,9 @@
ispec2D_670_top,ispec2D_670_bot
double precision r_moho,r_400,r_670
+ ! flags for transverse isotropic elements
+ logical, dimension(:), allocatable :: ispec_is_tiso
+
! create the name for the database of the current slide and region
call create_name_database(prname,myrank,iregion_code,LOCAL_PATH)
@@ -238,53 +241,64 @@
else
nspec_att = 1
end if
- allocate(Qmu_store(NGLLX,NGLLY,NGLLZ,nspec_att))
- allocate(tau_e_store(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att))
-
+ allocate(Qmu_store(NGLLX,NGLLY,NGLLZ,nspec_att), &
+ tau_e_store(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att),stat=ier)
+ if(ier /= 0) stop 'error in allocate 1'
+
! Gauss-Lobatto-Legendre points of integration
- allocate(xigll(NGLLX))
- allocate(yigll(NGLLY))
- allocate(zigll(NGLLZ))
+ allocate(xigll(NGLLX), &
+ yigll(NGLLY), &
+ zigll(NGLLZ),stat=ier)
+ if(ier /= 0) stop 'error in allocate 2'
! Gauss-Lobatto-Legendre weights of integration
- allocate(wxgll(NGLLX))
- allocate(wygll(NGLLY))
- allocate(wzgll(NGLLZ))
+ allocate(wxgll(NGLLX), &
+ wygll(NGLLY), &
+ wzgll(NGLLZ),stat=ier)
+ if(ier /= 0) stop 'error in allocate 3'
! 3D shape functions and their derivatives
- allocate(shape3D(NGNOD,NGLLX,NGLLY,NGLLZ))
- allocate(dershape3D(NDIM,NGNOD,NGLLX,NGLLY,NGLLZ))
+ allocate(shape3D(NGNOD,NGLLX,NGLLY,NGLLZ), &
+ dershape3D(NDIM,NGNOD,NGLLX,NGLLY,NGLLZ),stat=ier)
+ if(ier /= 0) stop 'error in allocat 4'
! 2D shape functions and their derivatives
- allocate(shape2D_x(NGNOD2D,NGLLY,NGLLZ))
- allocate(shape2D_y(NGNOD2D,NGLLX,NGLLZ))
- allocate(shape2D_bottom(NGNOD2D,NGLLX,NGLLY))
- allocate(shape2D_top(NGNOD2D,NGLLX,NGLLY))
- allocate(dershape2D_x(NDIM2D,NGNOD2D,NGLLY,NGLLZ))
- allocate(dershape2D_y(NDIM2D,NGNOD2D,NGLLX,NGLLZ))
- allocate(dershape2D_bottom(NDIM2D,NGNOD2D,NGLLX,NGLLY))
- allocate(dershape2D_top(NDIM2D,NGNOD2D,NGLLX,NGLLY))
+ allocate(shape2D_x(NGNOD2D,NGLLY,NGLLZ), &
+ shape2D_y(NGNOD2D,NGLLX,NGLLZ), &
+ shape2D_bottom(NGNOD2D,NGLLX,NGLLY), &
+ shape2D_top(NGNOD2D,NGLLX,NGLLY), &
+ dershape2D_x(NDIM2D,NGNOD2D,NGLLY,NGLLZ), &
+ dershape2D_y(NDIM2D,NGNOD2D,NGLLX,NGLLZ), &
+ dershape2D_bottom(NDIM2D,NGNOD2D,NGLLX,NGLLY), &
+ dershape2D_top(NDIM2D,NGNOD2D,NGLLX,NGLLY),stat=ier)
+ if(ier /= 0) stop 'error in allocate 5'
! array with model density
- allocate(rhostore(NGLLX,NGLLY,NGLLZ,nspec))
- allocate(dvpstore(NGLLX,NGLLY,NGLLZ,nspec))
+ allocate(rhostore(NGLLX,NGLLY,NGLLZ,nspec), &
+ dvpstore(NGLLX,NGLLY,NGLLZ,nspec),stat=ier)
+ if(ier /= 0) stop 'error in allocate 6'
! for anisotropy
- allocate(kappavstore(NGLLX,NGLLY,NGLLZ,nspec))
- allocate(muvstore(NGLLX,NGLLY,NGLLZ,nspec))
+ allocate(kappavstore(NGLLX,NGLLY,NGLLZ,nspec), &
+ muvstore(NGLLX,NGLLY,NGLLZ,nspec), &
+ kappahstore(NGLLX,NGLLY,NGLLZ,nspec), &
+ muhstore(NGLLX,NGLLY,NGLLZ,nspec), &
+ eta_anisostore(NGLLX,NGLLY,NGLLZ,nspec), &
+ ispec_is_tiso(nspec),stat=ier)
+ if(ier /= 0) stop 'error in allocate 7'
+
+ ! initializes flags for transverse isotropic elements
+ ispec_is_tiso(:) = .false.
- allocate(kappahstore(NGLLX,NGLLY,NGLLZ,nspec))
- allocate(muhstore(NGLLX,NGLLY,NGLLZ,nspec))
- allocate(eta_anisostore(NGLLX,NGLLY,NGLLZ,nspec))
-
! Stacey absorbing boundaries
if(NCHUNKS /= 6) then
nspec_stacey = nspec
else
nspec_stacey = 1
endif
- allocate(rho_vp(NGLLX,NGLLY,NGLLZ,nspec_stacey))
- allocate(rho_vs(NGLLX,NGLLY,NGLLZ,nspec_stacey))
+ allocate(rho_vp(NGLLX,NGLLY,NGLLZ,nspec_stacey), &
+ rho_vs(NGLLX,NGLLY,NGLLZ,nspec_stacey),stat=ier)
+ if(ier /= 0) stop 'error in allocate 8'
! anisotropy
if((ANISOTROPIC_INNER_CORE .and. iregion_code == IREGION_INNER_CORE) .or. &
@@ -293,65 +307,72 @@
else
nspec_ani = 1
endif
- allocate(c11store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c12store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c13store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c14store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c15store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c16store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c22store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c23store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c24store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c25store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c26store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c33store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c34store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c35store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c36store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c44store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c45store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c46store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c55store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c56store(NGLLX,NGLLY,NGLLZ,nspec_ani))
- allocate(c66store(NGLLX,NGLLY,NGLLZ,nspec_ani))
+ allocate(c11store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c12store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c13store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c14store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c15store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c16store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c22store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c23store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c24store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c25store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c26store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c33store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c34store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c35store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c36store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c44store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c45store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c46store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c55store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c56store(NGLLX,NGLLY,NGLLZ,nspec_ani), &
+ c66store(NGLLX,NGLLY,NGLLZ,nspec_ani),stat=ier)
+ if(ier /= 0) stop 'error in allocate 9'
! boundary locator
- allocate(iboun(6,nspec))
+ allocate(iboun(6,nspec),stat=ier)
+ if(ier /= 0) stop 'error in allocate 10'
! boundary parameters locator
- allocate(ibelm_xmin(NSPEC2DMAX_XMIN_XMAX))
- allocate(ibelm_xmax(NSPEC2DMAX_XMIN_XMAX))
- allocate(ibelm_ymin(NSPEC2DMAX_YMIN_YMAX))
- allocate(ibelm_ymax(NSPEC2DMAX_YMIN_YMAX))
- allocate(ibelm_bottom(NSPEC2D_BOTTOM))
- allocate(ibelm_top(NSPEC2D_TOP))
+ allocate(ibelm_xmin(NSPEC2DMAX_XMIN_XMAX), &
+ ibelm_xmax(NSPEC2DMAX_XMIN_XMAX), &
+ ibelm_ymin(NSPEC2DMAX_YMIN_YMAX), &
+ ibelm_ymax(NSPEC2DMAX_YMIN_YMAX), &
+ ibelm_bottom(NSPEC2D_BOTTOM), &
+ ibelm_top(NSPEC2D_TOP),stat=ier)
+ if(ier /= 0) stop 'error in allocate 11'
! 2-D jacobians and normals
- allocate(jacobian2D_xmin(NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX))
- allocate(jacobian2D_xmax(NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX))
- allocate(jacobian2D_ymin(NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX))
- allocate(jacobian2D_ymax(NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX))
- allocate(jacobian2D_bottom(NGLLX,NGLLY,NSPEC2D_BOTTOM))
- allocate(jacobian2D_top(NGLLX,NGLLY,NSPEC2D_TOP))
+ allocate(jacobian2D_xmin(NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX), &
+ jacobian2D_xmax(NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX), &
+ jacobian2D_ymin(NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX), &
+ jacobian2D_ymax(NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX), &
+ jacobian2D_bottom(NGLLX,NGLLY,NSPEC2D_BOTTOM), &
+ jacobian2D_top(NGLLX,NGLLY,NSPEC2D_TOP),stat=ier)
+ if(ier /= 0) stop 'error in allocate 12'
- allocate(normal_xmin(NDIM,NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX))
- allocate(normal_xmax(NDIM,NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX))
- allocate(normal_ymin(NDIM,NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX))
- allocate(normal_ymax(NDIM,NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX))
- allocate(normal_bottom(NDIM,NGLLX,NGLLY,NSPEC2D_BOTTOM))
- allocate(normal_top(NDIM,NGLLX,NGLLY,NSPEC2D_TOP))
+ allocate(normal_xmin(NDIM,NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX), &
+ normal_xmax(NDIM,NGLLY,NGLLZ,NSPEC2DMAX_XMIN_XMAX), &
+ normal_ymin(NDIM,NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX), &
+ normal_ymax(NDIM,NGLLX,NGLLZ,NSPEC2DMAX_YMIN_YMAX), &
+ normal_bottom(NDIM,NGLLX,NGLLY,NSPEC2D_BOTTOM), &
+ normal_top(NDIM,NGLLX,NGLLY,NSPEC2D_TOP),stat=ier)
+ if(ier /= 0) stop 'error in allocate 13'
! Stacey
- allocate(nimin(2,NSPEC2DMAX_YMIN_YMAX))
- allocate(nimax(2,NSPEC2DMAX_YMIN_YMAX))
- allocate(njmin(2,NSPEC2DMAX_XMIN_XMAX))
- allocate(njmax(2,NSPEC2DMAX_XMIN_XMAX))
- allocate(nkmin_xi(2,NSPEC2DMAX_XMIN_XMAX))
- allocate(nkmin_eta(2,NSPEC2DMAX_YMIN_YMAX))
+ allocate(nimin(2,NSPEC2DMAX_YMIN_YMAX), &
+ nimax(2,NSPEC2DMAX_YMIN_YMAX), &
+ njmin(2,NSPEC2DMAX_XMIN_XMAX), &
+ njmax(2,NSPEC2DMAX_XMIN_XMAX), &
+ nkmin_xi(2,NSPEC2DMAX_XMIN_XMAX), &
+ nkmin_eta(2,NSPEC2DMAX_YMIN_YMAX),stat=ier)
+ if(ier /= 0) stop 'error in allocate 14'
! MPI cut-planes parameters along xi and along eta
- allocate(iMPIcut_xi(2,nspec))
- allocate(iMPIcut_eta(2,nspec))
+ allocate(iMPIcut_xi(2,nspec), &
+ iMPIcut_eta(2,nspec),stat=ier)
+ if(ier /= 0) stop 'error in allocate 15'
! store and save the final arrays only in the second pass
! therefore in the first pass some arrays can be allocated with a dummy size
@@ -362,16 +383,17 @@
ACTUALLY_STORE_ARRAYS = .true.
nspec_actually = nspec
endif
- allocate(xixstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(xiystore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(xizstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(etaxstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(etaystore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(etazstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(gammaxstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(gammaystore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(gammazstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier); if(ier /= 0) stop 'error in allocate'
-
+ allocate(xixstore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ xiystore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ xizstore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ etaxstore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ etaystore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ etazstore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ gammaxstore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ gammaystore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
+ gammazstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier)
+ if(ier /= 0) stop 'error in allocate 16'
+
! boundary mesh
if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
NSPEC2D_MOHO = NSPEC2D_TOP
@@ -382,15 +404,16 @@
NSPEC2D_400 = 1
NSPEC2D_670 = 1
endif
- allocate(ibelm_moho_top(NSPEC2D_MOHO),ibelm_moho_bot(NSPEC2D_MOHO))
- allocate(ibelm_400_top(NSPEC2D_400),ibelm_400_bot(NSPEC2D_400))
- allocate(ibelm_670_top(NSPEC2D_670),ibelm_670_bot(NSPEC2D_670))
- allocate(normal_moho(NDIM,NGLLX,NGLLY,NSPEC2D_MOHO))
- allocate(normal_400(NDIM,NGLLX,NGLLY,NSPEC2D_400))
- allocate(normal_670(NDIM,NGLLX,NGLLY,NSPEC2D_670))
- allocate(jacobian2D_moho(NGLLX,NGLLY,NSPEC2D_MOHO))
- allocate(jacobian2D_400(NGLLX,NGLLY,NSPEC2D_400))
- allocate(jacobian2D_670(NGLLX,NGLLY,NSPEC2D_670))
+ allocate(ibelm_moho_top(NSPEC2D_MOHO),ibelm_moho_bot(NSPEC2D_MOHO), &
+ ibelm_400_top(NSPEC2D_400),ibelm_400_bot(NSPEC2D_400), &
+ ibelm_670_top(NSPEC2D_670),ibelm_670_bot(NSPEC2D_670), &
+ normal_moho(NDIM,NGLLX,NGLLY,NSPEC2D_MOHO), &
+ normal_400(NDIM,NGLLX,NGLLY,NSPEC2D_400), &
+ normal_670(NDIM,NGLLX,NGLLY,NSPEC2D_670), &
+ jacobian2D_moho(NGLLX,NGLLY,NSPEC2D_MOHO), &
+ jacobian2D_400(NGLLX,NGLLY,NSPEC2D_400), &
+ jacobian2D_670(NGLLX,NGLLY,NSPEC2D_670),stat=ier)
+ if(ier /= 0) stop 'error in allocate 17'
! initialize number of layers
call crm_initialize_layers(myrank,ipass,xigll,yigll,zigll,wxgll,wygll,wzgll, &
@@ -405,7 +428,8 @@
first_layer_aniso,last_layer_aniso,nb_layer_above_aniso,is_on_a_slice_edge)
! to consider anisotropic elements first and to build the mesh from the bottom to the top of the region
- allocate (perm_layer(ifirst_region:ilast_region))
+ allocate (perm_layer(ifirst_region:ilast_region),stat=ier)
+ if(ier /= 0) stop 'error in allocate 18'
perm_layer = (/ (i, i=ilast_region,ifirst_region,-1) /)
if(iregion_code == IREGION_CRUST_MANTLE) then
@@ -422,7 +446,8 @@
! crustal layer stretching: element layer's top and bottom radii will get stretched when in crust
! (number of element layers in crust can vary for different resolutions and 1chunk simulations)
- allocate(stretch_tab(2,ner(1)))
+ allocate(stretch_tab(2,ner(1)),stat=ier)
+ if(ier /= 0) stop 'error in allocate 19'
if (CASE_3D .and. iregion_code == IREGION_CRUST_MANTLE .and. .not. SUPPRESS_CRUSTAL_MESH) then
! stretching function determines top and bottom of each element layer in the
! crust region (between r_top(1) and r_bottom(1)), where ner(1) is the
@@ -514,7 +539,8 @@
ibelm_moho_top,ibelm_moho_bot,ibelm_400_top,ibelm_400_bot,ibelm_670_top,ibelm_670_bot, &
normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,&
- ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot)
+ ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot,&
+ ispec_is_tiso)
! mesh doubling elements
@@ -545,7 +571,8 @@
normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,&
ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot, &
- CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,offset_proc_xi,offset_proc_eta)
+ CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA,offset_proc_xi,offset_proc_eta, &
+ ispec_is_tiso)
enddo !ilayer_loop
@@ -566,7 +593,8 @@
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
- rho_vp,rho_vs,ABSORBING_CONDITIONS,ACTUALLY_STORE_ARRAYS,xigll,yigll,zigll)
+ rho_vp,rho_vs,ABSORBING_CONDITIONS,ACTUALLY_STORE_ARRAYS,xigll,yigll,zigll, &
+ ispec_is_tiso)
! check total number of spectral elements created
@@ -602,12 +630,13 @@
! NEX_XI,NCHUNKS,ABSORBING_CONDITIONS,PPM_V )
! allocate memory for arrays
- allocate(locval(npointot),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(ifseg(npointot),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(xp(npointot),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(yp(npointot),stat=ier); if(ier /= 0) stop 'error in allocate'
- allocate(zp(npointot),stat=ier); if(ier /= 0) stop 'error in allocate'
-
+ allocate(locval(npointot), &
+ ifseg(npointot), &
+ xp(npointot), &
+ yp(npointot), &
+ zp(npointot),stat=ier)
+ if(ier /= 0) stop 'error in allocate 20'
+
locval = 0
ifseg = .false.
xp = 0.d0
@@ -714,8 +743,11 @@
! count number of anisotropic elements in current region
! should be zero in all the regions except in the mantle
- nspec_tiso = count(idoubling(1:nspec) == IFLAG_220_80) + count(idoubling(1:nspec) == IFLAG_80_MOHO)
-
+ ! (used only for checks in meshfem3D() routine)
+ !nspec_tiso = count(idoubling(1:nspec) == IFLAG_220_80) + count(idoubling(1:nspec) == IFLAG_80_MOHO)
+ nspec_tiso = count(ispec_is_tiso(:))
+
+ ! precomputes jacobian for 2d absorbing boundary surfaces
call get_jacobian_boundaries(myrank,iboun,nspec,xstore,ystore,zstore, &
dershape2D_x,dershape2D_y,dershape2D_bottom,dershape2D_top, &
ibelm_xmin,ibelm_xmax,ibelm_ymin,ibelm_ymax,ibelm_bottom,ibelm_top, &
@@ -731,7 +763,8 @@
xigll,yigll,zigll)
! allocates mass matrix in this slice (will be fully assembled in the solver)
- allocate(rmass(nglob),stat=ier); if(ier /= 0) stop 'error in allocate'
+ allocate(rmass(nglob),stat=ier)
+ if(ier /= 0) stop 'error in allocate 21'
! allocates ocean load mass matrix as well if oceans
if(OCEANS .and. iregion_code == IREGION_CRUST_MANTLE) then
nglob_oceans = nglob
@@ -739,7 +772,8 @@
! allocate dummy array if no oceans
nglob_oceans = 1
endif
- allocate(rmass_ocean_load(nglob_oceans))
+ allocate(rmass_ocean_load(nglob_oceans),stat=ier)
+ if(ier /= 0) stop 'error in allocate 22'
! creating mass matrix in this slice (will be fully assembled in the solver)
call create_mass_matrices(myrank,nspec,idoubling,wxgll,wygll,wzgll,ibool, &
@@ -770,7 +804,7 @@
ANISOTROPIC_INNER_CORE,OCEANS, &
tau_s,tau_e_store,Qmu_store,T_c_source,ATTENUATION, &
size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4),size(tau_e_store,5),&
- ABSORBING_CONDITIONS,SAVE_MESH_FILES)
+ ABSORBING_CONDITIONS,SAVE_MESH_FILES,ispec_is_tiso)
deallocate(rmass,stat=ier); if(ier /= 0) stop 'error in deallocate'
deallocate(rmass_ocean_load,stat=ier); if(ier /= 0) stop 'error in deallocate'
@@ -826,6 +860,7 @@
deallocate(rhostore,dvpstore,kappavstore,kappahstore)
deallocate(muvstore,muhstore)
deallocate(eta_anisostore)
+ deallocate(ispec_is_tiso)
deallocate(c11store)
deallocate(c12store)
deallocate(c13store)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -51,8 +51,9 @@
NSPEC2D_MOHO,NSPEC2D_400,NSPEC2D_670,nex_eta_moho, &
ibelm_moho_top,ibelm_moho_bot,ibelm_400_top,ibelm_400_bot,ibelm_670_top,ibelm_670_bot, &
normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
- ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,&
- ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot)
+ ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top, &
+ ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot, &
+ ispec_is_tiso)
! adds a regular spectral element to the different regions of the mesh
@@ -160,6 +161,8 @@
integer ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top, &
ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot
+ logical, dimension(nspec) :: ispec_is_tiso
+
! local parameters
double precision, dimension(NGNOD) :: offset_x,offset_y,offset_z
double precision, dimension(NGNOD) :: xelm,yelm,zelm
@@ -263,7 +266,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
- xigll,yigll,zigll)
+ xigll,yigll,zigll,ispec_is_tiso)
! boundary mesh
if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -278,7 +278,7 @@
! correct number of spectral elements in each block depending on chunk type
- integer nspec_aniso,npointot
+ integer nspec_tiso,npointot
! parameters needed to store the radii of the grid points
! in the spherically symmetric Earth
@@ -659,7 +659,7 @@
call create_regions_mesh(iregion_code,ibool,idoubling,is_on_a_slice_edge, &
xstore,ystore,zstore,rmins,rmaxs, &
- iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_aniso, &
+ iproc_xi,iproc_eta,ichunk,NSPEC(iregion_code),nspec_tiso, &
volume_local,area_local_bottom,area_local_top, &
nglob(iregion_code),npointot, &
NEX_XI,NEX_PER_PROC_XI,NEX_PER_PROC_ETA, &
@@ -679,10 +679,10 @@
enddo
! store number of anisotropic elements found in the mantle
- if(nspec_aniso /= 0 .and. iregion_code /= IREGION_CRUST_MANTLE) &
+ if(nspec_tiso /= 0 .and. iregion_code /= IREGION_CRUST_MANTLE) &
call exit_MPI(myrank,'found anisotropic elements outside of the mantle')
- if(iregion_code == IREGION_CRUST_MANTLE .and. nspec_aniso == 0) &
+ if(iregion_code == IREGION_CRUST_MANTLE .and. nspec_tiso == 0) &
call exit_MPI(myrank,'found no anisotropic elements in the mantle')
! computes total area and volume
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -43,7 +43,7 @@
TRANSVERSE_ISOTROPY,HETEROGEN_3D_MANTLE,ANISOTROPIC_3D_MANTLE, &
ANISOTROPIC_INNER_CORE,OCEANS, &
tau_s,tau_e_store,Qmu_store,T_c_source,ATTENUATION,vx,vy,vz,vnspec, &
- ABSORBING_CONDITIONS,SAVE_MESH_FILES)
+ ABSORBING_CONDITIONS,SAVE_MESH_FILES,ispec_is_tiso)
implicit none
@@ -152,6 +152,8 @@
logical ABSORBING_CONDITIONS,SAVE_MESH_FILES
+ logical, dimension(nspec) :: ispec_is_tiso
+
! local parameters
integer i,j,k,ispec,iglob,nspec1, nglob1
real(kind=CUSTOM_REAL) scaleval1,scaleval2
@@ -338,6 +340,8 @@
write(27) is_on_a_slice_edge
+ write(27) ispec_is_tiso
+
close(27)
! absorbing boundary parameters
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -241,8 +241,19 @@
! ibool_outer_core
static_memory_size = static_memory_size + dble(NGLLX)*dble(NGLLY)*dble(NGLLZ)*NSPEC(IREGION_OUTER_CORE)*dble(SIZE_INTEGER)
-! idoubling_crust_mantle
- static_memory_size = static_memory_size + NSPEC(IREGION_CRUST_MANTLE)*dble(SIZE_INTEGER)
+! idoubling_crust_mantle (not needed anymore..)
+! static_memory_size = static_memory_size + NSPEC(IREGION_CRUST_MANTLE)*dble(SIZE_INTEGER)
+! idoubling_outer_core
+ static_memory_size = static_memory_size + NSPEC(IREGION_OUTER_CORE)*dble(SIZE_INTEGER)
+! idoubling_inner_core
+ static_memory_size = static_memory_size + NSPEC(IREGION_INNER_CORE)*dble(SIZE_INTEGER)
+
+! ispec_is_tiso_crust_mantle
+ static_memory_size = static_memory_size + NSPEC(IREGION_CRUST_MANTLE)*dble(SIZE_LOGICAL)
+! ispec_is_tiso_outer_core
+ static_memory_size = static_memory_size + NSPEC(IREGION_OUTER_CORE)*dble(SIZE_LOGICAL)
+! ispec_is_tiso_inner_core
+ static_memory_size = static_memory_size + NSPEC(IREGION_INNER_CORE)*dble(SIZE_LOGICAL)
! xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle,rmass_crust_mantle
static_memory_size = static_memory_size + NGLOB(IREGION_CRUST_MANTLE)*4.d0*dble(CUSTOM_REAL)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -26,7 +26,8 @@
!=====================================================================
subroutine compute_boundary_kernel(displ,accel,b_displ,nspec,iregion_code, &
- ystore,zstore,ibool,idoubling, &
+ ystore,zstore,ibool,ispec_is_tiso, &
+ !--- idoubling, &
xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz, &
hprime_xx,hprime_yy,hprime_zz, &
rhostore,kappavstore,muvstore,kappahstore,muhstore,eta_anisostore, &
@@ -42,7 +43,8 @@
real(kind=CUSTOM_REAL), dimension(NDIM,*) :: displ,accel,b_displ
integer nspec, iregion_code
integer, dimension(NGLLX,NGLLY,NGLLZ,nspec) :: ibool
- integer, dimension(*) :: idoubling
+! integer, dimension(*) :: idoubling
+ logical, dimension(*) :: ispec_is_tiso
real(kind=CUSTOM_REAL), dimension(*) :: ystore,zstore
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLX) :: hprime_xx
@@ -128,7 +130,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,idoubling)
+ ystore,zstore,ibool,ispec_is_tiso) !idoubling)
! ----- forward strain -------
temp1(:) = matmul(b_displl(:,:,j,k), hprime_xx(i,:))
@@ -153,7 +155,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,idoubling)
+ ystore,zstore,ibool,ispec_is_tiso) !-- idoubling)
! ---- precompute K_d for F-S boundaries ----
if (fluid_solid_boundary) then
@@ -236,7 +238,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,idoubling)
+ ystore,zstore,ibool,ispec_is_tiso) !--idoubling)
implicit none
@@ -254,8 +256,9 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store
real(kind=CUSTOM_REAL), dimension(*) :: ystore,zstore
integer, dimension(NGLLX,NGLLY,NGLLZ,*) :: ibool
- integer, dimension(*) :: idoubling
-
+! integer, dimension(*) :: idoubling
+ logical, dimension(*) :: ispec_is_tiso
+
! --- local variables ---
real(kind=CUSTOM_REAL) :: duxdxl_plus_duydyl,duxdxl_plus_duzdzl,duydyl_plus_duzdzl
real(kind=CUSTOM_REAL) :: duxdyl_plus_duydxl,duzdxl_plus_duxdzl,duzdyl_plus_duydzl
@@ -332,8 +335,11 @@
sigma(2,3) = c14*duxdxl + c46*duxdyl_plus_duydxl + c24*duydyl + &
c45*duzdxl_plus_duxdzl + c44*duzdyl_plus_duydzl + c34*duzdzl
- else if(.not. (TRANSVERSE_ISOTROPY_VAL .and. (idoubling(ispec) == IFLAG_80_MOHO .or. idoubling(ispec) == IFLAG_220_80))) then
+ 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
+
kappal = kappavstore(i,j,k,ispec)
mul = muvstore(i,j,k,ispec)
@@ -350,6 +356,8 @@
else
+ ! transverse isotropic elements
+
kappavl = kappavstore(i,j,k,ispec)
muvl = muvstore(i,j,k,ispec)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -53,7 +53,9 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- ibool,idoubling,R_memory,epsilondev,epsilon_trace_over_3,one_minus_sum_beta, &
+ ibool,ispec_is_tiso, &
+ !-- idoubling,
+ R_memory,epsilondev,epsilon_trace_over_3,one_minus_sum_beta, &
alphaval,betaval,gammaval,factor_common,vx,vy,vz,vnspec)
implicit none
@@ -85,7 +87,8 @@
! end type model_attenuation_variables
! array with the local to global mapping per slice
- integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling
+! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling
+ logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso
! displacement and acceleration
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) :: displ_crust_mantle,accel_crust_mantle
@@ -456,8 +459,8 @@
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. (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)
mul = muvstore(i,j,k,ispec)
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 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -53,7 +53,9 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- ibool,idoubling,R_memory,epsilondev,epsilon_trace_over_3,one_minus_sum_beta, &
+ ibool,ispec_is_tiso, &
+ ! --idoubling,
+ R_memory,epsilondev,epsilon_trace_over_3,one_minus_sum_beta, &
alphaval,betaval,gammaval,factor_common,vx,vy,vz,vnspec)
! this routine is optimized for NGLLX = NGLLY = NGLLZ = 5 using the Deville et al. (2002) inlined matrix-matrix products
@@ -113,7 +115,8 @@
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
! array with the local to global mapping per slice
- integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling
+! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling
+ logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso
! gravity
double precision, dimension(NRAD_GRAVITY) :: minus_gravity_table,density_table,minus_deriv_gravity_table
@@ -542,9 +545,10 @@
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. (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)
mul = muvstore(i,j,k,ispec)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -502,7 +502,8 @@
c22store_crust_mantle,c23store_crust_mantle, &
c33store_crust_mantle,c44store_crust_mantle, &
c55store_crust_mantle,c66store_crust_mantle, &
- muvstore_crust_mantle,muhstore_crust_mantle,idoubling_crust_mantle, &
+ muvstore_crust_mantle,muhstore_crust_mantle,ispec_is_tiso_crust_mantle, &
+ !----- idoubling_crust_mantle, &
muvstore_inner_core, &
SIMULATION_TYPE,MOVIE_VOLUME,muvstore_crust_mantle_3dmovie, &
c11store_inner_core,c12store_inner_core,c13store_inner_core, &
@@ -535,7 +536,8 @@
muvstore_crust_mantle
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPECMAX_TISO_MANTLE) :: &
muhstore_crust_mantle
- integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
+! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
+ logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso_crust_mantle
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: &
muvstore_inner_core
@@ -644,9 +646,13 @@
muvstore_crust_mantle_3dmovie(i,j,k,ispec)=muvstore_crust_mantle(i,j,k,ispec)
endif
muvstore_crust_mantle(i,j,k,ispec) = muvstore_crust_mantle(i,j,k,ispec) * scale_factor
- if(TRANSVERSE_ISOTROPY_VAL .and. (idoubling_crust_mantle(ispec) == IFLAG_220_80 &
- .or. idoubling_crust_mantle(ispec) == IFLAG_80_MOHO)) &
+
+ ! scales transverse isotropic values for mu_h
+ !if(TRANSVERSE_ISOTROPY_VAL .and. (idoubling_crust_mantle(ispec) == IFLAG_220_80 &
+ ! .or. idoubling_crust_mantle(ispec) == IFLAG_80_MOHO)) &
+ if( ispec_is_tiso_crust_mantle(ispec) ) then
muhstore_crust_mantle(i,j,k,ispec) = muhstore_crust_mantle(i,j,k,ispec) * scale_factor
+ endif
endif
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -35,7 +35,8 @@
c11store,c12store,c13store,c14store,c15store,c16store,c22store, &
c23store,c24store,c25store,c26store,c33store,c34store,c35store, &
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
- ibool,idoubling,is_on_a_slice_edge,rmass,rmass_ocean_load,nspec,nglob, &
+ ibool,idoubling,ispec_is_tiso, &
+ is_on_a_slice_edge,rmass,rmass_ocean_load,nspec,nglob, &
READ_KAPPA_MU,READ_TISO,TRANSVERSE_ISOTROPY, &
ANISOTROPIC_3D_MANTLE,ANISOTROPIC_INNER_CORE,OCEANS,LOCAL_PATH,ABSORBING_CONDITIONS)
@@ -92,6 +93,8 @@
! this for non blocking MPI
logical, dimension(nspec) :: is_on_a_slice_edge
+ logical, dimension(nspec) :: ispec_is_tiso
+
! processor identification
character(len=150) prname
@@ -191,6 +194,8 @@
read(IIN) is_on_a_slice_edge
+ read(IIN) ispec_is_tiso
+
close(IIN)
end subroutine read_arrays_solver
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -40,13 +40,16 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle,is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! -- idoubling_crust_mantle
+ is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
vp_outer_core,xstore_outer_core,ystore_outer_core,zstore_outer_core, &
xix_outer_core,xiy_outer_core,xiz_outer_core, &
etax_outer_core,etay_outer_core,etaz_outer_core, &
gammax_outer_core,gammay_outer_core,gammaz_outer_core, &
rhostore_outer_core,kappavstore_outer_core, &
- ibool_outer_core,idoubling_outer_core,is_on_a_slice_edge_outer_core,rmass_outer_core, &
+ ibool_outer_core,idoubling_outer_core,ispec_is_tiso_outer_core, &
+ is_on_a_slice_edge_outer_core,rmass_outer_core, &
xstore_inner_core,ystore_inner_core,zstore_inner_core, &
xix_inner_core,xiy_inner_core,xiz_inner_core, &
etax_inner_core,etay_inner_core,etaz_inner_core, &
@@ -54,7 +57,8 @@
rhostore_inner_core,kappavstore_inner_core,muvstore_inner_core, &
c11store_inner_core,c12store_inner_core,c13store_inner_core, &
c33store_inner_core,c44store_inner_core, &
- ibool_inner_core,idoubling_inner_core,is_on_a_slice_edge_inner_core,rmass_inner_core, &
+ ibool_inner_core,idoubling_inner_core,ispec_is_tiso_inner_core, &
+ is_on_a_slice_edge_inner_core,rmass_inner_core, &
ABSORBING_CONDITIONS,LOCAL_PATH)
implicit none
@@ -93,7 +97,10 @@
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
- integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
+
+! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
+ logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso_crust_mantle
+
! mass matrix
real(kind=CUSTOM_REAL), dimension(NGLOB_CRUST_MANTLE) :: rmass_crust_mantle
! additional mass matrix for ocean load
@@ -112,6 +119,8 @@
rhostore_outer_core,kappavstore_outer_core
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_OUTER_CORE) :: ibool_outer_core
integer, dimension(NSPEC_OUTER_CORE) :: idoubling_outer_core
+ logical, dimension(NSPEC_OUTER_CORE) :: ispec_is_tiso_outer_core
+
real(kind=CUSTOM_REAL), dimension(NGLOB_OUTER_CORE) :: rmass_outer_core
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: &
@@ -126,6 +135,8 @@
c13store_inner_core,c44store_inner_core
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: ibool_inner_core
integer, dimension(NSPEC_INNER_CORE) :: idoubling_inner_core
+ logical, dimension(NSPEC_INNER_CORE) :: ispec_is_tiso_inner_core
+
real(kind=CUSTOM_REAL), dimension(NGLOB_INNER_CORE) :: rmass_inner_core
logical ABSORBING_CONDITIONS
@@ -134,7 +145,8 @@
!local parameters
logical READ_KAPPA_MU,READ_TISO
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,1) :: dummy_array
-
+ integer, dimension(NSPEC_CRUST_MANTLE) :: dummy_i
+
! this for non blocking MPI
logical, dimension(NSPEC_CRUST_MANTLE) :: is_on_a_slice_edge_crust_mantle
logical, dimension(NSPEC_OUTER_CORE) :: is_on_a_slice_edge_outer_core
@@ -177,7 +189,10 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle,is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
+ ibool_crust_mantle,dummy_i, &
+ ! --idoubling_crust_mantle,
+ ispec_is_tiso_crust_mantle, &
+ is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
NSPEC_CRUST_MANTLE,NGLOB_CRUST_MANTLE, &
READ_KAPPA_MU,READ_TISO,TRANSVERSE_ISOTROPY_VAL,ANISOTROPIC_3D_MANTLE_VAL, &
ANISOTROPIC_INNER_CORE_VAL,OCEANS_VAL,LOCAL_PATH,ABSORBING_CONDITIONS)
@@ -206,7 +221,8 @@
dummy_array,dummy_array,dummy_array, &
dummy_array,dummy_array,dummy_array, &
dummy_array,dummy_array,dummy_array, &
- ibool_outer_core,idoubling_outer_core,is_on_a_slice_edge_outer_core,rmass_outer_core,rmass_ocean_load, &
+ ibool_outer_core,idoubling_outer_core,ispec_is_tiso_outer_core, &
+ is_on_a_slice_edge_outer_core,rmass_outer_core,rmass_ocean_load, &
NSPEC_OUTER_CORE,NGLOB_OUTER_CORE, &
READ_KAPPA_MU,READ_TISO,TRANSVERSE_ISOTROPY_VAL,ANISOTROPIC_3D_MANTLE_VAL, &
ANISOTROPIC_INNER_CORE_VAL,OCEANS_VAL,LOCAL_PATH,ABSORBING_CONDITIONS)
@@ -239,7 +255,8 @@
dummy_array,dummy_array,dummy_array, &
c44store_inner_core,dummy_array,dummy_array, &
dummy_array,dummy_array,dummy_array, &
- ibool_inner_core,idoubling_inner_core,is_on_a_slice_edge_inner_core,rmass_inner_core,rmass_ocean_load, &
+ ibool_inner_core,idoubling_inner_core,ispec_is_tiso_inner_core, &
+ is_on_a_slice_edge_inner_core,rmass_inner_core,rmass_ocean_load, &
NSPEC_INNER_CORE,NGLOB_INNER_CORE, &
READ_KAPPA_MU,READ_TISO,TRANSVERSE_ISOTROPY_VAL,ANISOTROPIC_3D_MANTLE_VAL, &
ANISOTROPIC_INNER_CORE_VAL,OCEANS_VAL,LOCAL_PATH,ABSORBING_CONDITIONS)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/save_kernels.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/save_kernels.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/save_kernels.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -33,7 +33,8 @@
rhostore_crust_mantle,muvstore_crust_mantle, &
kappavstore_crust_mantle,ibool_crust_mantle, &
kappahstore_crust_mantle,muhstore_crust_mantle, &
- eta_anisostore_crust_mantle,idoubling_crust_mantle, &
+ eta_anisostore_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
LOCAL_PATH)
implicit none
@@ -60,7 +61,8 @@
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
+! 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
@@ -154,9 +156,10 @@
! 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. (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
! A,C,L,N,F from isotropic model
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -107,7 +107,7 @@
integer :: irec,isource,nrec_tot_found,ier
integer :: icomp,itime,nadj_files_found,nadj_files_found_tot
character(len=3),dimension(NDIM) :: comp
- character(len=150) :: filename,adj_source_file,system_command,filename_new
+ character(len=256) :: filename,adj_source_file,system_command,filename_new
character(len=2) :: bic
! sources
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90 2011-04-04 03:07:55 UTC (rev 18166)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90 2011-04-04 04:19:03 UTC (rev 18167)
@@ -551,7 +551,8 @@
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle
! local to global mapping
- integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
+! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
+ logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso_crust_mantle
! mass matrix
real(kind=CUSTOM_REAL), dimension(NGLOB_CRUST_MANTLE) :: rmass_crust_mantle
@@ -577,6 +578,7 @@
! local to global mapping
integer, dimension(NSPEC_OUTER_CORE) :: idoubling_outer_core
+ logical, dimension(NSPEC_OUTER_CORE) :: ispec_is_tiso_outer_core ! only needed for compute_boundary_kernel()
! mass matrix
real(kind=CUSTOM_REAL), dimension(NGLOB_OUTER_CORE) :: rmass_outer_core
@@ -605,6 +607,7 @@
! local to global mapping
integer, dimension(NSPEC_INNER_CORE) :: idoubling_inner_core
+ logical, dimension(NSPEC_INNER_CORE) :: ispec_is_tiso_inner_core ! only needed for computer_boundary_kernel() routine
! mass matrix
real(kind=CUSTOM_REAL), dimension(NGLOB_INNER_CORE) :: rmass_inner_core
@@ -992,13 +995,16 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle,is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! -- idoubling_crust_mantle,
+ is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
vp_outer_core,xstore_outer_core,ystore_outer_core,zstore_outer_core, &
xix_outer_core,xiy_outer_core,xiz_outer_core, &
etax_outer_core,etay_outer_core,etaz_outer_core, &
gammax_outer_core,gammay_outer_core,gammaz_outer_core, &
rhostore_outer_core,kappavstore_outer_core, &
- ibool_outer_core,idoubling_outer_core,is_on_a_slice_edge_outer_core,rmass_outer_core, &
+ ibool_outer_core,idoubling_outer_core,ispec_is_tiso_outer_core, &
+ is_on_a_slice_edge_outer_core,rmass_outer_core, &
xstore_inner_core,ystore_inner_core,zstore_inner_core, &
xix_inner_core,xiy_inner_core,xiz_inner_core, &
etax_inner_core,etay_inner_core,etaz_inner_core, &
@@ -1006,7 +1012,8 @@
rhostore_inner_core,kappavstore_inner_core,muvstore_inner_core, &
c11store_inner_core,c12store_inner_core,c13store_inner_core, &
c33store_inner_core,c44store_inner_core, &
- ibool_inner_core,idoubling_inner_core,is_on_a_slice_edge_inner_core,rmass_inner_core, &
+ ibool_inner_core,idoubling_inner_core,ispec_is_tiso_inner_core, &
+ is_on_a_slice_edge_inner_core,rmass_inner_core, &
ABSORBING_CONDITIONS,LOCAL_PATH)
! read 2-D addressing for summation between slices with MPI
@@ -1686,7 +1693,7 @@
c22store_crust_mantle,c23store_crust_mantle, &
c33store_crust_mantle,c44store_crust_mantle, &
c55store_crust_mantle,c66store_crust_mantle, &
- muvstore_crust_mantle,muhstore_crust_mantle,idoubling_crust_mantle, &
+ muvstore_crust_mantle,muhstore_crust_mantle,ispec_is_tiso_crust_mantle, &
muvstore_inner_core, &
SIMULATION_TYPE,MOVIE_VOLUME,muvstore_crust_mantle_3dmovie, &
c11store_inner_core,c12store_inner_core,c13store_inner_core, &
@@ -2686,7 +2693,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
R_memory_crust_mantle,epsilondev_crust_mantle, &
eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
alphaval,betaval,gammaval,factor_common_crust_mantle, &
@@ -2729,7 +2737,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
R_memory_crust_mantle,epsilondev_crust_mantle, &
eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
alphaval,betaval,gammaval,factor_common_crust_mantle, &
@@ -2781,7 +2790,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! -- idoubling_crust_mantle, &
b_R_memory_crust_mantle,b_epsilondev_crust_mantle, &
b_eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
b_alphaval,b_betaval,b_gammaval,factor_common_crust_mantle, &
@@ -2824,7 +2834,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
b_R_memory_crust_mantle,b_epsilondev_crust_mantle, &
b_eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
b_alphaval,b_betaval,b_gammaval,factor_common_crust_mantle, &
@@ -3179,7 +3190,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ !---idoubling_crust_mantle, &
R_memory_crust_mantle,epsilondev_crust_mantle, &
eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
alphaval,betaval,gammaval,factor_common_crust_mantle, &
@@ -3222,7 +3234,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
R_memory_crust_mantle,epsilondev_crust_mantle, &
eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
alphaval,betaval,gammaval,factor_common_crust_mantle, &
@@ -3490,7 +3503,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
b_R_memory_crust_mantle,b_epsilondev_crust_mantle, &
b_eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
b_alphaval,b_betaval,b_gammaval,factor_common_crust_mantle, &
@@ -3533,7 +3547,8 @@
c34store_crust_mantle,c35store_crust_mantle,c36store_crust_mantle, &
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
- ibool_crust_mantle,idoubling_crust_mantle, &
+ ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ !--idoubling_crust_mantle, &
b_R_memory_crust_mantle,b_epsilondev_crust_mantle, &
b_eps_trace_over_3_crust_mantle,one_minus_sum_beta_crust_mantle, &
b_alphaval,b_betaval,b_gammaval,factor_common_crust_mantle, &
@@ -4016,7 +4031,8 @@
if (.not. SUPPRESS_CRUSTAL_MESH .and. HONOR_1D_SPHERICAL_MOHO) then
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! -- idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4032,7 +4048,8 @@
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4052,7 +4069,8 @@
! 400
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4068,7 +4086,8 @@
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4087,7 +4106,8 @@
! 670
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4103,7 +4123,8 @@
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4124,7 +4145,8 @@
iregion_code = IREGION_CRUST_MANTLE
call compute_boundary_kernel(displ_crust_mantle,accel_crust_mantle, &
b_displ_crust_mantle,nspec_crust_mantle,iregion_code, &
- ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,idoubling_crust_mantle, &
+ ystore_crust_mantle,zstore_crust_mantle,ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! -- idoubling_crust_mantle, &
xix_crust_mantle,xiy_crust_mantle,xiz_crust_mantle, &
etax_crust_mantle,etay_crust_mantle,etaz_crust_mantle,&
gammax_crust_mantle,gammay_crust_mantle,gammaz_crust_mantle,hprime_xx,hprime_yy,hprime_zz, &
@@ -4142,7 +4164,8 @@
iregion_code = IREGION_OUTER_CORE
call compute_boundary_kernel(vector_displ_outer_core,vector_accel_outer_core, &
b_vector_displ_outer_core,nspec_outer_core, &
- iregion_code,ystore_outer_core,zstore_outer_core,ibool_outer_core,idoubling_outer_core, &
+ iregion_code,ystore_outer_core,zstore_outer_core,ibool_outer_core,ispec_is_tiso_outer_core, &
+ ! --idoubling_outer_core, &
xix_outer_core,xiy_outer_core,xiz_outer_core, &
etax_outer_core,etay_outer_core,etaz_outer_core,&
gammax_outer_core,gammay_outer_core,gammaz_outer_core,hprime_xx,hprime_yy,hprime_zz, &
@@ -4163,7 +4186,8 @@
fluid_solid_boundary = .true.
call compute_boundary_kernel(vector_displ_outer_core,vector_accel_outer_core, &
b_vector_displ_outer_core,nspec_outer_core, &
- iregion_code,ystore_outer_core,zstore_outer_core,ibool_outer_core,idoubling_outer_core, &
+ iregion_code,ystore_outer_core,zstore_outer_core,ibool_outer_core,ispec_is_tiso_outer_core, &
+ ! --idoubling_outer_core, &
xix_outer_core,xiy_outer_core,xiz_outer_core, &
etax_outer_core,etay_outer_core,etaz_outer_core,&
gammax_outer_core,gammay_outer_core,gammaz_outer_core,hprime_xx,hprime_yy,hprime_zz, &
@@ -4181,7 +4205,8 @@
iregion_code = IREGION_INNER_CORE
call compute_boundary_kernel(displ_inner_core,accel_inner_core, &
b_displ_inner_core,nspec_inner_core,iregion_code, &
- ystore_inner_core,zstore_inner_core,ibool_inner_core,idoubling_inner_core, &
+ ystore_inner_core,zstore_inner_core,ibool_inner_core,ispec_is_tiso_inner_core, &
+ ! -- idoubling_inner_core, &
xix_inner_core,xiy_inner_core,xiz_inner_core, &
etax_inner_core,etay_inner_core,etaz_inner_core,&
gammax_inner_core,gammay_inner_core,gammaz_inner_core,hprime_xx,hprime_yy,hprime_zz, &
@@ -4399,7 +4424,8 @@
rhostore_crust_mantle,muvstore_crust_mantle, &
kappavstore_crust_mantle,ibool_crust_mantle, &
kappahstore_crust_mantle,muhstore_crust_mantle, &
- eta_anisostore_crust_mantle,idoubling_crust_mantle, &
+ eta_anisostore_crust_mantle,ispec_is_tiso_crust_mantle, &
+ ! --idoubling_crust_mantle, &
LOCAL_PATH)
!<YANGL
More information about the CIG-COMMITS
mailing list