[cig-commits] r21974 - seismo/3D/SPECFEM3D/trunk/src/specfem3D
xie.zhinan at geodynamics.org
xie.zhinan at geodynamics.org
Fri May 3 01:50:26 PDT 2013
Author: xie.zhinan
Date: 2013-05-03 01:50:25 -0700 (Fri, 03 May 2013)
New Revision: 21974
Modified:
seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_Dev.F90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_noDev.f90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_accel_contribution.f90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_memory_variables.f90
Log:
clean the pml code a little
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_Dev.F90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_Dev.F90 2013-05-03 07:46:41 UTC (rev 21973)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_Dev.F90 2013-05-03 08:50:25 UTC (rev 21974)
@@ -35,11 +35,11 @@
kappastore,mustore,jacobian,ibool, &
ATTENUATION,deltat, &
one_minus_sum_beta,factor_common,&
- one_minus_sum_beta_kappa,factor_common_kappa,& !ZN
- alphaval,betaval,gammaval,& !ZN
+ one_minus_sum_beta_kappa,factor_common_kappa,&
+ alphaval,betaval,gammaval,&
NSPEC_ATTENUATION_AB,NSPEC_ATTENUATION_AB_Kappa, &
- R_trace,R_xx,R_yy,R_xy,R_xz,R_yz, & !ZN
- epsilondev_trace,epsilondev_xx,epsilondev_yy,epsilondev_xy, & !ZN
+ R_trace,R_xx,R_yy,R_xy,R_xz,R_yz, &
+ epsilondev_trace,epsilondev_xx,epsilondev_yy,epsilondev_xy, &
epsilondev_xz,epsilondev_yz,epsilon_trace_over_3, &
ANISOTROPY,NSPEC_ANISO, &
c11store,c12store,c13store,c14store,c15store,c16store,&
@@ -94,17 +94,17 @@
integer :: NSPEC_ATTENUATION_AB,NSPEC_ATTENUATION_AB_Kappa
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB) :: one_minus_sum_beta
real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB) :: factor_common
- real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa) :: one_minus_sum_beta_kappa !ZN
- real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa) :: factor_common_kappa !ZN
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa) :: one_minus_sum_beta_kappa
+ real(kind=CUSTOM_REAL), dimension(N_SLS,NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa) :: factor_common_kappa
real(kind=CUSTOM_REAL), dimension(N_SLS) :: alphaval,betaval,gammaval
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB,N_SLS) :: &
R_xx,R_yy,R_xy,R_xz,R_yz
- real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa,N_SLS) :: R_trace !ZN
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa,N_SLS) :: R_trace
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_STRAIN_ONLY) :: &
epsilondev_xx,epsilondev_yy,epsilondev_xy,epsilondev_xz,epsilondev_yz
- real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa) :: epsilondev_trace !ZN
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ATTENUATION_AB_Kappa) :: epsilondev_trace
real(kind=CUSTOM_REAL),dimension(NGLLX,NGLLY,NGLLZ,NSPEC_ADJOINT) :: epsilon_trace_over_3
! anisotropy
@@ -198,10 +198,10 @@
equivalence(tempz3_att,C3_mxm_m2_m1_5points_att)
! local attenuation parameters
- real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: epsilondev_trace_loc,epsilondev_xx_loc, & !ZN
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: epsilondev_trace_loc,epsilondev_xx_loc, &
epsilondev_yy_loc, epsilondev_xy_loc, epsilondev_xz_loc, epsilondev_yz_loc
- real(kind=CUSTOM_REAL) R_xx_val1,R_yy_val1,R_xx_val2,R_yy_val2,R_xx_val3,R_yy_val3, & !ZN
- R_trace_val1,R_trace_val2,R_trace_val3 !ZN
+ real(kind=CUSTOM_REAL) R_xx_val1,R_yy_val1,R_xx_val2,R_yy_val2,R_xx_val3,R_yy_val3, &
+ R_trace_val1,R_trace_val2,R_trace_val3
real(kind=CUSTOM_REAL) factor_loc,alphaval_loc,betaval_loc,gammaval_loc
real(kind=CUSTOM_REAL) Sn,Snp1
real(kind=CUSTOM_REAL) templ
@@ -563,7 +563,7 @@
if(ATTENUATION) then
! use unrelaxed parameters if attenuation
mul = mul * one_minus_sum_beta(i,j,k,ispec)
- if(FULL_ATTENUATION_SOLID) kappal = kappal * one_minus_sum_beta_kappa(i,j,k,ispec) !ZN
+ if(FULL_ATTENUATION_SOLID) kappal = kappal * one_minus_sum_beta_kappa(i,j,k,ispec)
endif
! full anisotropic case, stress calculations
@@ -823,13 +823,13 @@
gammaval_loc = gammaval(i_sls)
if(FULL_ATTENUATION_SOLID) then
- ! term in trace !ZN
- factor_loc = kappastore(i,j,k,ispec) * factor_common_kappa(i_sls,i,j,k,ispec) !ZN
+ ! term in trace
+ factor_loc = kappastore(i,j,k,ispec) * factor_common_kappa(i_sls,i,j,k,ispec)
- Sn = factor_loc * epsilondev_trace(i,j,k,ispec) !ZN
- Snp1 = factor_loc * epsilondev_trace_loc(i,j,k) !ZN
- R_trace(i,j,k,ispec,i_sls) = alphaval_loc * R_trace(i,j,k,ispec,i_sls) + & !ZN
- betaval_loc * Sn + gammaval_loc * Snp1 !ZN
+ Sn = factor_loc * epsilondev_trace(i,j,k,ispec)
+ Snp1 = factor_loc * epsilondev_trace_loc(i,j,k)
+ R_trace(i,j,k,ispec,i_sls) = alphaval_loc * R_trace(i,j,k,ispec,i_sls) + &
+ betaval_loc * Sn + gammaval_loc * Snp1
endif
! term in xx yy zz xy xz yz
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_noDev.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_noDev.f90 2013-05-03 07:46:41 UTC (rev 21973)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_viscoelastic_noDev.f90 2013-05-03 08:50:25 UTC (rev 21974)
@@ -882,7 +882,7 @@
do iface=1,num_abs_boundary_faces
ispec = abs_boundary_ispec(iface)
if (ispec_is_inner(ispec) .eqv. phase_is_inner) then
- if( ispec_is_elastic(ispec) ) then
+ if( ispec_is_elastic(ispec) .and. is_CPML(ispec)) then
! reference gll points on boundary face
do igll = 1,NGLLSQUARE
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_accel_contribution.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_accel_contribution.f90 2013-05-03 07:46:41 UTC (rev 21973)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_accel_contribution.f90 2013-05-03 08:50:25 UTC (rev 21974)
@@ -35,13 +35,12 @@
! Anisotropic-Medium PML for Vector FETD With Modified Basis Functions,
! IEEE Transactions on Antennas and Propagation, vol. 54, no. 1, (2006)
- use specfem_par, only: ibool,wgllwgll_yz,wgllwgll_xz,wgllwgll_xy,jacobian,it,deltat,kappastore,rhostore
+ use specfem_par, only: it,deltat,wgll_cube,jacobian,ibool,rhostore
use specfem_par_elastic, only: displ,veloc
- use pml_par, only: NSPEC_CPML,CPML_regions,spec_to_CPML, &
- d_store_x,d_store_y,d_store_z,K_store_x,K_store_y,K_store_z,alpha_store,&
- rmemory_displ_elastic,accel_elastic_CPML
- use constants, only: CUSTOM_REAL,NDIM,NGLLX,NGLLY,NGLLZ, &
- CPML_X_ONLY,CPML_Y_ONLY,CPML_Z_ONLY,CPML_XY_ONLY,CPML_XZ_ONLY,CPML_YZ_ONLY,CPML_XYZ
+ use pml_par, only: CPML_regions,d_store_x,d_store_y,d_store_z,K_store_x,K_store_y,K_store_z,&
+ alpha_store,rmemory_displ_elastic,accel_elastic_CPML
+ use constants, only: CUSTOM_REAL,NGLLX,NGLLY,NGLLZ,CPML_X_ONLY,CPML_Y_ONLY,CPML_Z_ONLY, &
+ CPML_XY_ONLY,CPML_XZ_ONLY,CPML_YZ_ONLY,CPML_XYZ
implicit none
@@ -49,7 +48,7 @@
! local parameters
integer :: i,j,k,iglob
- real(kind=CUSTOM_REAL) :: fac1,fac2,fac3,fac4,rhol,jacobianl
+ real(kind=CUSTOM_REAL) :: wgllcube,rhol,jacobianl
real(kind=CUSTOM_REAL) :: bb,coef0_1,coef1_1,coef2_1,coef0_2,coef1_2,coef2_2,coef0_3,coef1_3,coef2_3
real(kind=CUSTOM_REAL) :: A0,A1,A2,A3,A4,A5,temp_A3! for convolution of acceleration
@@ -59,10 +58,7 @@
rhol = rhostore(i,j,k,ispec)
jacobianl = jacobian(i,j,k,ispec)
iglob = ibool(i,j,k,ispec)
- fac1 = wgllwgll_yz(j,k)
- fac2 = wgllwgll_xz(i,k)
- fac3 = wgllwgll_xy(i,j)
- fac4 = sqrt(fac1 * fac2 * fac3)
+ wgllcube = wgll_cube(i,j,k)
if( CPML_regions(ispec_CPML) == CPML_X_ONLY ) then
!------------------------------------------------------------------------------
@@ -438,21 +434,21 @@
endif
- accel_elastic_CPML(1,i,j,k,ispec_CPML) = fac4 * rhol * jacobianl * &
+ accel_elastic_CPML(1,i,j,k,ispec_CPML) = wgllcube * rhol * jacobianl * &
( A1 * veloc(1,iglob) + A2 * displ(1,iglob) + &
A3 * rmemory_displ_elastic(1,i,j,k,ispec_CPML,1) + &
A4 * rmemory_displ_elastic(1,i,j,k,ispec_CPML,2) + &
A5 * rmemory_displ_elastic(1,i,j,k,ispec_CPML,3) &
)
- accel_elastic_CPML(2,i,j,k,ispec_CPML) = fac4 * rhol * jacobianl * &
+ accel_elastic_CPML(2,i,j,k,ispec_CPML) = wgllcube * rhol * jacobianl * &
( A1 * veloc(2,iglob) + A2 * displ(2,iglob) + &
A3 * rmemory_displ_elastic(2,i,j,k,ispec_CPML,1) + &
A4 * rmemory_displ_elastic(2,i,j,k,ispec_CPML,2) + &
A5 * rmemory_displ_elastic(2,i,j,k,ispec_CPML,3) &
)
- accel_elastic_CPML(3,i,j,k,ispec_CPML) = fac4 * rhol * jacobianl * &
+ accel_elastic_CPML(3,i,j,k,ispec_CPML) = wgllcube * rhol * jacobianl * &
( A1 * veloc(3,iglob) + A2 * displ(3,iglob) + &
A3 * rmemory_displ_elastic(3,i,j,k,ispec_CPML,1) + &
A4 * rmemory_displ_elastic(3,i,j,k,ispec_CPML,2) + &
@@ -478,13 +474,12 @@
! Anisotropic-Medium PML for Vector FETD With Modified Basis Functions,
! IEEE Transactions on Antennas and Propagation, vol. 54, no. 1, (2006)
- use specfem_par, only: ibool,wgllwgll_yz,wgllwgll_xz,wgllwgll_xy,jacobian,it,deltat,kappastore,rhostore
- use specfem_par_acoustic, only: potential_acoustic,potential_dot_acoustic,potential_dot_dot_acoustic
- use pml_par, only: NSPEC_CPML,CPML_regions,spec_to_CPML, &
- d_store_x,d_store_y,d_store_z,K_store_x,K_store_y,K_store_z,alpha_store,&
- rmemory_potential_acoustic, potential_dot_dot_acoustic_CPML
- use constants, only: CUSTOM_REAL,NDIM,NGLLX,NGLLY,NGLLZ, &
- CPML_X_ONLY,CPML_Y_ONLY,CPML_Z_ONLY,CPML_XY_ONLY,CPML_XZ_ONLY,CPML_YZ_ONLY,CPML_XYZ
+ use specfem_par, only: it,deltat,wgll_cube,jacobian,ibool,kappastore
+ use specfem_par_acoustic, only: potential_acoustic,potential_dot_acoustic
+ use pml_par, only: CPML_regions,d_store_x,d_store_y,d_store_z,K_store_x,K_store_y,K_store_z,&
+ alpha_store,rmemory_potential_acoustic, potential_dot_dot_acoustic_CPML
+ use constants, only: CUSTOM_REAL,NGLLX,NGLLY,NGLLZ,CPML_X_ONLY,CPML_Y_ONLY,CPML_Z_ONLY, &
+ CPML_XY_ONLY,CPML_XZ_ONLY,CPML_YZ_ONLY,CPML_XYZ
implicit none
@@ -492,20 +487,17 @@
! local parameters
integer :: i,j,k,iglob
- real(kind=CUSTOM_REAL) :: fac1,fac2,fac3,fac4,kappal,jacobianl
+ real(kind=CUSTOM_REAL) :: wgllcube,kappal_inv,jacobianl
real(kind=CUSTOM_REAL) :: bb,coef0_1,coef1_1,coef2_1,coef0_2,coef1_2,coef2_2,coef0_3,coef1_3,coef2_3
real(kind=CUSTOM_REAL) :: A0,A1,A2,A3,A4,A5,temp_A3 ! for convolution of acceleration
do k=1,NGLLZ
do j=1,NGLLY
do i=1,NGLLX
- kappal = kappastore(i,j,k,ispec)
+ kappal_inv = 1.d0 / kappastore(i,j,k,ispec)
jacobianl = jacobian(i,j,k,ispec)
iglob = ibool(i,j,k,ispec)
- fac1 = wgllwgll_yz(j,k)
- fac2 = wgllwgll_xz(i,k)
- fac3 = wgllwgll_xy(i,j)
- fac4 = sqrt(fac1 * fac2 * fac3)
+ wgllcube = wgll_cube(i,j,k)
if( CPML_regions(ispec_CPML) == CPML_X_ONLY ) then
!------------------------------------------------------------------------------
@@ -799,7 +791,7 @@
endif
- potential_dot_dot_acoustic_CPML(i,j,k,ispec_CPML) = fac4 * 1.0/kappal *jacobianl * &
+ potential_dot_dot_acoustic_CPML(i,j,k,ispec_CPML) = wgllcube * kappal_inv *jacobianl * &
( A1 * potential_dot_acoustic(iglob) + A2 * potential_acoustic(iglob) + &
A3 * rmemory_potential_acoustic(i,j,k,ispec_CPML,1)+ &
A4 * rmemory_potential_acoustic(i,j,k,ispec_CPML,2)+ &
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_memory_variables.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_memory_variables.f90 2013-05-03 07:46:41 UTC (rev 21973)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/pml_compute_memory_variables.f90 2013-05-03 08:50:25 UTC (rev 21974)
@@ -37,7 +37,7 @@
use specfem_par, only: wgllwgll_xy,wgllwgll_xz,wgllwgll_yz,it,deltat, &
xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz,jacobian, &
- kappastore,mustore,rhostore
+ kappastore,mustore
use pml_par
use constants, only: NGLLX,NGLLY,NGLLZ,FOUR_THIRDS, &
CPML_X_ONLY,CPML_Y_ONLY,CPML_Z_ONLY,CPML_XY_ONLY,CPML_XZ_ONLY,CPML_YZ_ONLY,CPML_XYZ
More information about the CIG-COMMITS
mailing list