[cig-commits] [commit] devel: updates coupling routines to avoid implicit (hidden) memory copies, fixes #122; updates rule dependencies (ba1bf2c)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Fri Nov 28 14:14:23 PST 2014
Repository : https://github.com/geodynamics/specfem3d
On branch : devel
Link : https://github.com/geodynamics/specfem3d/compare/9a84d06c76e869f5276019e4f84affce23830a4d...8dce71b713c1fd0b510e7538fea0f2307c7b29e8
>---------------------------------------------------------------
commit ba1bf2c8590b29430afb6cc341cfdcf67b5a0232
Author: daniel peter <peterda at ethz.ch>
Date: Fri Nov 28 22:50:48 2014 +0100
updates coupling routines to avoid implicit (hidden) memory copies, fixes #122; updates rule dependencies
>---------------------------------------------------------------
ba1bf2c8590b29430afb6cc341cfdcf67b5a0232
src/auxiliaries/rules.mk | 4 +-
src/check_mesh_quality_CUBIT_Abaqus/rules.mk | 8 +-
src/specfem3D/compute_coupling_acoustic_el.f90 | 92 +++++++++++-----------
src/specfem3D/compute_coupling_viscoelastic_ac.f90 | 91 +++++++++++----------
.../compute_forces_acoustic_calling_routine.f90 | 6 +-
src/specfem3D/compute_forces_viscoelastic_Dev.F90 | 9 +--
...compute_forces_viscoelastic_calling_routine.F90 | 12 +--
src/specfem3D/fault_solver_dynamic.f90 | 71 +++++++++--------
src/specfem3D/pml_compute_memory_variables.f90 | 8 +-
9 files changed, 149 insertions(+), 152 deletions(-)
diff --git a/src/auxiliaries/rules.mk b/src/auxiliaries/rules.mk
index b4fd130..b0b7845 100644
--- a/src/auxiliaries/rules.mk
+++ b/src/auxiliaries/rules.mk
@@ -158,8 +158,8 @@ create_movie_shakemap_AVS_DX_GMT: xcreate_movie_shakemap_AVS_DX_GMT
xcreate_movie_shakemap_AVS_DX_GMT: $E/xcreate_movie_shakemap_AVS_DX_GMT
-$E/xconvolve_source_timefunction: $O/convolve_source_timefunction.aux.o
- ${FCCOMPILE_CHECK} -o ${E}/xconvolve_source_timefunction $O/convolve_source_timefunction.aux.o
+$E/xconvolve_source_timefunction: $O/convolve_source_timefunction.aux.o $O/shared_par.shared_module.o
+ ${FCCOMPILE_CHECK} -o ${E}/xconvolve_source_timefunction $O/convolve_source_timefunction.aux.o $O/shared_par.shared_module.o
$E/xcombine_surf_data: $(combine_surf_data_auxiliaries_OBJECTS) $(combine_surf_data_auxiliaries_SHARED_OBJECTS)
${FCLINK} -o $@ $+
diff --git a/src/check_mesh_quality_CUBIT_Abaqus/rules.mk b/src/check_mesh_quality_CUBIT_Abaqus/rules.mk
index 87ffa68..9512170 100644
--- a/src/check_mesh_quality_CUBIT_Abaqus/rules.mk
+++ b/src/check_mesh_quality_CUBIT_Abaqus/rules.mk
@@ -73,11 +73,11 @@ xmultiply_CUBIT_Abaqus_mesh_by_1000: $E/xmultiply_CUBIT_Abaqus_mesh_by_1000
# rules for the pure Fortran version
-$E/xcheck_mesh_quality_CUBIT_Abaqus: $O/check_mesh_quality_CUBIT_Abaqus.check.o
- ${FCLINK} -o $E/xcheck_mesh_quality_CUBIT_Abaqus $O/check_mesh_quality_CUBIT_Abaqus.check.o
+$E/xcheck_mesh_quality_CUBIT_Abaqus: $O/check_mesh_quality_CUBIT_Abaqus.check.o $O/shared_par.shared_module.o
+ ${FCLINK} -o $E/xcheck_mesh_quality_CUBIT_Abaqus $O/check_mesh_quality_CUBIT_Abaqus.check.o $O/shared_par.shared_module.o
-$E/xconvert_skewness_to_angle: $O/convert_skewness_to_angle.check.o
- ${FCLINK} -o $E/xconvert_skewness_to_angle $O/convert_skewness_to_angle.check.o
+$E/xconvert_skewness_to_angle: $O/convert_skewness_to_angle.check.o $O/shared_par.shared_module.o
+ ${FCLINK} -o $E/xconvert_skewness_to_angle $O/convert_skewness_to_angle.check.o $O/shared_par.shared_module.o
$E/xmultiply_CUBIT_Abaqus_mesh_by_1000: $O/multiply_CUBIT_Abaqus_mesh_by_1000.check.o
${FCLINK} -o $E/xmultiply_CUBIT_Abaqus_mesh_by_1000 $O/multiply_CUBIT_Abaqus_mesh_by_1000.check.o
diff --git a/src/specfem3D/compute_coupling_acoustic_el.f90 b/src/specfem3D/compute_coupling_acoustic_el.f90
index 76bc50d..4711b90 100644
--- a/src/specfem3D/compute_coupling_acoustic_el.f90
+++ b/src/specfem3D/compute_coupling_acoustic_el.f90
@@ -40,38 +40,40 @@
! returns the updated pressure array: potential_dot_dot_acoustic
- use constants
+ use constants,only: CUSTOM_REAL,NDIM,NGLLX,NGLLY,NGLLZ,NGLLSQUARE
+
use pml_par, only: NSPEC_CPML
implicit none
- integer :: NSPEC_AB,NGLOB_AB
+ integer,intent(in) :: NSPEC_AB,NGLOB_AB
! displacement and pressure
- real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_AB) :: displ
- real(kind=CUSTOM_REAL), dimension(NGLOB_AB) :: potential_dot_dot_acoustic
- integer :: SIMULATION_TYPE
- logical :: backward_simulation
+ real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_AB),intent(in) :: displ
+ real(kind=CUSTOM_REAL), dimension(NGLOB_AB),intent(inout) :: potential_dot_dot_acoustic
+ integer,intent(in) :: SIMULATION_TYPE
+ logical,intent(in) :: backward_simulation
! global indexing
- integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_AB) :: ibool
+ integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_AB),intent(in) :: ibool
! acoustic-elastic coupling surface
- integer :: num_coupling_ac_el_faces
- real(kind=CUSTOM_REAL) :: coupling_ac_el_normal(NDIM,NGLLSQUARE,num_coupling_ac_el_faces)
- real(kind=CUSTOM_REAL) :: coupling_ac_el_jacobian2Dw(NGLLSQUARE,num_coupling_ac_el_faces)
- integer :: coupling_ac_el_ijk(3,NGLLSQUARE,num_coupling_ac_el_faces)
- integer :: coupling_ac_el_ispec(num_coupling_ac_el_faces)
- real(kind=CUSTOM_REAL), dimension(NDIM,NGLLX,NGLLY,NGLLZ,num_coupling_ac_el_faces,2) :: rmemory_coupling_ac_el_displ
+ integer,intent(in) :: num_coupling_ac_el_faces
+ real(kind=CUSTOM_REAL),intent(in) :: coupling_ac_el_normal(NDIM,NGLLSQUARE,num_coupling_ac_el_faces)
+ real(kind=CUSTOM_REAL),intent(in) :: coupling_ac_el_jacobian2Dw(NGLLSQUARE,num_coupling_ac_el_faces)
+ integer,intent(in) :: coupling_ac_el_ijk(3,NGLLSQUARE,num_coupling_ac_el_faces)
+ integer,intent(in) :: coupling_ac_el_ispec(num_coupling_ac_el_faces)
! communication overlap
- logical, dimension(NSPEC_AB) :: ispec_is_inner
- logical :: phase_is_inner
+ logical, dimension(NSPEC_AB),intent(in) :: ispec_is_inner
+ logical,intent(in) :: phase_is_inner
! CPML
- logical :: PML_CONDITIONS
- integer :: spec_to_CPML(NSPEC_AB)
- logical :: is_CPML(NSPEC_AB)
+ logical,intent(in) :: PML_CONDITIONS
+ integer,intent(in) :: spec_to_CPML(NSPEC_AB)
+ logical,intent(in) :: is_CPML(NSPEC_AB)
+ real(kind=CUSTOM_REAL), dimension(NDIM,NGLLX,NGLLY,NGLLZ,num_coupling_ac_el_faces,2),intent(inout) :: &
+ rmemory_coupling_ac_el_displ
! local parameters
real(kind=CUSTOM_REAL) :: displ_x,displ_y,displ_z,displ_n
@@ -99,38 +101,34 @@
iglob = ibool(i,j,k,ispec)
! elastic displacement on global point
- if (PML_CONDITIONS .and. NSPEC_CPML > 0) then
- if (.not. backward_simulation) then
- if (is_CPML(ispec)) then
- if (SIMULATION_TYPE == 1) then
- ispec_CPML = spec_to_CPML(ispec)
- call pml_compute_memory_variables_acoustic_elastic(ispec_CPML,iface,iglob,i,j,k,&
- displ_x,displ_y,displ_z,displ,&
- num_coupling_ac_el_faces,rmemory_coupling_ac_el_displ)
- endif
-
- if (SIMULATION_TYPE == 3) then
- ! left blank for change
- endif
-
- else
- displ_x = displ(1,iglob)
- displ_y = displ(2,iglob)
- displ_z = displ(3,iglob)
+ displ_x = displ(1,iglob)
+ displ_y = displ(2,iglob)
+ displ_z = displ(3,iglob)
+
+ ! adjoint wavefield case
+ if (SIMULATION_TYPE /= 1 .and. (.not. backward_simulation)) then
+ ! handles adjoint runs coupling between adjoint potential and adjoint elastic wavefield
+ ! adjoint definition: \partial_t^2 \bfs^\dagger=-\frac{1}{\rho}\bfnabla\phi^\dagger
+ displ_x = - displ_x
+ displ_y = - displ_y
+ displ_z = - displ_z
+ endif
+
+ ! CPML overwrite cases
+ if (PML_CONDITIONS .and. (.not. backward_simulation) .and. NSPEC_CPML > 0) then
+ if (is_CPML(ispec)) then
+ if (SIMULATION_TYPE == 1) then
+ ispec_CPML = spec_to_CPML(ispec)
+ call pml_compute_memory_variables_acoustic_elastic(ispec_CPML,iface,iglob,i,j,k,&
+ displ_x,displ_y,displ_z,displ,&
+ num_coupling_ac_el_faces,rmemory_coupling_ac_el_displ)
endif
- else
- if (is_CPML(ispec)) then
- ! left blank, since no operation needed
- else
- displ_x = displ(1,iglob)
- displ_y = displ(2,iglob)
- displ_z = displ(3,iglob)
+
+ if (SIMULATION_TYPE == 3) then
+ ! safety stop
+ stop 'TODO: Coupling acoustic-elastic for CPML in compute_coupling_acoustic_el() not implemented yet...'
endif
endif
- else
- displ_x = displ(1,iglob)
- displ_y = displ(2,iglob)
- displ_z = displ(3,iglob)
endif
! gets associated normal on GLL point
diff --git a/src/specfem3D/compute_coupling_viscoelastic_ac.f90 b/src/specfem3D/compute_coupling_viscoelastic_ac.f90
index 58a4f4b..45fd9df 100644
--- a/src/specfem3D/compute_coupling_viscoelastic_ac.f90
+++ b/src/specfem3D/compute_coupling_viscoelastic_ac.f90
@@ -40,31 +40,32 @@
! returns the updated acceleration array: accel
- use constants
+ use constants,only: CUSTOM_REAL,NDIM,NGLLX,NGLLY,NGLLZ,NGLLSQUARE
+
use pml_par,only : rmemory_coupling_el_ac_potential_dot_dot,is_CPML,spec_to_CPML,&
potential_acoustic_old,NSPEC_CPML
implicit none
- integer :: NSPEC_AB,NGLOB_AB,SIMULATION_TYPE
- logical :: backward_simulation
+ integer,intent(in) :: NSPEC_AB,NGLOB_AB,SIMULATION_TYPE
+ logical,intent(in) :: backward_simulation
! displacement and pressure
- real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_AB) :: accel
- real(kind=CUSTOM_REAL), dimension(NGLOB_AB) :: potential_dot_dot_acoustic,potential_dot_acoustic,potential_acoustic
+ real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_AB),intent(inout) :: accel
+ real(kind=CUSTOM_REAL), dimension(NGLOB_AB),intent(in) :: potential_dot_dot_acoustic,potential_dot_acoustic,potential_acoustic
! global indexing
- integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_AB) :: ibool
+ integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_AB),intent(in) :: ibool
! acoustic-elastic coupling surface
- integer :: num_coupling_ac_el_faces
- real(kind=CUSTOM_REAL) :: coupling_ac_el_normal(NDIM,NGLLSQUARE,num_coupling_ac_el_faces)
- real(kind=CUSTOM_REAL) :: coupling_ac_el_jacobian2Dw(NGLLSQUARE,num_coupling_ac_el_faces)
- integer :: coupling_ac_el_ijk(3,NGLLSQUARE,num_coupling_ac_el_faces)
- integer :: coupling_ac_el_ispec(num_coupling_ac_el_faces)
+ integer,intent(in) :: num_coupling_ac_el_faces
+ real(kind=CUSTOM_REAL),intent(in) :: coupling_ac_el_normal(NDIM,NGLLSQUARE,num_coupling_ac_el_faces)
+ real(kind=CUSTOM_REAL),intent(in) :: coupling_ac_el_jacobian2Dw(NGLLSQUARE,num_coupling_ac_el_faces)
+ integer,intent(in) :: coupling_ac_el_ijk(3,NGLLSQUARE,num_coupling_ac_el_faces)
+ integer,intent(in) :: coupling_ac_el_ispec(num_coupling_ac_el_faces)
! communication overlap
- logical, dimension(NSPEC_AB) :: ispec_is_inner
- logical :: phase_is_inner
+ logical, dimension(NSPEC_AB),intent(in) :: ispec_is_inner
+ logical,intent(in) :: phase_is_inner
! local parameters
real(kind=CUSTOM_REAL) :: pressure
@@ -98,37 +99,35 @@
iglob = ibool(i,j,k,ispec)
! acoustic pressure on global point
- if (PML_CONDITIONS .and. NSPEC_CPML > 0)then
- if (.not. backward_simulation)then
- if (is_CPML(ispec))then
- if (SIMULATION_TYPE == 1)then
- ispec_CPML = spec_to_CPML(ispec)
- call pml_compute_memory_variables_elastic_acoustic(ispec_CPML,iface,iglob,i,j,k,&
- pressure,potential_acoustic,potential_acoustic_old,&
- potential_dot_acoustic,potential_dot_dot_acoustic, &
- num_coupling_ac_el_faces,rmemory_coupling_el_ac_potential_dot_dot)
- pressure = - pressure
- endif
-
- if (SIMULATION_TYPE == 3)then
- ispec_CPML = spec_to_CPML(ispec)
- call pml_compute_memory_variables_elastic_acoustic(ispec_CPML,iface,iglob,i,j,k,&
- pressure,potential_acoustic,potential_acoustic_old,&
- potential_dot_acoustic,potential_dot_dot_acoustic,&
- num_coupling_ac_el_faces,rmemory_coupling_el_ac_potential_dot_dot)
- endif
- else
- pressure = - potential_dot_dot_acoustic(iglob)
+ pressure = - potential_dot_dot_acoustic(iglob)
+
+ ! adjoint wavefield case
+ if (SIMULATION_TYPE /= 1 .and. (.not. backward_simulation)) then
+ ! handles adjoint runs coupling between adjoint potential and adjoint elastic wavefield
+ ! adjoint definition: pressure^\dagger = potential^\dagger
+ pressure = - pressure
+ endif
+
+ ! CPML overwrite cases
+ if (PML_CONDITIONS .and. (.not. backward_simulation) .and. NSPEC_CPML > 0) then
+ if (is_CPML(ispec))then
+ if (SIMULATION_TYPE == 1) then
+ ispec_CPML = spec_to_CPML(ispec)
+ call pml_compute_memory_variables_elastic_acoustic(ispec_CPML,iface,iglob,i,j,k,&
+ pressure,potential_acoustic,potential_acoustic_old,&
+ potential_dot_acoustic,potential_dot_dot_acoustic, &
+ num_coupling_ac_el_faces,rmemory_coupling_el_ac_potential_dot_dot)
+ pressure = - pressure
endif
- else
- if (is_CPML(ispec))then
-! left blank, since no operation needed
- else
- pressure = - potential_dot_dot_acoustic(iglob)
+
+ if (SIMULATION_TYPE == 3)then
+ ispec_CPML = spec_to_CPML(ispec)
+ call pml_compute_memory_variables_elastic_acoustic(ispec_CPML,iface,iglob,i,j,k,&
+ pressure,potential_acoustic,potential_acoustic_old,&
+ potential_dot_acoustic,potential_dot_dot_acoustic,&
+ num_coupling_ac_el_faces,rmemory_coupling_el_ac_potential_dot_dot)
endif
endif
- else
- pressure = - potential_dot_dot_acoustic(iglob)
endif
! gets associated normal on GLL point
@@ -180,7 +179,7 @@ end subroutine compute_coupling_viscoelastic_ac
implicit none
- integer :: NSPEC_AB,NGLOB_AB
+ integer,intent(in) :: NSPEC_AB,NGLOB_AB
real(kind=CUSTOM_REAL),dimension(NDIM,NGLOB_AB),intent(inout) :: accel
real(kind=CUSTOM_REAL),dimension(NGLOB_AB),intent(in) :: rmassx,rmassy,rmassz
@@ -189,10 +188,10 @@ end subroutine compute_coupling_viscoelastic_ac
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_AB),intent(in) :: ibool
! free surface
- integer :: num_free_surface_faces
- real(kind=CUSTOM_REAL) :: free_surface_normal(NDIM,NGLLSQUARE,num_free_surface_faces)
- integer :: free_surface_ijk(3,NGLLSQUARE,num_free_surface_faces)
- integer :: free_surface_ispec(num_free_surface_faces)
+ integer,intent(in) :: num_free_surface_faces
+ real(kind=CUSTOM_REAL),intent(in) :: free_surface_normal(NDIM,NGLLSQUARE,num_free_surface_faces)
+ integer,intent(in) :: free_surface_ijk(3,NGLLSQUARE,num_free_surface_faces)
+ integer,intent(in) :: free_surface_ispec(num_free_surface_faces)
! local parameters
real(kind=CUSTOM_REAL) :: nx,ny,nz
diff --git a/src/specfem3D/compute_forces_acoustic_calling_routine.f90 b/src/specfem3D/compute_forces_acoustic_calling_routine.f90
index 0bc9b6d..cb93058 100644
--- a/src/specfem3D/compute_forces_acoustic_calling_routine.f90
+++ b/src/specfem3D/compute_forces_acoustic_calling_routine.f90
@@ -136,11 +136,7 @@ subroutine compute_forces_acoustic()
! handles adjoint runs coupling between adjoint potential and adjoint elastic wavefield
! adjoint definition: \partial_t^2 \bfs^\dagger=-\frac{1}{\rho}\bfnabla\phi^\dagger
call compute_coupling_acoustic_el(NSPEC_AB,NGLOB_AB, &
-!! DK DK: beware: function or procedure arguments that contain a calculation produce a memory copy
-!! DK DK: created by the compiler; The code is fine, but the hidden copy may slow it down.
-!! DK DK: here is the warning from the Cray compiler:
-!! DK DK: ftn-1438 crayftn: This argument produces a copy in to a temporary variable.
- ibool,-accel,potential_dot_dot_acoustic, &
+ ibool,accel,potential_dot_dot_acoustic, &
num_coupling_ac_el_faces, &
coupling_ac_el_ispec,coupling_ac_el_ijk, &
coupling_ac_el_normal, &
diff --git a/src/specfem3D/compute_forces_viscoelastic_Dev.F90 b/src/specfem3D/compute_forces_viscoelastic_Dev.F90
index eded6af..5c1a9f0 100644
--- a/src/specfem3D/compute_forces_viscoelastic_Dev.F90
+++ b/src/specfem3D/compute_forces_viscoelastic_Dev.F90
@@ -37,7 +37,7 @@
xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz, &
hprime_xx,hprime_xxT, &
hprimewgll_xx,hprimewgll_xxT, &
- wgllwgll_xy,wgllwgll_xz,wgllwgll_yz, &
+ wgllwgll_xy_3D,wgllwgll_xz_3D,wgllwgll_yz_3D, &
kappastore,mustore,jacobian,ibool, &
ATTENUATION,deltat,PML_CONDITIONS, &
one_minus_sum_beta,factor_common,&
@@ -69,7 +69,6 @@
use fault_solver_dynamic, only : Kelvin_Voigt_eta
use specfem_par, only: FULL_ATTENUATION_SOLID
- use specfem_par, only: wgllwgll_xy_3D,wgllwgll_xz_3D,wgllwgll_yz_3D
use pml_par, only: is_CPML, spec_to_CPML, accel_elastic_CPML,NSPEC_CPML, &
PML_dux_dxl, PML_dux_dyl, PML_dux_dzl, PML_duy_dxl, PML_duy_dyl, PML_duy_dzl, &
@@ -105,9 +104,9 @@
! array with derivatives of Lagrange polynomials and precalculated products
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLX) :: hprime_xx,hprimewgll_xxT
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLX) :: hprime_xxT,hprimewgll_xx
- real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY) :: wgllwgll_xy
- real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ) :: wgllwgll_xz
- real(kind=CUSTOM_REAL), dimension(NGLLY,NGLLZ) :: wgllwgll_yz
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: wgllwgll_xy_3D
+ real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLZ,NGLLZ) :: wgllwgll_xz_3D
+ real(kind=CUSTOM_REAL), dimension(NGLLY,NGLLZ,NGLLZ) :: wgllwgll_yz_3D
! memory variables and standard linear solids for attenuation
logical :: ATTENUATION
diff --git a/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90 b/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90
index b897e34..277c656 100644
--- a/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90
+++ b/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90
@@ -141,13 +141,9 @@ subroutine compute_forces_viscoelastic()
else
! handles adjoint runs coupling between adjoint potential and adjoint elastic wavefield
- ! adoint definition: pressure^\dagger=potential^\dagger
+ ! adjoint definition: pressure^\dagger = potential^\dagger
call compute_coupling_viscoelastic_ac(NSPEC_AB,NGLOB_AB, &
-!! DK DK: beware: function or procedure arguments that contain a calculation produce a memory copy
-!! DK DK: created by the compiler; The code is fine, but the hidden copy may slow it down.
-!! DK DK: here is the warning from the Cray compiler:
-!! DK DK: ftn-1438 crayftn: This argument produces a copy in to a temporary variable.
- ibool,accel,-potential_acoustic, &
+ ibool,accel,potential_acoustic, &
num_coupling_ac_el_faces, &
coupling_ac_el_ispec,coupling_ac_el_ijk, &
coupling_ac_el_normal, &
@@ -543,7 +539,7 @@ subroutine compute_forces_viscoelastic_Dev_sim1(iphase)
call compute_forces_viscoelastic_Dev_5p(iphase, NSPEC_AB,NGLOB_AB,displ,veloc,accel, &
xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz, &
hprime_xx,hprime_xxT,hprimewgll_xx,hprimewgll_xxT, &
- wgllwgll_xy,wgllwgll_xz,wgllwgll_yz, &
+ wgllwgll_xy_3D,wgllwgll_xz_3D,wgllwgll_yz_3D, &
kappastore,mustore,jacobian,ibool, &
ATTENUATION,deltat,PML_CONDITIONS, &
one_minus_sum_beta,factor_common, &
@@ -599,7 +595,7 @@ subroutine compute_forces_viscoelastic_Dev_sim3(iphase)
b_displ,b_veloc,b_accel, &
xix,xiy,xiz,etax,etay,etaz,gammax,gammay,gammaz, &
hprime_xx,hprime_xxT,hprimewgll_xx,hprimewgll_xxT, &
- wgllwgll_xy,wgllwgll_xz,wgllwgll_yz, &
+ wgllwgll_xy_3D,wgllwgll_xz_3D,wgllwgll_yz_3D, &
kappastore,mustore,jacobian,ibool, &
ATTENUATION,deltat,PML_CONDITIONS, &
one_minus_sum_beta,factor_common, &
diff --git a/src/specfem3D/fault_solver_dynamic.f90 b/src/specfem3D/fault_solver_dynamic.f90
index 31792f8..c94b810 100644
--- a/src/specfem3D/fault_solver_dynamic.f90
+++ b/src/specfem3D/fault_solver_dynamic.f90
@@ -311,6 +311,8 @@ subroutine init_2d_distribution(a,coord,iin,n)
character(len=MAX_STRING_LEN) :: shapeval
real(kind=CUSTOM_REAL) :: val,valh, xc, yc, zc, r, l, lx,ly,lz
real(kind=CUSTOM_REAL) :: r1(size(a))
+ real(kind=CUSTOM_REAL) :: tmp1(size(a)),tmp2(size(a)),tmp3(size(a))
+
integer :: i
real(kind=CUSTOM_REAL) :: SMALLVAL
@@ -334,46 +336,48 @@ subroutine init_2d_distribution(a,coord,iin,n)
lz = 0e0_CUSTOM_REAL
read(iin,DIST2D)
+
select case (shapeval)
case ('circle')
-!! DK DK: beware: function or procedure arguments that contain a calculation produce a memory copy
-!! DK DK: created by the compiler; The code is fine, but the hidden copy may slow it down.
-!! DK DK: here is the warning from the Cray compiler:
-!! DK DK: ftn-1438 crayftn: CAUTION INIT_2D_DISTRIBUTION, File = src/specfem3D/fault_solver_dynamic.f90, Line = 355, Column = 24
-!! DK DK: This argument produces a copy in to a temporary variable.
- b = heaviside( r - sqrt((coord(1,:)-xc)**2 + (coord(2,:)-yc)**2 + (coord(3,:)-zc)**2) ) *val
+ tmp1 = r - sqrt((coord(1,:)-xc)**2 + (coord(2,:)-yc)**2 + (coord(3,:)-zc)**2)
+ b = heaviside( tmp1 ) * val
+
case ('circle-exp')
r1 = sqrt((coord(1,:)-xc)**2 + (coord(2,:)-yc)**2 + (coord(3,:)-zc)**2)
where(r1<r)
- b =exp(r1**2/(r1**2 - r**2) ) *val + valh
+ b = exp(r1**2/(r1**2 - r**2) ) * val + valh
elsewhere
- b =0._CUSTOM_REAL
+ b = 0._CUSTOM_REAL
endwhere
+
case ('ellipse')
-!! DK DK: beware: function or procedure arguments that contain a calculation produce a memory copy
-!! DK DK: created by the compiler; The code is fine, but the hidden copy may slow it down.
- b = heaviside( 1e0_CUSTOM_REAL - sqrt( (coord(1,:)-xc)**2/lx**2 + (coord(2,:)-yc)**2/ly**2 + (coord(3,:)-zc)**2/lz**2) ) *val
+ tmp1 = 1.0_CUSTOM_REAL - sqrt( (coord(1,:)-xc)**2/lx**2 + (coord(2,:)-yc)**2/ly**2 + (coord(3,:)-zc)**2/lz**2)
+ b = heaviside( tmp1 ) * val
+
case ('square')
-!! DK DK: beware: function or procedure arguments that contain a calculation produce a memory copy
-!! DK DK: created by the compiler; The code is fine, but the hidden copy may slow it down.
- b = heaviside((l/2._CUSTOM_REAL)-abs(coord(1,:)-xc)+SMALLVAL) * &
- heaviside((l/2._CUSTOM_REAL)-abs(coord(2,:)-yc)+SMALLVAL) * &
- heaviside((l/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL) * &
- val
+ tmp1 = (l/2._CUSTOM_REAL)-abs(coord(1,:)-xc)+SMALLVAL
+ tmp2 = (l/2._CUSTOM_REAL)-abs(coord(2,:)-yc)+SMALLVAL
+ tmp3 = (l/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL
+ b = heaviside( tmp1 ) * heaviside( tmp2 ) * heaviside( tmp3) * val
+
case ('cylinder')
- b = heaviside(r - sqrt((coord(1,:)-xc)**2 + (coord(2,:)-yc)**2)) * &
- heaviside((lz/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL) * &
- val
+ tmp1 = r - sqrt((coord(1,:)-xc)**2 + (coord(2,:)-yc)**2)
+ tmp2 = (lz/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL
+ b = heaviside( tmp1 ) * heaviside( tmp2 ) * val
+
case ('rectangle')
- b = heaviside((lx/2._CUSTOM_REAL)-abs(coord(1,:)-xc)+SMALLVAL) * &
- heaviside((ly/2._CUSTOM_REAL)-abs(coord(2,:)-yc)+SMALLVAL) * &
- heaviside((lz/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL) * &
- val
+ tmp1 = (lx/2._CUSTOM_REAL)-abs(coord(1,:)-xc)+SMALLVAL
+ tmp2 = (ly/2._CUSTOM_REAL)-abs(coord(2,:)-yc)+SMALLVAL
+ tmp3 = (lz/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL
+ b = heaviside( tmp1 ) * heaviside( tmp2 ) * heaviside( tmp3 ) * val
+
case ('rectangle-taper')
- b = heaviside((lx/2._CUSTOM_REAL)-abs(coord(1,:)-xc)+SMALLVAL) * &
- heaviside((ly/2._CUSTOM_REAL)-abs(coord(2,:)-yc)+SMALLVAL) * &
- heaviside((lz/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL) * &
- (val + ( coord(3,:) - zc + lz/2._CUSTOM_REAL ) * (valh-val)/lz )
+ tmp1 = (lx/2._CUSTOM_REAL)-abs(coord(1,:)-xc)+SMALLVAL
+ tmp2 = (ly/2._CUSTOM_REAL)-abs(coord(2,:)-yc)+SMALLVAL
+ tmp3 = (lz/2._CUSTOM_REAL)-abs(coord(3,:)-zc)+SMALLVAL
+ b = heaviside( tmp1 ) * heaviside( tmp2 ) * heaviside( tmp3 ) &
+ * (val + ( coord(3,:) - zc + lz/2._CUSTOM_REAL ) * (valh-val)/lz )
+
case default
stop 'bc_dynflt_3d::init_2d_distribution:: unknown shape'
end select
@@ -428,10 +432,10 @@ subroutine BC_DYNFLT_set3d(bc,MxA,V,D,iflt)
real(kind=CUSTOM_REAL), intent(in) :: V(:,:),D(:,:)
integer, intent(in) :: iflt
+ ! local parameters
real(kind=CUSTOM_REAL), dimension(3,bc%nglob) :: T,dD,dV,dA
real(kind=CUSTOM_REAL), dimension(bc%nglob) :: strength,tStick,tnew, &
- theta_old, theta_new, dc, &
- Vf_old,Vf_new,TxExt
+ theta_old, theta_new, dc, Vf_old, Vf_new, TxExt, tmp_Vf
real(kind=CUSTOM_REAL) :: half_dt,TLoad,DTau0,GLoad,timeval
integer :: i
@@ -521,7 +525,8 @@ subroutine BC_DYNFLT_set3d(bc,MxA,V,D,iflt)
! second pass
bc%rsf%theta = theta_old
- call rsf_update_state(0.5_CUSTOM_REAL*(Vf_old + Vf_new),bc%dt,bc%rsf)
+ tmp_Vf(:) = 0.5_CUSTOM_REAL*(Vf_old(:) + Vf_new(:))
+ call rsf_update_state(tmp_Vf,bc%dt,bc%rsf)
do i=1,bc%nglob
Vf_new(i)=rtsafe(funcd,0.0_CUSTOM_REAL,Vf_old(i)+5.0_CUSTOM_REAL,1e-5_CUSTOM_REAL,tStick(i),-T(3,i),bc%Z(i),bc%rsf%f0(i), &
bc%rsf%V0(i),bc%rsf%a(i),bc%rsf%b(i),bc%rsf%L(i),bc%rsf%theta(i),bc%rsf%StateLaw)
@@ -1022,6 +1027,10 @@ subroutine init_dataXZ(dataXZ,bc)
allocate(bc%dataXZ_all%sta(npoin_all))
endif
+!note: crayftn compiler warns about possible copy which may slow down the code for dataXZ%npoin,dataXZ%xcoord,..
+!ftn-1438 crayftn: CAUTION INIT_DATAXZ, File = src/specfem3D/fault_solver_dynamic.f90, Line = 1036, Column = 45
+! This argument produces a possible copy in and out to a temporary variable.
+
allocate(bc%npoin_perproc(NPROC))
bc%npoin_perproc=0
call gather_all_singlei(dataXZ%npoin,bc%npoin_perproc,NPROC)
diff --git a/src/specfem3D/pml_compute_memory_variables.f90 b/src/specfem3D/pml_compute_memory_variables.f90
index 46846d3..c6c0247 100644
--- a/src/specfem3D/pml_compute_memory_variables.f90
+++ b/src/specfem3D/pml_compute_memory_variables.f90
@@ -519,10 +519,10 @@ subroutine pml_compute_memory_variables_acoustic_elastic(ispec_CPML,iface,iglob,
implicit none
integer, intent(in) :: ispec_CPML,iface,iglob,num_coupling_ac_el_faces
- real(kind=CUSTOM_REAL) :: displ_x,displ_y,displ_z
+ real(kind=CUSTOM_REAL),intent(out) :: displ_x,displ_y,displ_z
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_AB), intent(in) :: displ
- real(kind=CUSTOM_REAL), dimension(3,NGLLX,NGLLY,NGLLZ,num_coupling_ac_el_faces,2) :: &
- rmemory_coupling_ac_el_displ
+ real(kind=CUSTOM_REAL), dimension(3,NGLLX,NGLLY,NGLLZ,num_coupling_ac_el_faces,2),intent(inout) :: &
+ rmemory_coupling_ac_el_displ
! local parameters
integer :: i,j,k,CPML_region_local,singularity_type_2
@@ -603,7 +603,7 @@ subroutine pml_compute_memory_variables_elastic_acoustic(ispec_CPML,iface,iglob,
implicit none
integer, intent(in) :: ispec_CPML,iface,iglob,num_coupling_ac_el_faces
- real(kind=CUSTOM_REAL) :: pressure
+ real(kind=CUSTOM_REAL),intent(out) :: pressure
real(kind=CUSTOM_REAL), dimension(NGLOB_AB), intent(in) :: potential_acoustic,potential_acoustic_old, &
potential_dot_acoustic,potential_dot_dot_acoustic
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,num_coupling_ac_el_faces,2) :: &
More information about the CIG-COMMITS
mailing list