[cig-commits] r21225 - seismo/2D/SPECFEM2D/trunk/src/specfem2D
xie.zhinan at geodynamics.org
xie.zhinan at geodynamics.org
Mon Jan 14 01:14:41 PST 2013
Author: xie.zhinan
Date: 2013-01-14 01:14:40 -0800 (Mon, 14 Jan 2013)
New Revision: 21225
Modified:
seismo/2D/SPECFEM2D/trunk/src/specfem2D/attenuation_model.f90
seismo/2D/SPECFEM2D/trunk/src/specfem2D/compute_forces_viscoelastic.F90
seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_external_model.f90
seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90
Log:
solve data type inconsistency when set CUSTOM_REAL=real in viscoelastic simulation
Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/attenuation_model.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/attenuation_model.f90 2013-01-12 21:24:23 UTC (rev 21224)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/attenuation_model.f90 2013-01-14 09:14:40 UTC (rev 21225)
@@ -53,8 +53,8 @@
integer :: N_SLS
double precision :: QKappa_attenuation,Qmu_attenuation,f0_attenuation
- double precision, dimension(N_SLS) :: inv_tau_sigma_nu1,phi_nu1,inv_tau_sigma_nu2,phi_nu2
- double precision :: Mu_nu1,Mu_nu2
+ real(kind=CUSTOM_REAL), dimension(N_SLS) :: inv_tau_sigma_nu1,phi_nu1,inv_tau_sigma_nu2,phi_nu2
+ real(kind=CUSTOM_REAL) :: Mu_nu1,Mu_nu2
integer :: i_sls
@@ -134,19 +134,35 @@
!
!--- other constants computed from the parameters above, do not modify
!
- inv_tau_sigma_nu1(:) = ONE / tau_sigma_nu1(:)
- inv_tau_sigma_nu2(:) = ONE / tau_sigma_nu2(:)
+ if(CUSTOM_REAL == SIZE_REAL) then
+ inv_tau_sigma_nu1(:) = sngl(dble(ONE) / tau_sigma_nu1(:))
+ inv_tau_sigma_nu2(:) = sngl(dble(ONE) / tau_sigma_nu2(:))
- phi_nu1(:) = (ONE - tau_epsilon_nu1(:)/tau_sigma_nu1(:)) / tau_sigma_nu1(:)
- phi_nu2(:) = (ONE - tau_epsilon_nu2(:)/tau_sigma_nu2(:)) / tau_sigma_nu2(:)
+ phi_nu1(:) = sngl((dble(ONE) - tau_epsilon_nu1(:)/tau_sigma_nu1(:)) / tau_sigma_nu1(:))
+ phi_nu2(:) = sngl((dble(ONE) - tau_epsilon_nu2(:)/tau_sigma_nu2(:)) / tau_sigma_nu2(:))
- Mu_nu1 = ONE
- Mu_nu2 = ONE
+ Mu_nu1 = dble(ONE)
+ Mu_nu2 = dble(ONE)
- do i_sls = 1,N_SLS
- Mu_nu1 = Mu_nu1 - (ONE - tau_epsilon_nu1(i_sls)/tau_sigma_nu1(i_sls))
- Mu_nu2 = Mu_nu2 - (ONE - tau_epsilon_nu2(i_sls)/tau_sigma_nu2(i_sls))
- enddo
+ do i_sls = 1,N_SLS
+ Mu_nu1 = sngl(dble(Mu_nu1) - (dble(ONE) - tau_epsilon_nu1(i_sls)/tau_sigma_nu1(i_sls)))
+ Mu_nu2 = sngl(dble(Mu_nu2) - (dble(ONE) - tau_epsilon_nu2(i_sls)/tau_sigma_nu2(i_sls)))
+ enddo
+ else
+ inv_tau_sigma_nu1(:) = ONE / tau_sigma_nu1(:)
+ inv_tau_sigma_nu2(:) = ONE / tau_sigma_nu2(:)
+ phi_nu1(:) = (ONE - tau_epsilon_nu1(:)/tau_sigma_nu1(:)) / tau_sigma_nu1(:)
+ phi_nu2(:) = (ONE - tau_epsilon_nu2(:)/tau_sigma_nu2(:)) / tau_sigma_nu2(:)
+
+ Mu_nu1 = ONE
+ Mu_nu2 = ONE
+
+ do i_sls = 1,N_SLS
+ Mu_nu1 = Mu_nu1 - (ONE - tau_epsilon_nu1(i_sls)/tau_sigma_nu1(i_sls))
+ Mu_nu2 = Mu_nu2 - (ONE - tau_epsilon_nu2(i_sls)/tau_sigma_nu2(i_sls))
+ enddo
+ endif
+
end subroutine attenuation_model
Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/compute_forces_viscoelastic.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/compute_forces_viscoelastic.F90 2013-01-12 21:24:23 UTC (rev 21224)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/compute_forces_viscoelastic.F90 2013-01-14 09:14:40 UTC (rev 21225)
@@ -109,7 +109,7 @@
logical, dimension(4,nelemabs) :: codeabs
real(kind=CUSTOM_REAL), dimension(3,nglob) :: accel_elastic,veloc_elastic,displ_elastic
- real(kind=CUSTOM_REAL), dimension(3,nglob) :: temp_displ_elastic
+ real(kind=CUSTOM_REAL), dimension(3,nglob) :: displ_elastic_np1
double precision, dimension(2,numat) :: density
double precision, dimension(4,3,numat) :: poroelastcoef
double precision, dimension(6,numat) :: anisotropy
@@ -133,8 +133,8 @@
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec,N_SLS) :: e1_LDDRK,e11_LDDRK,e13_LDDRK
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec,N_SLS) :: e1_initial_rk,e11_initial_rk,e13_initial_rk
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec,N_SLS,stage_time_scheme) :: e1_force_RK, e11_force_RK, e13_force_RK
- double precision, dimension(NGLLX,NGLLZ,nspec,N_SLS) :: inv_tau_sigma_nu1,phi_nu1,inv_tau_sigma_nu2,phi_nu2
- double precision, dimension(NGLLX,NGLLZ,nspec) :: Mu_nu1,Mu_nu2
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec,N_SLS) :: inv_tau_sigma_nu1,phi_nu1,inv_tau_sigma_nu2,phi_nu2
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec) :: Mu_nu1,Mu_nu2
real(kind=CUSTOM_REAL) :: e1_sum,e11_sum,e13_sum
integer :: i_sls
@@ -271,14 +271,6 @@
displ_elastic_new = displ_elastic + deltat * veloc_elastic
endif
-
- ! compute Grad(displ_elastic) at time step n for attenuation
- if(ATTENUATION_VISCOELASTIC_SOLID) then
- temp_displ_elastic = displ_elastic + deltat * veloc_elastic + 0.5 * deltat**2 * accel_elastic
- call compute_gradient_attenuation(displ_elastic,dux_dxl_n,duz_dxl_n, &
- dux_dzl_n,duz_dzl_n,xix,xiz,gammax,gammaz,ibool,elastic,hprime_xx,hprime_zz,nspec,nglob)
- endif
-
ifirstelem = 1
ilastelem = nspec
@@ -2070,10 +2062,17 @@
! implement attenuation
if(ATTENUATION_VISCOELASTIC_SOLID) then
+ ! compute Grad(displ_elastic) at time step n for attenuation
+ call compute_gradient_attenuation(displ_elastic,dux_dxl_n,duz_dxl_n, &
+ dux_dzl_n,duz_dzl_n,xix,xiz,gammax,gammaz,ibool,elastic,hprime_xx,hprime_zz,nspec,nglob)
+
+ displ_elastic_np1 = displ_elastic + deltat * veloc_elastic
+
! compute Grad(displ_elastic) at time step n+1 for attenuation
- call compute_gradient_attenuation(temp_displ_elastic,dux_dxl_np1,duz_dxl_np1, &
+ call compute_gradient_attenuation(displ_elastic_np1,dux_dxl_np1,duz_dxl_np1, &
dux_dzl_np1,duz_dzl_np1,xix,xiz,gammax,gammaz,ibool,elastic,hprime_xx,hprime_zz,nspec,nglob)
+
! update memory variables with fourth-order Runge-Kutta time scheme for attenuation
! loop over spectral elements
do ispec = 1,nspec
Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_external_model.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_external_model.f90 2013-01-12 21:24:23 UTC (rev 21224)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_external_model.f90 2013-01-14 09:14:40 UTC (rev 21225)
@@ -71,12 +71,12 @@
! for attenuation
integer :: N_SLS
- double precision :: Mu_nu1_sent,Mu_nu2_sent
- double precision, dimension(N_SLS) :: inv_tau_sigma_nu1_sent,phi_nu1_sent, &
+ real(kind=CUSTOM_REAL) :: Mu_nu1_sent,Mu_nu2_sent
+ real(kind=CUSTOM_REAL), dimension(N_SLS) :: inv_tau_sigma_nu1_sent,phi_nu1_sent, &
inv_tau_sigma_nu2_sent,phi_nu2_sent
- double precision, dimension(NGLLX,NGLLZ,nspec,N_SLS) :: inv_tau_sigma_nu1,phi_nu1, &
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec,N_SLS) :: inv_tau_sigma_nu1,phi_nu1, &
inv_tau_sigma_nu2,phi_nu2
- double precision, dimension(NGLLX,NGLLZ,nspec) :: Mu_nu1,Mu_nu2
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,nspec) :: Mu_nu1,Mu_nu2
double precision, dimension(NGLLX,NGLLZ,nspec) :: QKappa_attenuationext,Qmu_attenuationext
! for anisotropy
Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90 2013-01-12 21:24:23 UTC (rev 21224)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90 2013-01-14 09:14:40 UTC (rev 21225)
@@ -592,10 +592,10 @@
real(kind=CUSTOM_REAL), dimension(:,:,:,:), allocatable :: e1_LDDRK,e11_LDDRK,e13_LDDRK
real(kind=CUSTOM_REAL), dimension(:,:,:,:), allocatable :: e1_initial_rk,e11_initial_rk,e13_initial_rk
real(kind=CUSTOM_REAL), dimension(:,:,:,:,:), allocatable :: e1_force_rk,e11_force_rk,e13_force_rk
- double precision, dimension(:,:,:,:), allocatable :: inv_tau_sigma_nu1,phi_nu1,inv_tau_sigma_nu2,phi_nu2
- double precision, dimension(:), allocatable :: inv_tau_sigma_nu1_sent,phi_nu1_sent,inv_tau_sigma_nu2_sent,phi_nu2_sent
- double precision, dimension(:,:,:) , allocatable :: Mu_nu1,Mu_nu2
- double precision :: Mu_nu1_sent,Mu_nu2_sent
+ real(kind=CUSTOM_REAL), dimension(:,:,:,:), allocatable :: inv_tau_sigma_nu1,phi_nu1,inv_tau_sigma_nu2,phi_nu2
+ real(kind=CUSTOM_REAL), dimension(:), allocatable :: inv_tau_sigma_nu1_sent,phi_nu1_sent,inv_tau_sigma_nu2_sent,phi_nu2_sent
+ real(kind=CUSTOM_REAL), dimension(:,:,:) , allocatable :: Mu_nu1,Mu_nu2
+ real(kind=CUSTOM_REAL) :: Mu_nu1_sent,Mu_nu2_sent
real(kind=CUSTOM_REAL), dimension(:,:,:), allocatable :: &
dux_dxl_n,duz_dzl_n,duz_dxl_n,dux_dzl_n,dux_dxl_np1,duz_dzl_np1,duz_dxl_np1,dux_dzl_np1
More information about the CIG-COMMITS
mailing list