[cig-commits] r20531 - in seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src: create_header_file meshfem3D shared specfem3D
joseph.charles at geodynamics.org
joseph.charles at geodynamics.org
Fri Jul 20 08:52:52 PDT 2012
Author: joseph.charles
Date: 2012-07-20 08:52:51 -0700 (Fri, 20 Jul 2012)
New Revision: 20531
Modified:
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/create_header_file/create_header_file.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/compute_element_properties.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_central_cube.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_doubling_elements.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regions_mesh.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regular_elements.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/finalize_mesher.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/get_model.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/memory_eval.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_forward_arrays.f90
seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_forward_arrays.f90
Log:
sets automatically the size of attenuation arrays according to the model defined in Par_file
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/create_header_file/create_header_file.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/create_header_file/create_header_file.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/create_header_file/create_header_file.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -157,7 +157,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,TOPOGRAPHY, &
- ONE_CRUST,doubling_index,this_region_has_a_doubling, &
+ ONE_CRUST,doubling_index,this_region_has_a_doubling,NCHUNKS, &
ner,NEX_PER_PROC_XI,NEX_PER_PROC_ETA, &
ratio_sampling_array,NPROCTOT, &
NSPEC,nglob,SIMULATION_TYPE,MOVIE_VOLUME,SAVE_FORWARD, &
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/compute_element_properties.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/compute_element_properties.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/compute_element_properties.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -38,7 +38,7 @@
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,ACTUALLY_STORE_ARRAYS,&
+ vx,vy,vz,rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
xigll,yigll,zigll,ispec_is_tiso)
use meshfem3D_models_par
@@ -97,10 +97,10 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,nspec_stacey) :: rho_vp,rho_vs
! attenuation
- integer nspec_att
- double precision, dimension(NGLLX,NGLLY,NGLLZ,nspec_att) :: Qmu_store
- double precision, dimension(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att) :: tau_e_store
- double precision, dimension(N_SLS) :: tau_s
+ integer vx,vy,vz,nspec_att
+ double precision, dimension(vx,vy,vz,nspec_att) :: Qmu_store
+ double precision, dimension(N_SLS,vx,vy,vz,nspec_att) :: tau_e_store
+ double precision, dimension(N_SLS) :: tau_s
double precision T_c_source
! Parameters used to calculate Jacobian based upon 125 GLL points
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_central_cube.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_central_cube.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_central_cube.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -39,7 +39,7 @@
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, &
+ nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,vx,vy,vz, &
rho_vp,rho_vs,ABSORBING_CONDITIONS,ACTUALLY_STORE_ARRAYS,xigll,yigll,zigll, &
ispec_is_tiso)
@@ -113,10 +113,10 @@
integer iproc_xi,iproc_eta,ichunk,ipass
! attenuation
- integer nspec_att
- double precision, dimension(NGLLX,NGLLY,NGLLZ,nspec_att) :: Qmu_store
- double precision, dimension(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att) :: tau_e_store
- double precision, dimension(N_SLS) :: tau_s
+ integer vx,vy,vz,nspec_att
+ double precision, dimension(vx,vy,vz,nspec_att) :: Qmu_store
+ double precision, dimension(N_SLS,vx,vy,vz,nspec_att) :: tau_e_store
+ double precision, dimension(N_SLS) :: tau_s
double precision T_c_source
logical :: ACTUALLY_STORE_ARRAYS,ABSORBING_CONDITIONS
@@ -268,8 +268,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,ACTUALLY_STORE_ARRAYS,&
+ nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4), &
+ rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS, &
xigll,yigll,zigll,ispec_is_tiso)
enddo
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_doubling_elements.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_doubling_elements.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_doubling_elements.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -44,7 +44,7 @@
gammaxstore,gammaystore,gammazstore,&
nspec_stacey,rho_vp,rho_vs,iboun,iMPIcut_xi,iMPIcut_eta, &
ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD,iproc_xi,iproc_eta, &
- nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,vx,vy,vz, &
rotation_matrix,idoubling,doubling_index,USE_ONE_LAYER_SB,ACTUALLY_STORE_ARRAYS, &
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, &
@@ -124,10 +124,10 @@
integer iproc_xi,iproc_eta
! attenuation
- integer nspec_att
- double precision, dimension(NGLLX,NGLLY,NGLLZ,nspec_att) :: Qmu_store
- double precision, dimension(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att) :: tau_e_store
- double precision, dimension(N_SLS) :: tau_s
+ integer vx,vy,vz,nspec_att
+ double precision, dimension(vx,vy,vz,nspec_att) :: Qmu_store
+ double precision, dimension(N_SLS,vx,vy,vz,nspec_att) :: tau_e_store
+ double precision, dimension(N_SLS) :: tau_s
double precision T_c_source
! rotation matrix from Euler angles
@@ -345,8 +345,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,ACTUALLY_STORE_ARRAYS,&
+ nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4), &
+ rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS, &
xigll,yigll,zigll,ispec_is_tiso)
! boundary mesh
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regions_mesh.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regions_mesh.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regions_mesh.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -251,7 +251,7 @@
! create the name for the database of the current slide and region
call create_name_database(prname,myrank,iregion_code,LOCAL_PATH)
-! initializes arrays
+ ! initializes arrays
call sync_all()
if( myrank == 0) then
write(IMAIN,*) ' ...allocating arrays '
@@ -555,8 +555,8 @@
NCHUNKS
use meshfem3D_models_par,only: &
- ATTENUATION,ANISOTROPIC_INNER_CORE,ANISOTROPIC_3D_MANTLE,SAVE_BOUNDARY_MESH, &
- AM_V
+ ATTENUATION,ATTENUATION_3D,ANISOTROPIC_INNER_CORE,ANISOTROPIC_3D_MANTLE, &
+ SAVE_BOUNDARY_MESH,AM_V
use create_regions_mesh_par2
@@ -579,8 +579,18 @@
else
nspec_att = 1
end if
- allocate(Qmu_store(NGLLX,NGLLY,NGLLZ,nspec_att), &
- tau_e_store(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att),stat=ier)
+
+ if (ATTENUATION) then
+ if (ATTENUATION_3D) then
+ ! attenuation arrays are fully 3D
+ allocate(Qmu_store(NGLLX,NGLLY,NGLLZ,nspec_att), &
+ tau_e_store(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att),stat=ier)
+ else if (.not. ATTENUATION_3D) then
+ ! save some memory in case of 1D attenuation models
+ allocate(Qmu_store(1,1,1,nspec_att), &
+ tau_e_store(N_SLS,1,1,1,nspec_att),stat=ier)
+ endif
+ endif
if(ier /= 0) stop 'error in allocate 1'
! array with model density
@@ -958,6 +968,7 @@
nspec_stacey,rho_vp,rho_vs,iboun,iMPIcut_xi,iMPIcut_eta, &
ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD,iproc_xi,iproc_eta, &
nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4), &
rotation_matrix,idoubling,doubling_index,USE_ONE_LAYER_SB, &
stretch_tab,ACTUALLY_STORE_ARRAYS, &
NSPEC2D_MOHO,NSPEC2D_400,NSPEC2D_670,nex_eta_moho, &
@@ -990,6 +1001,7 @@
nspec_stacey,rho_vp,rho_vs,iboun,iMPIcut_xi,iMPIcut_eta, &
ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD,iproc_xi,iproc_eta, &
nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4), &
rotation_matrix,idoubling,doubling_index,USE_ONE_LAYER_SB,ACTUALLY_STORE_ARRAYS, &
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, &
@@ -1036,6 +1048,7 @@
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,&
+ size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4), &
rho_vp,rho_vs,ABSORBING_CONDITIONS,ACTUALLY_STORE_ARRAYS,xigll,yigll,zigll, &
ispec_is_tiso)
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regular_elements.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regular_elements.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/create_regular_elements.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -45,7 +45,7 @@
gammaxstore,gammaystore,gammazstore,&
nspec_stacey,rho_vp,rho_vs,iboun,iMPIcut_xi,iMPIcut_eta, &
ANGULAR_WIDTH_XI_RAD,ANGULAR_WIDTH_ETA_RAD,iproc_xi,iproc_eta, &
- nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,vx,vy,vz, &
rotation_matrix,idoubling,doubling_index,USE_ONE_LAYER_SB, &
stretch_tab,ACTUALLY_STORE_ARRAYS, &
NSPEC2D_MOHO,NSPEC2D_400,NSPEC2D_670,nex_eta_moho, &
@@ -129,10 +129,10 @@
integer iproc_xi,iproc_eta
! attenuation
- integer nspec_att
- double precision, dimension(NGLLX,NGLLY,NGLLZ,nspec_att) :: Qmu_store
- double precision, dimension(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att) :: tau_e_store
- double precision, dimension(N_SLS) :: tau_s
+ integer vx,vy,vz,nspec_att
+ double precision, dimension(vx,vy,vz,nspec_att) :: Qmu_store
+ double precision, dimension(N_SLS,vx,vy,vz,nspec_att) :: tau_e_store
+ double precision, dimension(N_SLS) :: tau_s
double precision T_c_source
! rotation matrix from Euler angles
@@ -264,7 +264,8 @@
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,&
+ nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source, &
+ size(tau_e_store,2),size(tau_e_store,3),size(tau_e_store,4), &
rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
xigll,yigll,zigll,ispec_is_tiso)
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/finalize_mesher.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/finalize_mesher.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/finalize_mesher.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -100,7 +100,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,TOPOGRAPHY, &
- ONE_CRUST,doubling_index,this_region_has_a_doubling, &
+ ONE_CRUST,doubling_index,this_region_has_a_doubling,NCHUNKS, &
ner,NEX_PER_PROC_XI,NEX_PER_PROC_ETA, &
ratio_sampling_array,NPROCTOT, &
NSPEC,nglob,SIMULATION_TYPE,MOVIE_VOLUME,SAVE_FORWARD, &
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/get_model.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/get_model.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/get_model.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -312,9 +312,7 @@
endif !CUSTOM_REAL
- !> Hejun
- ! No matter 1D or 3D attenuation, we save all gll point values
- if(ATTENUATION) then
+ if(ATTENUATION .and. ATTENUATION_3D) then
tau_e_store(:,i,j,k,ispec) = tau_e(:)
Qmu_store(i,j,k,ispec) = Qmu
endif
@@ -323,6 +321,11 @@
enddo
enddo
+ if (ATTENUATION .and. .not. ATTENUATION_3D) then
+ tau_e_store(:,1,1,1,ispec) = tau_e(:)
+ Qmu_store(1,1,1,ispec) = Qmu
+ endif
+
end subroutine get_model
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/meshfem3D/save_arrays_solver.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -370,7 +370,6 @@
close(27)
- ! No matter 1D or 3D Attenuation, we save value for gll points
if(ATTENUATION) then
open(unit=27, file=prname(1:len_trim(prname))//'attenuation.bin', &
status='unknown', form='unformatted',action='write',iostat=ier)
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/memory_eval.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/memory_eval.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/memory_eval.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -27,9 +27,9 @@
! compute the approximate amount of static memory needed to run the solver
- subroutine memory_eval(OCEANS,ABSORBING_CONDITIONS,ATTENUATION,ANISOTROPIC_3D_MANTLE,&
- TRANSVERSE_ISOTROPY,ANISOTROPIC_INNER_CORE,ROTATION,TOPOGRAPHY,&
- ONE_CRUST,doubling_index,this_region_has_a_doubling,&
+ subroutine memory_eval(OCEANS,ABSORBING_CONDITIONS,ATTENUATION,ANISOTROPIC_3D_MANTLE, &
+ TRANSVERSE_ISOTROPY,ANISOTROPIC_INNER_CORE,ROTATION,TOPOGRAPHY, &
+ ONE_CRUST,doubling_index,this_region_has_a_doubling,NCHUNKS, &
ner,NEX_PER_PROC_XI,NEX_PER_PROC_ETA, &
ratio_sampling_array,NPROCTOT, &
NSPEC,nglob,SIMULATION_TYPE,MOVIE_VOLUME,SAVE_FORWARD, &
@@ -54,7 +54,7 @@
! input
logical, intent(in) :: TRANSVERSE_ISOTROPY,ANISOTROPIC_3D_MANTLE,ANISOTROPIC_INNER_CORE, &
- ROTATION,TOPOGRAPHY, &
+ ROTATION,TOPOGRAPHY,NCHUNKS, &
ATTENUATION,ONE_CRUST,OCEANS,ABSORBING_CONDITIONS, &
MOVIE_VOLUME,SAVE_FORWARD
integer, dimension(MAX_NUM_REGIONS), intent(in) :: NSPEC, nglob, &
@@ -236,8 +236,15 @@
9.d0*dble(NGLLX)*dble(NGLLY)*dble(NGLLZ)*NSPEC(IREGION_CRUST_MANTLE)*dble(CUSTOM_REAL)
! xstore_crust_mantle,ystore_crust_mantle,zstore_crust_mantle,rmass_crust_mantle
- static_memory_size = static_memory_size + &
- 4.d0*nglob(IREGION_CRUST_MANTLE)*dble(CUSTOM_REAL)
+ if(NCHUNKS /= 6 .and. ABSORBING_CONDITIONS) then
+ ! three mass matrices for the crust and mantle region: rmassx, rmassy and rmassz
+ static_memory_size = static_memory_size + &
+ 6.d0*nglob(IREGION_CRUST_MANTLE)*dble(CUSTOM_REAL)
+ else
+ ! one only keeps one mass matrix for the calculations: rmassz
+ static_memory_size = static_memory_size + &
+ 4.d0*nglob(IREGION_CRUST_MANTLE)*dble(CUSTOM_REAL)
+ endif
! rhostore_crust_mantle,kappavstore_crust_mantle,muvstore_crust_mantle
static_memory_size = static_memory_size + &
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/shared/save_header_file.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -450,17 +450,23 @@
write(IOUT,*) 'integer, parameter :: NCORNERSCHUNKS_VAL = ',NCORNERSCHUNKS
if(ATTENUATION) then
- att1 = NGLLX
- att2 = NGLLY
- att3 = NGLLZ
- att4 = NSPEC(IREGION_CRUST_MANTLE)
- att5 = NSPEC(IREGION_INNER_CORE)
+ if(ATTENUATION_3D) then
+ att1 = NGLLX
+ att2 = NGLLY
+ att3 = NGLLZ
+ else
+ att1 = 1
+ att2 = 1
+ att3 = 1
+ endif
+ att4 = NSPEC(IREGION_CRUST_MANTLE)
+ att5 = NSPEC(IREGION_INNER_CORE)
else
- att1 = 1
- att2 = 1
- att3 = 1
- att4 = 1
- att5 = 1
+ att1 = 1
+ att2 = 1
+ att3 = 1
+ att4 = 1
+ att5 = 1
endif
write(IOUT,*) 'integer, parameter :: ATT1 = ',att1
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_element.F90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -80,20 +80,19 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPECMAX_ISO_MANTLE) :: &
kappavstore,muvstore
+ ! variable sized array variables
+ integer :: vx,vy,vz,vnspec
! attenuation
! memory variables for attenuation
! memory variables R_ij are stored at the local rather than global level
! to allow for optimization of cache access by compiler
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_memory
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_xx,R_yy,R_xy,R_xz,R_yz
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
-
real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilon_trace_over_3
- integer :: vx,vy,vz,vnspec
- real(kind=CUSTOM_REAL), dimension(vx, vy, vz, vnspec) :: one_minus_sum_beta
+ real(kind=CUSTOM_REAL), dimension(vx,vy,vz,vnspec) :: one_minus_sum_beta
! gravity
double precision, dimension(NRAD_GRAVITY) :: minus_gravity_table,density_table,minus_deriv_gravity_table
@@ -235,8 +234,10 @@
endif
! precompute terms for attenuation if needed
- if(ATTENUATION_VAL) then
- one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(1,1,1,ispec)
endif
!
@@ -463,21 +464,20 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPECMAX_ISO_MANTLE) :: &
kappavstore,muvstore
-
! attenuation
! memory variables for attenuation
! memory variables R_ij are stored at the local rather than global level
! to allow for optimization of cache access by compiler
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_memory
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_xx,R_yy,R_xy,R_xz,R_yz
real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilon_trace_over_3
- integer vx,vy,vz,vnspec
+ ! variable sized array variables
+ integer :: vx,vy,vz,vnspec
! [alpha,beta,gamma]val reduced to N_SLS to N_SLS*NUM_NODES
- real(kind=CUSTOM_REAL), dimension(vx, vy, vz, vnspec) :: one_minus_sum_beta
+ real(kind=CUSTOM_REAL), dimension(vx,vy,vz,vnspec) :: one_minus_sum_beta
! gravity
@@ -634,10 +634,13 @@
endif
! precompute terms for attenuation if needed
- if(ATTENUATION_VAL) then
- one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(1,1,1,ispec)
endif
+
!
! compute either isotropic or anisotropic elements
!
@@ -1059,15 +1062,15 @@
! memory variables R_ij are stored at the local rather than global level
! to allow for optimization of cache access by compiler
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_memory
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_xx,R_yy,R_xy,R_xz,R_yz
real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilon_trace_over_3
- integer vx,vy,vz,vnspec
+ ! variable sized array variables
+ integer :: vx,vy,vz,vnspec
! [alpha,beta,gamma]val reduced to N_SLS to N_SLS*NUM_NODES
- real(kind=CUSTOM_REAL), dimension(vx, vy, vz, vnspec) :: one_minus_sum_beta
+ real(kind=CUSTOM_REAL), dimension(vx,vy,vz,vnspec) :: one_minus_sum_beta
! gravity
double precision, dimension(NRAD_GRAVITY) :: minus_gravity_table,density_table,minus_deriv_gravity_table
@@ -1211,9 +1214,12 @@
endif
! precompute terms for attenuation if needed
- if(ATTENUATION_VAL) then
- one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
- minus_sum_beta = one_minus_sum_beta_use - 1.0
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
+ minus_sum_beta = one_minus_sum_beta_use - 1.0
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(1,1,1,ispec)
+ minus_sum_beta = one_minus_sum_beta_use - 1.0
endif
!
@@ -1461,62 +1467,62 @@
imodulo_N_SLS = mod(N_SLS,3)
if(imodulo_N_SLS >= 1) then
- do i_SLS = 1,imodulo_N_SLS
- R_xx_val1 = R_xx_loc(i_SLS)
- R_yy_val1 = R_yy_loc(i_SLS)
- sigma_xx = sigma_xx - R_xx_val1
- sigma_yy = sigma_yy - R_yy_val1
- sigma_zz = sigma_zz + R_xx_val1 + R_yy_val1
- sigma_xy = sigma_xy - R_xy_loc(i_SLS)
- sigma_xz = sigma_xz - R_xz_loc(i_SLS)
- sigma_yz = sigma_yz - R_yz_loc(i_SLS)
- enddo
+ do i_SLS = 1,imodulo_N_SLS
+ R_xx_val1 = R_xx_loc(i_SLS)
+ R_yy_val1 = R_yy_loc(i_SLS)
+ sigma_xx = sigma_xx - R_xx_val1
+ sigma_yy = sigma_yy - R_yy_val1
+ sigma_zz = sigma_zz + R_xx_val1 + R_yy_val1
+ sigma_xy = sigma_xy - R_xy_loc(i_SLS)
+ sigma_xz = sigma_xz - R_xz_loc(i_SLS)
+ sigma_yz = sigma_yz - R_yz_loc(i_SLS)
+ enddo
endif
if(N_SLS >= imodulo_N_SLS+1) then
- ! note: another possibility would be using a reduction example for this loop; was tested but it does not improve,
- ! probably since N_SLS == 3 is too small for a loop benefit
- do i_SLS = imodulo_N_SLS+1,N_SLS,3
- R_xx_val1 = R_xx_loc(i_SLS)
- R_yy_val1 = R_yy_loc(i_SLS)
+ ! note: another possibility would be using a reduction example for this loop; was tested but it does not improve,
+ ! probably since N_SLS == 3 is too small for a loop benefit
+ do i_SLS = imodulo_N_SLS+1,N_SLS,3
+ R_xx_val1 = R_xx_loc(i_SLS)
+ R_yy_val1 = R_yy_loc(i_SLS)
- i_SLS1=i_SLS+1
- R_xx_val2 = R_xx_loc(i_SLS1)
- R_yy_val2 = R_yy_loc(i_SLS1)
+ i_SLS1=i_SLS+1
+ R_xx_val2 = R_xx_loc(i_SLS1)
+ R_yy_val2 = R_yy_loc(i_SLS1)
- i_SLS2 =i_SLS+2
- R_xx_val3 = R_xx_loc(i_SLS2)
- R_yy_val3 = R_yy_loc(i_SLS2)
+ i_SLS2 =i_SLS+2
+ R_xx_val3 = R_xx_loc(i_SLS2)
+ R_yy_val3 = R_yy_loc(i_SLS2)
- sigma_xx = sigma_xx - R_xx_val1 - R_xx_val2 - R_xx_val3
- sigma_yy = sigma_yy - R_yy_val1 - R_yy_val2 - R_yy_val3
- sigma_zz = sigma_zz + R_xx_val1 + R_yy_val1 &
- + R_xx_val2 + R_yy_val2 &
- + R_xx_val3 + R_yy_val3
+ sigma_xx = sigma_xx - R_xx_val1 - R_xx_val2 - R_xx_val3
+ sigma_yy = sigma_yy - R_yy_val1 - R_yy_val2 - R_yy_val3
+ sigma_zz = sigma_zz + R_xx_val1 + R_yy_val1 &
+ + R_xx_val2 + R_yy_val2 &
+ + R_xx_val3 + R_yy_val3
- sigma_xy = sigma_xy - R_xy_loc(i_SLS)
- sigma_xz = sigma_xz - R_xz_loc(i_SLS)
- sigma_yz = sigma_yz - R_yz_loc(i_SLS)
+ sigma_xy = sigma_xy - R_xy_loc(i_SLS)
+ sigma_xz = sigma_xz - R_xz_loc(i_SLS)
+ sigma_yz = sigma_yz - R_yz_loc(i_SLS)
- sigma_xy = sigma_xy - R_xy_loc(i_SLS1)
- sigma_xz = sigma_xz - R_xz_loc(i_SLS1)
- sigma_yz = sigma_yz - R_yz_loc(i_SLS1)
+ sigma_xy = sigma_xy - R_xy_loc(i_SLS1)
+ sigma_xz = sigma_xz - R_xz_loc(i_SLS1)
+ sigma_yz = sigma_yz - R_yz_loc(i_SLS1)
- sigma_xy = sigma_xy - R_xy_loc(i_SLS2)
- sigma_xz = sigma_xz - R_xz_loc(i_SLS2)
- sigma_yz = sigma_yz - R_yz_loc(i_SLS2)
- enddo
+ sigma_xy = sigma_xy - R_xy_loc(i_SLS2)
+ sigma_xz = sigma_xz - R_xz_loc(i_SLS2)
+ sigma_yz = sigma_yz - R_yz_loc(i_SLS2)
+ enddo
endif
#else
! way 1:
do i_SLS = 1,N_SLS
- R_xx_val1 = R_xx_loc(i_SLS) ! R_memory(1,i_SLS,i,j,k,ispec)
- R_yy_val1 = R_yy_loc(i_SLS) ! R_memory(2,i_SLS,i,j,k,ispec)
- sigma_xx = sigma_xx - R_xx_val1
- sigma_yy = sigma_yy - R_yy_val1
- sigma_zz = sigma_zz + R_xx_val1 + R_yy_val1
- sigma_xy = sigma_xy - R_xy_loc(i_SLS) ! R_memory(3,i_SLS,i,j,k,ispec)
- sigma_xz = sigma_xz - R_xz_loc(i_SLS) ! R_memory(4,i_SLS,i,j,k,ispec)
- sigma_yz = sigma_yz - R_yz_loc(i_SLS) ! R_memory(5,i_SLS,i,j,k,ispec)
+ R_xx_val1 = R_xx_loc(i_SLS) ! R_memory(1,i_SLS,i,j,k,ispec)
+ R_yy_val1 = R_yy_loc(i_SLS) ! R_memory(2,i_SLS,i,j,k,ispec)
+ sigma_xx = sigma_xx - R_xx_val1
+ sigma_yy = sigma_yy - R_yy_val1
+ sigma_zz = sigma_zz + R_xx_val1 + R_yy_val1
+ sigma_xy = sigma_xy - R_xy_loc(i_SLS) ! R_memory(3,i_SLS,i,j,k,ispec)
+ sigma_xz = sigma_xz - R_xz_loc(i_SLS) ! R_memory(4,i_SLS,i,j,k,ispec)
+ sigma_yz = sigma_yz - R_yz_loc(i_SLS) ! R_memory(5,i_SLS,i,j,k,ispec)
enddo
#endif
@@ -1562,12 +1568,17 @@
! element id
integer :: ispec
+ ! attenuation
+ ! memory variables for attenuation
+ ! memory variables R_ij are stored at the local rather than global level
+ ! to allow for optimization of cache access by compiler
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_memory
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_xx,R_yy,R_xy,R_xz,R_yz
+ ! variable sized array variables
integer :: vx,vy,vz,vnspec
- real(kind=CUSTOM_REAL), dimension(N_SLS, vx, vy, vz, vnspec) :: factor_common
+
+ real(kind=CUSTOM_REAL), dimension(N_SLS,vx,vy,vz,vnspec) :: factor_common
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPECMAX_ANISO_MANTLE) :: c44store
@@ -1580,7 +1591,7 @@
real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ) :: epsilondev_loc
! local parameters
- real(kind=CUSTOM_REAL), dimension(NGLLX, NGLLY, NGLLZ) :: factor_common_c44_muv
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: factor_common_c44_muv
integer :: i_SLS
#ifdef _HANDOPT_ATT
@@ -1603,86 +1614,86 @@
gammal = gammaval(i_SLS)
! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
- factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec)
-
- if(ANISOTROPIC_3D_MANTLE_VAL) then
- factor_common_c44_muv(:,:,:) = factor_common_c44_muv(:,:,:) * c44store(:,:,:,ispec)
- else
- factor_common_c44_muv(:,:,:) = factor_common_c44_muv(:,:,:) * muvstore(:,:,:,ispec)
+ if (ATTENUATION_3D_VAL) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * c44store(:,:,:,ispec)
+ else
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * muvstore(:,:,:,ispec)
+ endif
+ else if (.not. ATTENUATION_3D_VAL) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * c44store(:,:,:,ispec)
+ else
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * muvstore(:,:,:,ispec)
+ endif
endif
! this helps to vectorize the inner most loop
do k=1,NGLLZ
- do j=1,NGLLY
- do i=1,NGLLX
-! R_memory(:,i_SLS,i,j,k,ispec) = alphal * R_memory(:,i_SLS,i,j,k,ispec) &
-! + factor_common_c44_muv(i,j,k) &
-! *( betal * epsilondev(:,i,j,k,ispec) + gammal * epsilondev_loc(:,i,j,k))
+ do j=1,NGLLY
+ do i=1,NGLLX
+ ! R_memory(:,i_SLS,i,j,k,ispec) = alphal * R_memory(:,i_SLS,i,j,k,ispec) &
+ ! + factor_common_c44_muv(i,j,k) &
+ ! *( betal * epsilondev(:,i,j,k,ispec) + gammal * epsilondev_loc(:,i,j,k))
- R_xx(i_SLS,i,j,k,ispec) = alphal * R_xx(i_SLS,i,j,k,ispec) &
- + factor_common_c44_muv(i,j,k) &
- *( betal * epsilondev_xx(i,j,k,ispec) + gammal * epsilondev_loc(1,i,j,k))
+ R_xx(i_SLS,i,j,k,ispec) = alphal * R_xx(i_SLS,i,j,k,ispec) + factor_common_c44_muv(i,j,k) * &
+ (betal * epsilondev_xx(i,j,k,ispec) + gammal * epsilondev_loc(1,i,j,k))
- R_yy(i_SLS,i,j,k,ispec) = alphal * R_yy(i_SLS,i,j,k,ispec) &
- + factor_common_c44_muv(i,j,k) &
- *( betal * epsilondev_yy(i,j,k,ispec) + gammal * epsilondev_loc(2,i,j,k))
+ R_yy(i_SLS,i,j,k,ispec) = alphal * R_yy(i_SLS,i,j,k,ispec) + factor_common_c44_muv(i,j,k) * &
+ (betal * epsilondev_yy(i,j,k,ispec) + gammal * epsilondev_loc(2,i,j,k))
- R_xy(i_SLS,i,j,k,ispec) = alphal * R_xy(i_SLS,i,j,k,ispec) &
- + factor_common_c44_muv(i,j,k) &
- *( betal * epsilondev_xy(i,j,k,ispec) + gammal * epsilondev_loc(3,i,j,k))
+ R_xy(i_SLS,i,j,k,ispec) = alphal * R_xy(i_SLS,i,j,k,ispec) + factor_common_c44_muv(i,j,k) * &
+ (betal * epsilondev_xy(i,j,k,ispec) + gammal * epsilondev_loc(3,i,j,k))
- R_xz(i_SLS,i,j,k,ispec) = alphal * R_xz(i_SLS,i,j,k,ispec) &
- + factor_common_c44_muv(i,j,k) &
- *( betal * epsilondev_xz(i,j,k,ispec) + gammal * epsilondev_loc(4,i,j,k))
+ R_xz(i_SLS,i,j,k,ispec) = alphal * R_xz(i_SLS,i,j,k,ispec) + factor_common_c44_muv(i,j,k) * &
+ (betal * epsilondev_xz(i,j,k,ispec) + gammal * epsilondev_loc(4,i,j,k))
- R_yz(i_SLS,i,j,k,ispec) = alphal * R_yz(i_SLS,i,j,k,ispec) &
- + factor_common_c44_muv(i,j,k) &
- *( betal * epsilondev_yz(i,j,k,ispec) + gammal * epsilondev_loc(5,i,j,k))
+ R_yz(i_SLS,i,j,k,ispec) = alphal * R_yz(i_SLS,i,j,k,ispec) + factor_common_c44_muv(i,j,k) * &
+ (betal * epsilondev_yz(i,j,k,ispec) + gammal * epsilondev_loc(5,i,j,k))
-
- enddo
- enddo
+ enddo
+ enddo
enddo
enddo ! i_SLS
#else
! way 1:
do i_SLS = 1,N_SLS
- ! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
- factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec)
+ ! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
+ if (ATTENUATION_3D_VAL) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * c44store(:,:,:,ispec)
+ else
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * muvstore(:,:,:,ispec)
+ endif
+ else if (.not. ATTENUATION_3D_VAL) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * c44store(:,:,:,ispec)
+ else
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * muvstore(:,:,:,ispec)
+ endif
+ endif
- if(ANISOTROPIC_3D_MANTLE_VAL) then
- factor_common_c44_muv(:,:,:) = factor_common_c44_muv(:,:,:) * c44store(:,:,:,ispec)
- else
- factor_common_c44_muv(:,:,:) = factor_common_c44_muv(:,:,:) * muvstore(:,:,:,ispec)
- endif
+ ! do i_memory = 1,5
+ ! R_memory(i_memory,i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_memory(i_memory,i_SLS,:,:,:,ispec) &
+ ! + factor_common_c44_muv(:,:,:) &
+ ! * (betaval(i_SLS) * epsilondev(i_memory,:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(i_memory,:,:,:))
+ ! enddo
-! do i_memory = 1,5
-! R_memory(i_memory,i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_memory(i_memory,i_SLS,:,:,:,ispec) &
-! + factor_common_c44_muv(:,:,:) &
-! * (betaval(i_SLS) * epsilondev(i_memory,:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(i_memory,:,:,:))
-! enddo
+ R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
- R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) &
- + factor_common_c44_muv(:,:,:) &
- * (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
+ R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
- R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) &
- + factor_common_c44_muv(:,:,:) &
- * (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
+ R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
- R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) &
- + factor_common_c44_muv(:,:,:) &
- * (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
+ R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
- R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) &
- + factor_common_c44_muv(:,:,:) &
- * (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
+ R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
- R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) &
- + factor_common_c44_muv(:,:,:) &
- * (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
-
-
enddo ! i_SLS
#endif
@@ -1733,8 +1744,10 @@
real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE_ATTENUATION) :: &
R_xx,R_yy,R_xy,R_xz,R_yz
+ ! variable sized array variables
integer :: vx,vy,vz,vnspec
- real(kind=CUSTOM_REAL), dimension(N_SLS, vx, vy, vz, vnspec) :: factor_common
+
+ real(kind=CUSTOM_REAL), dimension(N_SLS,vx,vy,vz,vnspec) :: factor_common
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: muvstore
@@ -1746,7 +1759,7 @@
real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ) :: epsilondev_loc
! local parameters
- real(kind=CUSTOM_REAL), dimension(NGLLX, NGLLY, NGLLZ) :: factor_common_use
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: factor_common_use
integer :: i_SLS
@@ -1770,68 +1783,67 @@
gammal = gammaval(i_SLS)
! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
- factor_common_use(:,:,:) = factor_common(i_SLS,:,:,:,ispec)
+ if (ATTENUATION_3D_VAL) then
+ factor_common_use(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * muvstore(:,:,:,ispec)
+ else if (.not. ATTENUATION_3D_VAL) then
+ factor_common_use(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * muvstore(:,:,:,ispec)
+ endif
- factor_common_use(:,:,:) = factor_common_use(:,:,:) * muvstore(:,:,:,ispec)
-
! this helps to vectorize the inner most loop
do k=1,NGLLZ
- do j=1,NGLLY
- do i=1,NGLLX
-! R_memory(:,i_SLS,i,j,k,ispec) = alphal * R_memory(:,i_SLS,i,j,k,ispec) &
-! + factor_common_use(i,j,k) &
-! *( betal * epsilondev(:,i,j,k,ispec) + gammal * epsilondev_loc(:,i,j,k))
+ do j=1,NGLLY
+ do i=1,NGLLX
+ ! R_memory(:,i_SLS,i,j,k,ispec) = alphal * R_memory(:,i_SLS,i,j,k,ispec) &
+ ! + factor_common_use(i,j,k) &
+ ! *( betal * epsilondev(:,i,j,k,ispec) + gammal * epsilondev_loc(:,i,j,k))
- R_xx(i_SLS,i,j,k,ispec) = alphal * R_xx(i_SLS,i,j,k,ispec) &
- + factor_common_use(i,j,k) &
- *( betal * epsilondev_xx(i,j,k,ispec) + gammal * epsilondev_loc(1,i,j,k))
- R_yy(i_SLS,i,j,k,ispec) = alphal * R_yy(i_SLS,i,j,k,ispec) &
- + factor_common_use(i,j,k) &
- *( betal * epsilondev_yy(i,j,k,ispec) + gammal * epsilondev_loc(2,i,j,k))
- R_xy(i_SLS,i,j,k,ispec) = alphal * R_xy(i_SLS,i,j,k,ispec) &
- + factor_common_use(i,j,k) &
- *( betal * epsilondev_xy(i,j,k,ispec) + gammal * epsilondev_loc(3,i,j,k))
- R_xz(i_SLS,i,j,k,ispec) = alphal * R_xz(i_SLS,i,j,k,ispec) &
- + factor_common_use(i,j,k) &
- *( betal * epsilondev_xz(i,j,k,ispec) + gammal * epsilondev_loc(4,i,j,k))
- R_yz(i_SLS,i,j,k,ispec) = alphal * R_yz(i_SLS,i,j,k,ispec) &
- + factor_common_use(i,j,k) &
- *( betal * epsilondev_yz(i,j,k,ispec) + gammal * epsilondev_loc(5,i,j,k))
+ R_xx(i_SLS,i,j,k,ispec) = alphal * R_xx(i_SLS,i,j,k,ispec) + factor_common_use(i,j,k) * &
+ (betal * epsilondev_xx(i,j,k,ispec) + gammal * epsilondev_loc(1,i,j,k))
+ R_yy(i_SLS,i,j,k,ispec) = alphal * R_yy(i_SLS,i,j,k,ispec) + factor_common_use(i,j,k) * &
+ (betal * epsilondev_yy(i,j,k,ispec) + gammal * epsilondev_loc(2,i,j,k))
+ R_xy(i_SLS,i,j,k,ispec) = alphal * R_xy(i_SLS,i,j,k,ispec) + factor_common_use(i,j,k) * &
+ (betal * epsilondev_xy(i,j,k,ispec) + gammal * epsilondev_loc(3,i,j,k))
+ R_xz(i_SLS,i,j,k,ispec) = alphal * R_xz(i_SLS,i,j,k,ispec) + factor_common_use(i,j,k) * &
+ (betal * epsilondev_xz(i,j,k,ispec) + gammal * epsilondev_loc(4,i,j,k))
+ R_yz(i_SLS,i,j,k,ispec) = alphal * R_yz(i_SLS,i,j,k,ispec) + factor_common_use(i,j,k) * &
+ (betal * epsilondev_yz(i,j,k,ispec) + gammal * epsilondev_loc(5,i,j,k))
- enddo
- enddo
+ enddo
+ enddo
enddo
enddo ! i_SLS
#else
! way 1:
do i_SLS = 1,N_SLS
- factor_common_use(:,:,:) = factor_common(i_SLS,:,:,:,ispec)
+
+ ! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
+ if (ATTENUATION_3D_VAL) then
+ factor_common_use(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * muvstore(:,:,:,ispec)
+ else if (.not. ATTENUATION_3D_VAL) then
+ factor_common_use(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * muvstore(:,:,:,ispec)
+ endif
+
! do i_memory = 1,5
! R_memory(i_memory,i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_memory(i_memory,i_SLS,:,:,:,ispec) &
! + muvstore(:,:,:,ispec) * factor_common_use(:,:,:) * &
! (betaval(i_SLS) * epsilondev(i_memory,:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(i_memory,:,:,:))
! enddo
- R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use(:,:,:) * &
- (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
+ R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
- R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use(:,:,:) * &
- (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
+ R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
- R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use(:,:,:) * &
- (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
+ R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
- R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use(:,:,:) * &
- (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
+ R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
- R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use(:,:,:) * &
- (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
+ R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
enddo
#endif
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -79,24 +79,23 @@
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB) :: veloc_crust_mantle
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB) :: accel_crust_mantle
+ ! variable sized array variables
+ integer :: vx,vy,vz,vnspec
+
! memory variables for attenuation
! memory variables R_ij are stored at the local rather than global level
! to allow for optimization of cache access by compiler
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_memory
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATT) :: R_xx,R_yy,R_xy,R_xz,R_yz
-
! real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilondev
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC) :: &
epsilondev_xx,epsilondev_yy,epsilondev_xy,epsilondev_xz,epsilondev_yz
real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC) :: epsilon_trace_over_3
- ! variable sized array variables for one_minus_sum_beta and factor_common
- integer vx, vy, vz, vnspec
! [alpha,beta,gamma]val reduced to N_SLS and factor_common to N_SLS*NUM_NODES
- real(kind=CUSTOM_REAL), dimension(N_SLS, vx, vy, vz, vnspec) :: factor_common
+ real(kind=CUSTOM_REAL), dimension(N_SLS,vx,vy,vz,vnspec) :: factor_common
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
! inner/outer element run flag
@@ -107,7 +106,7 @@
! for attenuation
real(kind=CUSTOM_REAL) one_minus_sum_beta_use,minus_sum_beta
real(kind=CUSTOM_REAL) R_xx_val,R_yy_val
- real(kind=CUSTOM_REAL), dimension(NGLLX, NGLLY, NGLLZ) :: factor_common_c44_muv
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: factor_common_c44_muv
real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ) :: epsilondev_loc
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: &
@@ -363,14 +362,16 @@
endif
! precompute terms for attenuation if needed
- if(ATTENUATION_VAL) then
- one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
- minus_sum_beta = one_minus_sum_beta_use - 1.0
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ one_minus_sum_beta_use = one_minus_sum_beta(1,1,1,ispec)
endif
+ minus_sum_beta = one_minus_sum_beta_use - 1.0
- !
- ! compute either isotropic or anisotropic elements
- !
+ !
+ ! compute either isotropic or anisotropic elements
+ !
if(ANISOTROPIC_3D_MANTLE_VAL) then
@@ -663,16 +664,16 @@
! subtract memory variables if attenuation
if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. ) ) then
- do i_SLS = 1,N_SLS
- R_xx_val = R_xx(i_SLS,i,j,k,ispec)
- R_yy_val = R_yy(i_SLS,i,j,k,ispec)
- sigma_xx = sigma_xx - R_xx_val
- sigma_yy = sigma_yy - R_yy_val
- sigma_zz = sigma_zz + R_xx_val + R_yy_val
- sigma_xy = sigma_xy - R_xy(i_SLS,i,j,k,ispec)
- sigma_xz = sigma_xz - R_xz(i_SLS,i,j,k,ispec)
- sigma_yz = sigma_yz - R_yz(i_SLS,i,j,k,ispec)
- enddo
+ do i_SLS = 1,N_SLS
+ R_xx_val = R_xx(i_SLS,i,j,k,ispec)
+ R_yy_val = R_yy(i_SLS,i,j,k,ispec)
+ sigma_xx = sigma_xx - R_xx_val
+ sigma_yy = sigma_yy - R_yy_val
+ sigma_zz = sigma_zz + R_xx_val + R_yy_val
+ sigma_xy = sigma_xy - R_xy(i_SLS,i,j,k,ispec)
+ sigma_xz = sigma_xz - R_xz(i_SLS,i,j,k,ispec)
+ sigma_yz = sigma_yz - R_yz(i_SLS,i,j,k,ispec)
+ enddo
endif
! define symmetric components of sigma for gravity
@@ -891,12 +892,19 @@
! IMPROVE we should probably use an average value instead
! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
- factor_common_c44_muv = factor_common(i_SLS,:,:,:,ispec)
- if(ANISOTROPIC_3D_MANTLE_VAL) then
- factor_common_c44_muv = factor_common_c44_muv * c44store(:,:,:,ispec)
- else
- factor_common_c44_muv = factor_common_c44_muv * muvstore(:,:,:,ispec)
- endif
+ if (ATTENUATION_3D_VAL) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * c44store(:,:,:,ispec)
+ else
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * muvstore(:,:,:,ispec)
+ endif
+ else if (.not. ATTENUATION_3D_VAL) then
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * c44store(:,:,:,ispec)
+ else
+ factor_common_c44_muv(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * muvstore(:,:,:,ispec)
+ endif
+ endif
! do i_memory = 1,5
! R_memory(i_memory,i_SLS,:,:,:,ispec) = alphaval(i_SLS) * &
@@ -906,23 +914,22 @@
! gammaval(i_SLS) * epsilondev_loc(i_memory,:,:,:))
! enddo
- R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) + factor_common_c44_muv * &
- (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
+ R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
- R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) + factor_common_c44_muv * &
- (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
+ R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
- R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) + factor_common_c44_muv * &
- (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
+ R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
- R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) + factor_common_c44_muv * &
- (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
+ R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
- R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) + factor_common_c44_muv * &
- (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
+ R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) + factor_common_c44_muv(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
enddo
-
endif
! save deviatoric strain for Runge-Kutta scheme
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_crust_mantle_Dev.F90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -89,13 +89,14 @@
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB) :: veloc_crust_mantle
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB) :: accel_crust_mantle
- ! attenuation
+ ! variable sized array variables
+ integer :: vx,vy,vz,vnspec
+
! memory variables for attenuation
! memory variables R_ij are stored at the local rather than global level
! to allow for optimization of cache access by compiler
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE_ATTENUAT) :: R_memory
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATT) :: R_xx,R_yy,R_xy,R_xz,R_yz
! real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: epsilondev
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC) :: &
@@ -103,10 +104,8 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC) :: epsilon_trace_over_3
- integer :: vx,vy,vz,vnspec
-
! [alpha,beta,gamma]val reduced to N_SLS and factor_common to N_SLS*NUM_NODES
- real(kind=CUSTOM_REAL), dimension(N_SLS, vx, vy, vz, vnspec) :: factor_common
+ real(kind=CUSTOM_REAL), dimension(N_SLS,vx,vy,vz,vnspec) :: factor_common
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
! inner/outer element run flag
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -76,13 +76,15 @@
! memory variables R_ij are stored at the local rather than global level
! to allow for optimization of cache access by compiler
! variable lengths for factor_common and one_minus_sum_beta
+
+ ! variable sized array variables
integer vx, vy, vz, vnspec
- real(kind=CUSTOM_REAL), dimension(N_SLS, vx, vy, vz, vnspec) :: factor_common
+
+ real(kind=CUSTOM_REAL), dimension(N_SLS,vx,vy,vz,vnspec) :: factor_common
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
! real(kind=CUSTOM_REAL), dimension(5,N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE_ATTENUATION) :: R_memory
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATT) :: &
- R_xx,R_yy,R_xy,R_xz,R_yz
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATT) :: R_xx,R_yy,R_xy,R_xz,R_yz
! real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: epsilondev
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC) :: &
@@ -96,7 +98,7 @@
! local parameters
real(kind=CUSTOM_REAL), dimension(5,NGLLX,NGLLY,NGLLZ) :: epsilondev_loc
- real(kind=CUSTOM_REAL), dimension(NGLLX, NGLLY, NGLLZ) :: factor_common_use
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: factor_common_use
real(kind=CUSTOM_REAL) R_xx_val,R_yy_val
integer i_SLS
@@ -343,8 +345,11 @@
endif
endif
- if(ATTENUATION_VAL) then
- minus_sum_beta = one_minus_sum_beta(i,j,k,ispec) - 1.0
+ ! precompute terms for attenuation if needed
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ minus_sum_beta = one_minus_sum_beta(i,j,k,ispec) - 1.0
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ minus_sum_beta = one_minus_sum_beta(1,1,1,ispec) - 1.0
endif
if(ANISOTROPIC_INNER_CORE_VAL) then
@@ -396,8 +401,10 @@
mul = muvstore(i,j,k,ispec)
! use unrelaxed parameters if attenuation
- if(ATTENUATION_VAL) then
- mul = mul * one_minus_sum_beta(i,j,k,ispec)
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ mul = mul * one_minus_sum_beta(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ mul = mul * one_minus_sum_beta(1,1,1,ispec)
endif
lambdalplus2mul = kappal + FOUR_THIRDS * mul
@@ -644,7 +651,14 @@
if(ATTENUATION_VAL .and. ( USE_ATTENUATION_MIMIC .eqv. .false. )) then
do i_SLS = 1,N_SLS
- factor_common_use = factor_common(i_SLS,:,:,:,ispec)
+
+ ! reformatted R_memory to handle large factor_common and reduced [alpha,beta,gamma]val
+ if (ATTENUATION_3D_VAL) then
+ factor_common_use(:,:,:) = factor_common(i_SLS,:,:,:,ispec) * muvstore(:,:,:,ispec)
+ else if (.not. ATTENUATION_3D_VAL) then
+ factor_common_use(:,:,:) = factor_common(i_SLS,1,1,1,ispec) * muvstore(:,:,:,ispec)
+ endif
+
! do i_memory = 1,5
! R_memory(i_memory,i_SLS,:,:,:,ispec) = &
! alphaval(i_SLS) * &
@@ -653,28 +667,22 @@
! (betaval(i_SLS) * &
! epsilondev(i_memory,:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(i_memory,:,:,:))
! enddo
+
+ R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
- R_xx(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xx(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use * &
- (betaval(i_SLS) * epsilondev_xx(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(1,:,:,:))
+ R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
- R_yy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yy(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use * &
- (betaval(i_SLS) * epsilondev_yy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(2,:,:,:))
+ R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
- R_xy(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xy(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use * &
- (betaval(i_SLS) * epsilondev_xy(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(3,:,:,:))
+ R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
- R_xz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_xz(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use * &
- (betaval(i_SLS) * epsilondev_xz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(4,:,:,:))
+ R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) + factor_common_use(:,:,:) * &
+ (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
- R_yz(i_SLS,:,:,:,ispec) = alphaval(i_SLS) * R_yz(i_SLS,:,:,:,ispec) &
- + muvstore(:,:,:,ispec) * factor_common_use * &
- (betaval(i_SLS) * epsilondev_yz(:,:,:,ispec) + gammaval(i_SLS) * epsilondev_loc(5,:,:,:))
-
-
enddo
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/compute_forces_inner_core_Dev.F90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -616,8 +616,10 @@
endif
endif
- if(ATTENUATION_VAL) then
- minus_sum_beta = one_minus_sum_beta(i,j,k,ispec) - 1.0
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ minus_sum_beta = one_minus_sum_beta(i,j,k,ispec) - 1.0
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ minus_sum_beta = one_minus_sum_beta(1,1,1,ispec) - 1.0
endif
if(ANISOTROPIC_INNER_CORE_VAL) then
@@ -667,8 +669,10 @@
mul = muvstore(i,j,k,ispec)
! use unrelaxed parameters if attenuation
- if(ATTENUATION_VAL) then
- mul = mul * one_minus_sum_beta(i,j,k,ispec)
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ mul = mul * one_minus_sum_beta(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ mul = mul * one_minus_sum_beta(1,1,1,ispec)
endif
lambdalplus2mul = kappal + FOUR_THIRDS * mul
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/get_attenuation.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -26,18 +26,20 @@
!=====================================================================
- subroutine get_attenuation_model_3D(myrank, prname, one_minus_sum_beta, &
- factor_common, scale_factor, tau_s, vnspec)
+ subroutine get_attenuation_model_3D_or_1D(myrank, prname, one_minus_sum_beta, &
+ factor_common, scale_factor, tau_s, vx, vy, vz, vnspec)
+ use specfem_par,only: ATTENUATION_VAL, ATTENUATION_3D_VAL
+
implicit none
include 'constants.h'
- integer myrank, vnspec
+ integer myrank,vx,vy,vz,vnspec
character(len=150) prname
- double precision, dimension(NGLLX,NGLLY,NGLLZ,vnspec) :: one_minus_sum_beta, scale_factor
- double precision, dimension(N_SLS,NGLLX,NGLLY,NGLLZ,vnspec) :: factor_common
- double precision, dimension(N_SLS) :: tau_s
+ double precision, dimension(vx,vy,vz,vnspec) :: one_minus_sum_beta, scale_factor
+ double precision, dimension(N_SLS,vx,vy,vz,vnspec) :: factor_common
+ double precision, dimension(N_SLS) :: tau_s
integer i,j,k,ispec
@@ -61,28 +63,45 @@
T_c_source = 1000.0d0 / T_c_source
T_c_source = T_c_source / scale_t
- do ispec = 1, vnspec
- do k = 1, NGLLZ
- do j = 1, NGLLY
- do i = 1, NGLLX
- tau_e(:) = factor_common(:,i,j,k,ispec)
- Q_mu = scale_factor(i,j,k,ispec)
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ do ispec = 1, vnspec
+ do k = 1, NGLLZ
+ do j = 1, NGLLY
+ do i = 1, NGLLX
+ tau_e(:) = factor_common(:,i,j,k,ispec)
+ Q_mu = scale_factor(i,j,k,ispec)
- ! Determine the factor_common and one_minus_sum_beta from tau_s and tau_e
- call get_attenuation_property_values(tau_s, tau_e, fc, omsb)
+ ! Determine the factor_common and one_minus_sum_beta from tau_s and tau_e
+ call get_attenuation_property_values(tau_s, tau_e, fc, omsb)
- factor_common(:,i,j,k,ispec) = fc(:)
- one_minus_sum_beta(i,j,k,ispec) = omsb
+ factor_common(:,i,j,k,ispec) = fc(:)
+ one_minus_sum_beta(i,j,k,ispec) = omsb
- ! Determine the "scale_factor" from tau_s, tau_e, central source frequency, and Q
- call get_attenuation_scale_factor(myrank, T_c_source, tau_e, tau_s, Q_mu, sf)
- scale_factor(i,j,k,ispec) = sf
+ ! Determine the "scale_factor" from tau_s, tau_e, central source frequency, and Q
+ call get_attenuation_scale_factor(myrank, T_c_source, tau_e, tau_s, Q_mu, sf)
+ scale_factor(i,j,k,ispec) = sf
+ enddo
enddo
enddo
enddo
- enddo
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ do ispec = 1, vnspec
+ tau_e(:) = factor_common(:,1,1,1,ispec)
+ Q_mu = scale_factor(1,1,1,ispec)
+
+ ! Determine the factor_common and one_minus_sum_beta from tau_s and tau_e
+ call get_attenuation_property_values(tau_s, tau_e, fc, omsb)
+
+ factor_common(:,1,1,1,ispec) = fc(:)
+ one_minus_sum_beta(1,1,1,ispec) = omsb
+
+ ! Determine the "scale_factor" from tau_s, tau_e, central source frequency, and Q
+ call get_attenuation_scale_factor(myrank, T_c_source, tau_e, tau_s, Q_mu, sf)
+ scale_factor(1,1,1,ispec) = sf
+ enddo
+ endif
- end subroutine get_attenuation_model_3D
+ end subroutine get_attenuation_model_3D_or_1D
!
!-------------------------------------------------------------------------------------------------
@@ -104,7 +123,7 @@
beta(:) = 1.0d0 - tau_e(:) / tau_s(:)
one_minus_sum_beta = 1.0d0
- do i = 1,N_SLS
+ do i = 1, N_SLS
one_minus_sum_beta = one_minus_sum_beta - beta(i)
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/prepare_timerun.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -705,13 +705,15 @@
! CRUST_MANTLE ATTENUATION
call create_name_database(prnamel, myrank, IREGION_CRUST_MANTLE, LOCAL_PATH)
- call get_attenuation_model_3D(myrank, prnamel, omsb_crust_mantle_dble, &
- factor_common_crust_mantle_dble,factor_scale_crust_mantle_dble,tau_sigma_dble,NSPEC_CRUST_MANTLE)
+ call get_attenuation_model_3D_or_1D(myrank, prnamel, omsb_crust_mantle_dble, &
+ factor_common_crust_mantle_dble,factor_scale_crust_mantle_dble,tau_sigma_dble, &
+ ATT1,ATT2,ATT3,ATT4,NSPEC_CRUST_MANTLE)
! INNER_CORE ATTENUATION
call create_name_database(prnamel, myrank, IREGION_INNER_CORE, LOCAL_PATH)
- call get_attenuation_model_3D(myrank, prnamel, omsb_inner_core_dble, &
- factor_common_inner_core_dble,factor_scale_inner_core_dble,tau_sigma_dble,NSPEC_INNER_CORE)
+ call get_attenuation_model_3D_or_1D(myrank, prnamel, omsb_inner_core_dble, &
+ factor_common_inner_core_dble,factor_scale_inner_core_dble,tau_sigma_dble, &
+ ATT1,ATT2,ATT3,ATT5,NSPEC_INNER_CORE)
if(CUSTOM_REAL == SIZE_REAL) then
factor_scale_crust_mantle = sngl(factor_scale_crust_mantle_dble)
@@ -748,48 +750,52 @@
! rescale in crust and mantle
do ispec = 1,NSPEC_CRUST_MANTLE
- do k=1,NGLLZ
- do j=1,NGLLY
- do i=1,NGLLX
- scale_factor = factor_scale_crust_mantle(i,j,k,ispec)
+ do k=1,NGLLZ
+ do j=1,NGLLY
+ do i=1,NGLLX
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ scale_factor = factor_scale_crust_mantle(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ scale_factor = factor_scale_crust_mantle(1,1,1,ispec)
+ endif
+
+ if(ANISOTROPIC_3D_MANTLE_VAL) then
+ scale_factor_minus_one = scale_factor - 1.
+ mul = c44store_crust_mantle(i,j,k,ispec)
+ c11store_crust_mantle(i,j,k,ispec) = c11store_crust_mantle(i,j,k,ispec) &
+ + FOUR_THIRDS * scale_factor_minus_one * mul
+ c12store_crust_mantle(i,j,k,ispec) = c12store_crust_mantle(i,j,k,ispec) &
+ - TWO_THIRDS * scale_factor_minus_one * mul
+ c13store_crust_mantle(i,j,k,ispec) = c13store_crust_mantle(i,j,k,ispec) &
+ - TWO_THIRDS * scale_factor_minus_one * mul
+ c22store_crust_mantle(i,j,k,ispec) = c22store_crust_mantle(i,j,k,ispec) &
+ + FOUR_THIRDS * scale_factor_minus_one * mul
+ c23store_crust_mantle(i,j,k,ispec) = c23store_crust_mantle(i,j,k,ispec) &
+ - TWO_THIRDS * scale_factor_minus_one * mul
+ c33store_crust_mantle(i,j,k,ispec) = c33store_crust_mantle(i,j,k,ispec) &
+ + FOUR_THIRDS * scale_factor_minus_one * mul
+ c44store_crust_mantle(i,j,k,ispec) = c44store_crust_mantle(i,j,k,ispec) &
+ + scale_factor_minus_one * mul
+ c55store_crust_mantle(i,j,k,ispec) = c55store_crust_mantle(i,j,k,ispec) &
+ + scale_factor_minus_one * mul
+ c66store_crust_mantle(i,j,k,ispec) = c66store_crust_mantle(i,j,k,ispec) &
+ + scale_factor_minus_one * mul
+ else
+ if(MOVIE_VOLUME .and. SIMULATION_TYPE==3) then
+ ! store the original value of \mu to comput \mu*\eps
+ 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(ANISOTROPIC_3D_MANTLE_VAL) then
- scale_factor_minus_one = scale_factor - 1.
- mul = c44store_crust_mantle(i,j,k,ispec)
- c11store_crust_mantle(i,j,k,ispec) = c11store_crust_mantle(i,j,k,ispec) &
- + FOUR_THIRDS * scale_factor_minus_one * mul
- c12store_crust_mantle(i,j,k,ispec) = c12store_crust_mantle(i,j,k,ispec) &
- - TWO_THIRDS * scale_factor_minus_one * mul
- c13store_crust_mantle(i,j,k,ispec) = c13store_crust_mantle(i,j,k,ispec) &
- - TWO_THIRDS * scale_factor_minus_one * mul
- c22store_crust_mantle(i,j,k,ispec) = c22store_crust_mantle(i,j,k,ispec) &
- + FOUR_THIRDS * scale_factor_minus_one * mul
- c23store_crust_mantle(i,j,k,ispec) = c23store_crust_mantle(i,j,k,ispec) &
- - TWO_THIRDS * scale_factor_minus_one * mul
- c33store_crust_mantle(i,j,k,ispec) = c33store_crust_mantle(i,j,k,ispec) &
- + FOUR_THIRDS * scale_factor_minus_one * mul
- c44store_crust_mantle(i,j,k,ispec) = c44store_crust_mantle(i,j,k,ispec) &
- + scale_factor_minus_one * mul
- c55store_crust_mantle(i,j,k,ispec) = c55store_crust_mantle(i,j,k,ispec) &
- + scale_factor_minus_one * mul
- c66store_crust_mantle(i,j,k,ispec) = c66store_crust_mantle(i,j,k,ispec) &
- + scale_factor_minus_one * mul
- else
- if(MOVIE_VOLUME .and. SIMULATION_TYPE==3) then
- ! store the original value of \mu to comput \mu*\eps
- 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
+ ! scales transverse isotropic values for mu_h
+ 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
- ! scales transverse isotropic values for mu_h
- 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
enddo
- enddo
- enddo
+ enddo
enddo ! END DO CRUST MANTLE
! rescale in inner core
@@ -798,7 +804,11 @@
do k=1,NGLLZ
do j=1,NGLLY
do i=1,NGLLX
- scale_factor_minus_one = factor_scale_inner_core(i,j,k,ispec) - 1.0
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ scale_factor_minus_one = factor_scale_inner_core(i,j,k,ispec) - 1.0
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ scale_factor_minus_one = factor_scale_inner_core(1,1,1,ispec) - 1.0
+ endif
if(ANISOTROPIC_INNER_CORE_VAL) then
mul = muvstore_inner_core(i,j,k,ispec)
@@ -814,7 +824,11 @@
+ scale_factor_minus_one * mul
endif
- muvstore_inner_core(i,j,k,ispec) = muvstore_inner_core(i,j,k,ispec) * factor_scale_inner_core(i,j,k,ispec)
+ if (ATTENUATION_VAL .and. ATTENUATION_3D_VAL) then
+ muvstore_inner_core(i,j,k,ispec) = muvstore_inner_core(i,j,k,ispec) * factor_scale_inner_core(i,j,k,ispec)
+ else if (ATTENUATION_VAL .and. .not. ATTENUATION_3D_VAL) then
+ muvstore_inner_core(i,j,k,ispec) = muvstore_inner_core(i,j,k,ispec) * factor_scale_inner_core(1,1,1,ispec)
+ endif
enddo
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_forward_arrays.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_forward_arrays.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/read_forward_arrays.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -245,29 +245,28 @@
endif
if (ATTENUATION_VAL) then
- read(55) b_R_xx_crust_mantle
- read(55) b_R_yy_crust_mantle
- read(55) b_R_xy_crust_mantle
- read(55) b_R_xz_crust_mantle
- read(55) b_R_yz_crust_mantle
-
- read(55) b_R_xx_inner_core
- read(55) b_R_yy_inner_core
- read(55) b_R_xy_inner_core
- read(55) b_R_xz_inner_core
- read(55) b_R_yz_inner_core
-
- ! note: for kernel simulations (SIMULATION_TYPE == 3), attenuation is
- ! only mimicking effects on phase shifts, but not on amplitudes.
- ! flag USE_ATTENUATION_MIMIC will have to be set to true in this case.
- !
- ! arrays b_R_xx, ... are not used when USE_ATTENUATION_MIMIC is set,
- ! therefore no need to transfer arrays onto GPU
- !if(GPU_MODE) then
- !endif
-
+ read(55) b_R_xx_crust_mantle
+ read(55) b_R_yy_crust_mantle
+ read(55) b_R_xy_crust_mantle
+ read(55) b_R_xz_crust_mantle
+ read(55) b_R_yz_crust_mantle
+
+ read(55) b_R_xx_inner_core
+ read(55) b_R_yy_inner_core
+ read(55) b_R_xy_inner_core
+ read(55) b_R_xz_inner_core
+ read(55) b_R_yz_inner_core
+
+ ! note: for kernel simulations (SIMULATION_TYPE == 3), attenuation is
+ ! only mimicking effects on phase shifts, but not on amplitudes.
+ ! flag USE_ATTENUATION_MIMIC will have to be set to true in this case.
+ !
+ ! arrays b_R_xx, ... are not used when USE_ATTENUATION_MIMIC is set,
+ ! therefore no need to transfer arrays onto GPU
+ !if(GPU_MODE) then
+ !endif
+
endif
-
close(55)
end subroutine read_forward_arrays
Modified: seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_forward_arrays.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_forward_arrays.f90 2012-07-19 22:05:23 UTC (rev 20530)
+++ seismo/3D/SPECFEM3D_GLOBE/branches/SPECFEM3D_GLOBE_SUNFLOWER/src/specfem3D/save_forward_arrays.f90 2012-07-20 15:52:51 UTC (rev 20531)
@@ -114,17 +114,17 @@
write(55) B_array_rotation
endif
if (ATTENUATION_VAL) then
- write(55) R_xx_crust_mantle
- write(55) R_yy_crust_mantle
- write(55) R_xy_crust_mantle
- write(55) R_xz_crust_mantle
- write(55) R_yz_crust_mantle
+ write(55) R_xx_crust_mantle
+ write(55) R_yy_crust_mantle
+ write(55) R_xy_crust_mantle
+ write(55) R_xz_crust_mantle
+ write(55) R_yz_crust_mantle
- write(55) R_xx_inner_core
- write(55) R_yy_inner_core
- write(55) R_xy_inner_core
- write(55) R_xz_inner_core
- write(55) R_yz_inner_core
+ write(55) R_xx_inner_core
+ write(55) R_yy_inner_core
+ write(55) R_xy_inner_core
+ write(55) R_xz_inner_core
+ write(55) R_yz_inner_core
endif
close(55)
More information about the CIG-COMMITS
mailing list