[commit] devel: merged Vadim's support for coupling with DSM, cleaned and improved by Clément Durochat (01063d6)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Tue Jul 1 17:28:09 PDT 2014
Repository : https://github.com/geodynamics/specfem3d
On branch : devel
Link : https://github.com/geodynamics/specfem3d/compare/35ff4e4b7f0c54c55b1a26cdbe5c6cd7d5dd3b55...271c1ec41d525aef92e78a9894fdef484d9927e8
>---------------------------------------------------------------
commit 01063d67149b5eae628c9e9e1eb7c29984ef037a
Author: Clément Durochat <durochat at lma.cnrs-mrs.fr>
Date: Tue Jul 1 19:22:28 2014 +0200
merged Vadim's support for coupling with DSM, cleaned and improved by Clément Durochat
>---------------------------------------------------------------
01063d67149b5eae628c9e9e1eb7c29984ef037a
setup/constants.h.in | 15 +++--
src/generate_databases/get_absorbing_boundary.f90 | 64 +++++++++++++++++++++-
src/generate_databases/save_arrays_solver.f90 | 35 ++++++++++++
src/shared/read_parameter_file.f90 | 2 +-
src/specfem3D/compute_add_sources_viscoelastic.f90 | 15 +++--
...compute_forces_viscoelastic_calling_routine.F90 | 22 +++++---
src/specfem3D/compute_stacey_viscoelastic.f90 | 46 ++++++++++++----
src/specfem3D/initialize_simulation.f90 | 8 ++-
src/specfem3D/prepare_timerun.F90 | 6 +-
src/specfem3D/read_mesh_databases.F90 | 5 +-
src/specfem3D/read_mesh_databases_adios.F90 | 5 +-
src/specfem3D/specfem3D_par.f90 | 2 +-
12 files changed, 185 insertions(+), 40 deletions(-)
diff --git a/setup/constants.h.in b/setup/constants.h.in
index a98f322..001d1a0 100644
--- a/setup/constants.h.in
+++ b/setup/constants.h.in
@@ -157,6 +157,17 @@
!!-----------------------------------------------------------
!!
+!! For coupling with DSM by VM
+!!
+!!-----------------------------------------------------------
+
+ logical, parameter :: COUPLE_WITH_DSM = .false. !!! .true.
+
+! some old tests (currently unstable; do not remove them though, we might fix this one day)
+ integer, parameter :: IIN_veloc_dsm = 51, IIN_tract_dsm = 52, Ntime_step_dsm = 100
+
+!!-----------------------------------------------------------
+!!
!! source/receiver setup
!!
!!-----------------------------------------------------------
@@ -477,7 +488,3 @@
integer, parameter :: IMODEL_IPATI = 10
integer, parameter :: IMODEL_IPATI_WATER = 11
-! some old tests (currently unstable; do not remove them though, we might fix this one day)
- logical, parameter :: OLD_TEST_TO_FIX_ONE_DAY = .false.
- integer, parameter :: IIN_veloc_dsm = 51, IIN_tract_dsm = 52, Ntime_step_dsm = 100
-
diff --git a/src/generate_databases/get_absorbing_boundary.f90 b/src/generate_databases/get_absorbing_boundary.f90
index 40b8ba9..8b2c3fc 100644
--- a/src/generate_databases/get_absorbing_boundary.f90
+++ b/src/generate_databases/get_absorbing_boundary.f90
@@ -42,7 +42,7 @@
implicit none
! number of spectral elements in each block
- integer :: myrank,nspec
+ integer :: myrank,ier,nspec
! arrays with the mesh
integer, dimension(NGLLX,NGLLY,NGLLZ,nspec) :: ibool
@@ -81,6 +81,34 @@
real(kind=CUSTOM_REAL),dimension(NGNOD2D_FOUR_CORNERS) :: xcoord,ycoord,zcoord
integer :: ispec,ispec2D,icorner,itop,iabsval,iface,igll,i,j,igllfree,ifree
+ !! CD modif. : begin (implemented by VM) !! additonal local parameters For coupling with DSM
+
+ logical, dimension(:,:),allocatable :: iboun ! pll
+
+ ! corner locations for faces
+ real(kind=CUSTOM_REAL), dimension(:,:,:),allocatable :: xcoord_iboun,ycoord_iboun,zcoord_iboun
+ character(len=27) namefile
+
+ ! sets flag in array iboun for elements with an absorbing boundary faces
+ if(COUPLE_WITH_DSM) then
+
+ ! allocate temporary flag array
+ allocate(iboun(6,nspec), &
+ xcoord_iboun(NGNOD2D,6,nspec), &
+ ycoord_iboun(NGNOD2D,6,nspec), &
+ zcoord_iboun(NGNOD2D,6,nspec),stat=ier)
+ if(ier /= 0) stop 'not enough memory to allocate arrays'
+
+ iboun(:,:) = .false.
+
+ write(namefile,'(a17,i6.6,a4)') 'xmin_gll_for_dsm_',myrank,'.txt'
+ open(123,file=namefile)
+ write(123,*) nspec2D_xmin
+
+ endif
+
+ !! CD modif. : end
+
! abs face counter
iabsval = 0
@@ -106,6 +134,10 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
iface)
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) iboun(iface,ispec) = .true.
+ !! CD modif. : end
+
! ijk indices of GLL points for face id
call get_element_face_gll_indices(iface,ijk_face,NGLLX,NGLLZ)
@@ -126,6 +158,12 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
lnormal )
normal_face(:,i,j) = lnormal(:)
+
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) write(123,'(i10,3f20.10)') ispec,xstore_dummy(ibool(i,j,1,ispec)),&
+ ystore_dummy(ibool(i,j,1,ispec)),zstore_dummy(ibool(i,j,1,ispec))
+ !! CD modif. : end
+
enddo
enddo
@@ -146,6 +184,10 @@
enddo ! nspec2D_xmin
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) close(123)
+ !! CD modif. : end
+
! xmax
ijk_face(:,:,:) = 0
normal_face(:,:,:) = 0.0_CUSTOM_REAL
@@ -168,6 +210,10 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
iface )
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) iboun(iface,ispec) = .true.
+ !! CD modif. : end
+
! ijk indices of GLL points on face
call get_element_face_gll_indices(iface,ijk_face,NGLLX,NGLLZ)
@@ -230,6 +276,10 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
iface )
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) iboun(iface,ispec) = .true.
+ !! CD modif. : end
+
! ijk indices of GLL points on face
call get_element_face_gll_indices(iface,ijk_face,NGLLY,NGLLZ)
@@ -292,6 +342,10 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
iface )
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) iboun(iface,ispec) = .true.
+ !! CD modif. : end
+
! ijk indices of GLL points on face
call get_element_face_gll_indices(iface,ijk_face,NGLLY,NGLLZ)
@@ -354,6 +408,10 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
iface )
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) iboun(iface,ispec) = .true.
+ !! CD modif. : end
+
! ijk indices of GLL points on face
call get_element_face_gll_indices(iface,ijk_face,NGLLX,NGLLY)
@@ -419,6 +477,10 @@
xstore_dummy,ystore_dummy,zstore_dummy, &
iface )
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(COUPLE_WITH_DSM) iboun(iface,ispec) = .true.
+ !! CD modif. : end
+
! ijk indices of GLL points on face
call get_element_face_gll_indices(iface,ijk_face,NGLLX,NGLLY)
diff --git a/src/generate_databases/save_arrays_solver.f90 b/src/generate_databases/save_arrays_solver.f90
index 95a69f3..8cc2e5f 100644
--- a/src/generate_databases/save_arrays_solver.f90
+++ b/src/generate_databases/save_arrays_solver.f90
@@ -571,6 +571,41 @@
deallocate(iglob_tmp,v_tmp_i)
endif
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
+ !if (num_abs_boundary_faces > 0) then
+ filename = prname(1:len_trim(prname))//'absorb_dsm'
+ open(27,file=filename(1:len_trim(filename)),status='unknown',form='unformatted',iostat=ier)
+ if( ier /= 0 ) stop 'error opening file absorb_dsm'
+ write(27) num_abs_boundary_faces
+ write(27) abs_boundary_ispec
+ write(27) abs_boundary_ijk
+ write(27) abs_boundary_jacobian2Dw
+ write(27) abs_boundary_normal
+ close(27)
+
+ filename = prname(1:len_trim(prname))//'inner'
+ open(27,file=filename(1:len_trim(filename)),status='unknown',form='unformatted',iostat=ier)
+ write(27) ispec_is_inner
+ write(27) ispec_is_elastic
+ close(27)
+
+ !end if
+
+ !! Don't delete this comment for the moment
+ !!
+ !! saves 1. MPI interface
+ !!
+ !!if( num_interfaces_ext_mesh >= 1 ) then
+ !! filename = prname(1:len_trim(prname))//'MPI_1_points'
+ !! call write_VTK_data_points(nglob, xstore_dummy,ystore_dummy,zstore_dummy, &
+ !! ibool_interfaces_ext_mesh_dummy(1:nibool_interfaces_ext_mesh(1),1), &
+ !! nibool_interfaces_ext_mesh(1), filename)
+ !!endif
+
+ endif
+ !! CD modif. : end
+
! acoustic-poroelastic domains
if( ACOUSTIC_SIMULATION .and. POROELASTIC_SIMULATION ) then
! saves points on acoustic-poroelastic coupling interface
diff --git a/src/shared/read_parameter_file.f90 b/src/shared/read_parameter_file.f90
index 758cd67..d257bd4 100644
--- a/src/shared/read_parameter_file.f90
+++ b/src/shared/read_parameter_file.f90
@@ -170,7 +170,7 @@
if (ierr /= 0) return
!! read the traction path directory
- if (OLD_TEST_TO_FIX_ONE_DAY) then
+ if (COUPLE_WITH_DSM) then
call read_value_string(TRAC_PATH, 'TRAC_PATH', ierr)
if (ierr /= 0) return
endif
diff --git a/src/specfem3D/compute_add_sources_viscoelastic.f90 b/src/specfem3D/compute_add_sources_viscoelastic.f90
index 0b0a25c..67de68f 100644
--- a/src/specfem3D/compute_add_sources_viscoelastic.f90
+++ b/src/specfem3D/compute_add_sources_viscoelastic.f90
@@ -103,8 +103,9 @@
!equivalence (i2head,i4head,r4head) ! share the same 240-byte-memory
double precision :: hxir(NGLLX),hpxir(NGLLX),hetar(NGLLY),hpetar(NGLLY),hgammar(NGLLZ),hpgammar(NGLLZ)
- ! some old tests (currently unstable; do not remove them though, we might fix this one day)
- if (OLD_TEST_TO_FIX_ONE_DAY) return
+ !! CD modif (implemented by VM) : !! For coupling with DSM
+ if (COUPLE_WITH_DSM) return
+ !! CD modif : end
! plotting source time function
if(PRINT_SOURCE_TIME_FUNCTION .and. .not. phase_is_inner ) then
@@ -446,8 +447,9 @@
real(kind=CUSTOM_REAL) stf_used,stf_used_total_all,time_source
integer :: isource,iglob,i,j,k,ispec
- ! some old tests (currently unstable; do not remove them though, we might fix this one day)
- if (OLD_TEST_TO_FIX_ONE_DAY) return
+ !! CD modif (implemented by VM) : !! For coupling with DSM
+ if (COUPLE_WITH_DSM) return
+ !! CD modif : end
! plotting source time function
if(PRINT_SOURCE_TIME_FUNCTION .and. .not. phase_is_inner ) then
@@ -666,8 +668,9 @@
!equivalence (i2head,i4head,r4head) ! share the same 240-byte-memory
double precision :: hxir(NGLLX),hpxir(NGLLX),hetar(NGLLY),hpetar(NGLLY),hgammar(NGLLZ),hpgammar(NGLLZ)
- ! some old tests (currently unstable; do not remove them though, we might fix this one day)
- if (OLD_TEST_TO_FIX_ONE_DAY) return
+ !! CD modif (implemented by VM) : !! For coupling with DSM
+ if (COUPLE_WITH_DSM) return
+ !! CD modif : end
! plotting source time function
if(PRINT_SOURCE_TIME_FUNCTION .and. .not. phase_is_inner ) then
diff --git a/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90 b/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90
index 2fbfec5..ee90a80 100644
--- a/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90
+++ b/src/specfem3D/compute_forces_viscoelastic_calling_routine.F90
@@ -170,15 +170,19 @@ subroutine compute_forces_viscoelastic()
ispec_is_inner,phase_is_inner)
endif
-! adds source term (single-force/moment-tensor solution)
- call compute_add_sources_viscoelastic( NSPEC_AB,NGLOB_AB,accel, &
- ibool,ispec_is_inner,phase_is_inner, &
- NSOURCES,myrank,it,islice_selected_source,ispec_selected_source,&
- hdur,hdur_gaussian,tshift_src,dt,t0,sourcearrays, &
- ispec_is_elastic,SIMULATION_TYPE,NSTEP, &
- nrec,islice_selected_rec,ispec_selected_rec, &
- nadj_rec_local,adj_sourcearrays, &
- NTSTEP_BETWEEN_READ_ADJSRC,NOISE_TOMOGRAPHY)
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if(.not. COUPLE_WITH_DSM) then
+ ! adds source term (single-force/moment-tensor solution)
+ call compute_add_sources_viscoelastic(NSPEC_AB,NGLOB_AB,accel, &
+ ibool,ispec_is_inner,phase_is_inner, &
+ NSOURCES,myrank,it,islice_selected_source,ispec_selected_source,&
+ hdur,hdur_gaussian,tshift_src,dt,t0,sourcearrays, &
+ ispec_is_elastic,SIMULATION_TYPE,NSTEP, &
+ nrec,islice_selected_rec,ispec_selected_rec, &
+ nadj_rec_local,adj_sourcearrays, &
+ NTSTEP_BETWEEN_READ_ADJSRC,NOISE_TOMOGRAPHY)
+ endif
+ !! CD modif. : end
! assemble all the contributions between slices using MPI
if( phase_is_inner .eqv. .false. ) then
diff --git a/src/specfem3D/compute_stacey_viscoelastic.f90 b/src/specfem3D/compute_stacey_viscoelastic.f90
index 9f4109d..88ad4b0 100644
--- a/src/specfem3D/compute_stacey_viscoelastic.f90
+++ b/src/specfem3D/compute_stacey_viscoelastic.f90
@@ -79,13 +79,15 @@
real(kind=CUSTOM_REAL) vx,vy,vz,nx,ny,nz,tx,ty,tz,vn,jacobianw
integer :: ispec,iglob,i,j,k,iface,igll
-! for new method
+ !! CD modif (implemented by VM) : begin !! For coupling with DSM
+
+ ! See also DSM parameters in setup/constants.h.in
real(kind=CUSTOM_REAL) :: Veloc_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces)
real(kind=CUSTOM_REAL) :: Tract_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces)
integer :: it_dsm
- if (OLD_TEST_TO_FIX_ONE_DAY) then
+ if (COUPLE_WITH_DSM) then
if (phase_is_inner .eqv. .false.) then
if (mod(it_dsm,Ntime_step_dsm+1) == 0 .or. it == 1) then
call read_dsm_file(Veloc_dsm_boundary,Tract_dsm_boundary,num_abs_boundary_faces,it_dsm)
@@ -93,6 +95,8 @@
endif
endif
+ !! CD modif : end
+
! checks if anything to do
if( num_abs_boundary_faces == 0 ) return
@@ -118,11 +122,15 @@
vx=veloc(1,iglob)
vy=veloc(2,iglob)
vz=veloc(3,iglob)
- if (OLD_TEST_TO_FIX_ONE_DAY) then
+
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
vx = vx - Veloc_dsm_boundary(1,it_dsm,igll,iface)
vy = vy - Veloc_dsm_boundary(2,it_dsm,igll,iface)
vz = vz - Veloc_dsm_boundary(3,it_dsm,igll,iface)
endif
+ !! CD modif. : end
+
! gets associated normal
nx = abs_boundary_normal(1,igll,iface)
ny = abs_boundary_normal(2,igll,iface)
@@ -136,11 +144,13 @@
ty = rho_vp(i,j,k,ispec)*vn*ny + rho_vs(i,j,k,ispec)*(vy-vn*ny)
tz = rho_vp(i,j,k,ispec)*vn*nz + rho_vs(i,j,k,ispec)*(vz-vn*nz)
- if (OLD_TEST_TO_FIX_ONE_DAY) then
- tx = tx -Tract_dsm_boundary(1,it_dsm,igll,iface)
- ty = ty -Tract_dsm_boundary(2,it_dsm,igll,iface)
- tz = tz -Tract_dsm_boundary(3,it_dsm,igll,iface)
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
+ tx = tx - Tract_dsm_boundary(1,it_dsm,igll,iface)
+ ty = ty - Tract_dsm_boundary(2,it_dsm,igll,iface)
+ tz = tz - Tract_dsm_boundary(3,it_dsm,igll,iface)
endif
+ !! CD modif. : end
! gets associated, weighted jacobian
jacobianw = abs_boundary_jacobian2Dw(igll,iface)
@@ -170,11 +180,13 @@
endif
endif
- if (OLD_TEST_TO_FIX_ONE_DAY) then
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
if (phase_is_inner .eqv. .true.) then
it_dsm = it_dsm + 1
endif
endif
+ !! CD modif. : end
end subroutine compute_stacey_viscoelastic
!
@@ -264,7 +276,12 @@
end subroutine compute_stacey_viscoelastic_bpwf
-!---------------------------------------------------------------------------------------
+!=============================================================================
+!
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+!
+!-----------------------------------------------------------------------------
+
subroutine read_dsm_file(Veloc_dsm_boundary,Tract_dsm_boundary,num_abs_boundary_faces,it_dsm)
@@ -298,6 +315,7 @@
end subroutine read_dsm_file
+ !! CD modif. : end
!
!=====================================================================
! for elastic solver on GPU
@@ -330,13 +348,13 @@
! GPU_MODE variables
integer(kind=8) :: Mesh_pointer
- ! for new method
+ !! CD modif (implemented by VM) : begin !! For coupling with DSM
real(kind=CUSTOM_REAL) :: Veloc_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces)
real(kind=CUSTOM_REAL) :: Tract_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces)
integer :: it_dsm
- if (OLD_TEST_TO_FIX_ONE_DAY) then
+ if (COUPLE_WITH_DSM) then
if (phase_is_inner .eqv. .false.) then
if (mod(it_dsm,Ntime_step_dsm+1) == 0 .or. it == 1) then
call read_dsm_file(Veloc_dsm_boundary,Tract_dsm_boundary,num_abs_boundary_faces,it_dsm)
@@ -344,6 +362,8 @@
endif
endif
+ !! CD modif. : end
+
! checks if anything to do
if( num_abs_boundary_faces == 0 ) return
@@ -366,11 +386,13 @@
endif
endif
- if (OLD_TEST_TO_FIX_ONE_DAY) then
+ !! CD modif. (implemented by VM) : begin !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
if (phase_is_inner .eqv. .true.) then
it_dsm = it_dsm + 1
endif
endif
+ !! CD modif. : end
end subroutine compute_stacey_viscoelastic_GPU
diff --git a/src/specfem3D/initialize_simulation.f90 b/src/specfem3D/initialize_simulation.f90
index 5aa3217..d09b491 100644
--- a/src/specfem3D/initialize_simulation.f90
+++ b/src/specfem3D/initialize_simulation.f90
@@ -66,6 +66,10 @@
! GPU_MODE is in par_file
call read_gpu_mode(GPU_MODE,GRAVITY)
+!! CD modif. : begin !! For coupling with DSM
+ if(GPU_MODE .and. COUPLE_WITH_DSM) stop 'Coupling with DSM currently not implemented for GPUs'
+!! CD modif. : end
+
! myrank is the rank of each process, between 0 and NPROC-1.
! as usual in MPI, process 0 is in charge of coordinating everything
! and also takes care of the main output
@@ -142,7 +146,9 @@
! reads in numbers of spectral elements and points for the part of the mesh handled by this process
call create_name_database(prname,myrank,LOCAL_PATH)
- if (OLD_TEST_TO_FIX_ONE_DAY) call create_name_database(dsmname,myrank,TRAC_PATH)
+
+! for coupling with DSM
+ if (COUPLE_WITH_DSM) call create_name_database(dsmname,myrank,TRAC_PATH)
if (ADIOS_FOR_MESH) then
call read_mesh_for_init(NSPEC_AB, NGLOB_AB)
diff --git a/src/specfem3D/prepare_timerun.F90 b/src/specfem3D/prepare_timerun.F90
index c27c310..0e23106 100644
--- a/src/specfem3D/prepare_timerun.F90
+++ b/src/specfem3D/prepare_timerun.F90
@@ -275,7 +275,9 @@
if(ELASTIC_SIMULATION) then
! switches to three-component mass matrix
- if( STACEY_ABSORBING_CONDITIONS ) then
+
+ !! CD modif. : begin
+ if( STACEY_ABSORBING_CONDITIONS .and. .not. COUPLE_WITH_DSM ) then
! adds boundary contributions
rmassx(:) = rmass(:) + rmassx(:)
rmassy(:) = rmass(:) + rmassy(:)
@@ -285,6 +287,8 @@
rmassy(:) = rmass(:)
rmassz(:) = rmass(:)
endif
+ !! CD modif. : end
+
! not needed anymore
deallocate(rmass)
diff --git a/src/specfem3D/read_mesh_databases.F90 b/src/specfem3D/read_mesh_databases.F90
index 6f2cd96..e042ea8 100644
--- a/src/specfem3D/read_mesh_databases.F90
+++ b/src/specfem3D/read_mesh_databases.F90
@@ -396,8 +396,8 @@
abs_boundary_normal(NDIM,NGLLSQUARE,num_abs_boundary_faces),stat=ier)
if( ier /= 0 ) stop 'error allocating array abs_boundary_ispec etc.'
- if (OLD_TEST_TO_FIX_ONE_DAY) then
- ! for new method
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
allocate(Veloc_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces))
allocate(Tract_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces))
open(unit=IIN_veloc_dsm,file=dsmname(1:len_trim(dsmname))//'vel.bin',status='old', &
@@ -408,6 +408,7 @@
allocate(Veloc_dsm_boundary(1,1,1,1))
allocate(Tract_dsm_boundary(1,1,1,1))
endif
+ !! CD modif. : end
if(PML_CONDITIONS)then
if (num_abs_boundary_faces > 0) then
diff --git a/src/specfem3D/read_mesh_databases_adios.F90 b/src/specfem3D/read_mesh_databases_adios.F90
index 3c7b0a1..62c07b7 100644
--- a/src/specfem3D/read_mesh_databases_adios.F90
+++ b/src/specfem3D/read_mesh_databases_adios.F90
@@ -853,8 +853,8 @@ subroutine read_mesh_databases_adios()
abs_boundary_normal(NDIM,NGLLSQUARE,num_abs_boundary_faces),stat=ier)
if( ier /= 0 ) stop 'error allocating array abs_boundary_ispec etc.'
- if (OLD_TEST_TO_FIX_ONE_DAY) then
- ! for new method
+ !! CD modif. : begin (implemented by VM) !! For coupling with DSM
+ if (COUPLE_WITH_DSM) then
allocate(Veloc_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces))
allocate(Tract_dsm_boundary(3,Ntime_step_dsm,NGLLSQUARE,num_abs_boundary_faces))
open(unit=IIN_veloc_dsm,file=dsmname(1:len_trim(dsmname))//'vel.bin',status='old', &
@@ -865,6 +865,7 @@ subroutine read_mesh_databases_adios()
allocate(Veloc_dsm_boundary(1,1,1,1))
allocate(Tract_dsm_boundary(1,1,1,1))
endif
+ !! CD modif. : end
allocate(ibelm_xmin(nspec2D_xmin),ibelm_xmax(nspec2D_xmax), &
ibelm_ymin(nspec2D_ymin),ibelm_ymax(nspec2D_ymax), &
diff --git a/src/specfem3D/specfem3D_par.f90 b/src/specfem3D/specfem3D_par.f90
index daf0b83..4cd207e 100644
--- a/src/specfem3D/specfem3D_par.f90
+++ b/src/specfem3D/specfem3D_par.f90
@@ -95,7 +95,7 @@ module specfem_par
! time loop step
integer :: it
-! for new method
+! for coupling with DSM
integer :: it_dsm
! parameters for the source
More information about the CIG-COMMITS
mailing list