[cig-commits] r18959 - seismo/3D/SPECFEM3D/trunk/src/specfem3D
yangl at geodynamics.org
yangl at geodynamics.org
Wed Sep 21 12:02:03 PDT 2011
Author: yangl
Date: 2011-09-21 12:02:02 -0700 (Wed, 21 Sep 2011)
New Revision: 18959
Modified:
seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_acoustic.f90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_elastic.f90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/iterate_time.f90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/read_mesh_databases.f90
seismo/3D/SPECFEM3D/trunk/src/specfem3D/specfem3D_par.f90
Log:
(SPECFEM3D_SESAME) new adjoint potential/displacement definition and coupling terms in acoustic-elastic-coupled adjoint simulations
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_acoustic.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_acoustic.f90 2011-09-21 16:21:40 UTC (rev 18958)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_acoustic.f90 2011-09-21 19:02:02 UTC (rev 18959)
@@ -172,13 +172,25 @@
! elastic coupling
if(ELASTIC_SIMULATION ) then
- call compute_coupling_acoustic_el(NSPEC_AB,NGLOB_AB, &
- ibool,displ,potential_dot_dot_acoustic, &
- num_coupling_ac_el_faces, &
- coupling_ac_el_ispec,coupling_ac_el_ijk, &
- coupling_ac_el_normal, &
- coupling_ac_el_jacobian2Dw, &
- ispec_is_inner,phase_is_inner)
+ if( SIMULATION_TYPE == 1 ) then
+ ! forward definition: \bfs=\frac{1}{\rho}\bfnabla\phi
+ call compute_coupling_acoustic_el(NSPEC_AB,NGLOB_AB, &
+ ibool,displ,potential_dot_dot_acoustic, &
+ num_coupling_ac_el_faces, &
+ coupling_ac_el_ispec,coupling_ac_el_ijk, &
+ coupling_ac_el_normal, &
+ coupling_ac_el_jacobian2Dw, &
+ ispec_is_inner,phase_is_inner)
+ else
+ ! adjoint definition: \partial_t^2 \bfs^\dagger=-\frac{1}{\rho}\bfnabla\phi^\dagger
+ call compute_coupling_acoustic_el(NSPEC_AB,NGLOB_AB, &
+ ibool,-accel_adj_coupling,potential_dot_dot_acoustic, &
+ num_coupling_ac_el_faces, &
+ coupling_ac_el_ispec,coupling_ac_el_ijk, &
+ coupling_ac_el_normal, &
+ coupling_ac_el_jacobian2Dw, &
+ ispec_is_inner,phase_is_inner)
+ endif
! adjoint simulations
if( SIMULATION_TYPE == 3 ) &
call compute_coupling_acoustic_el(NSPEC_ADJOINT,NGLOB_ADJOINT, &
@@ -308,6 +320,10 @@
ibool,free_surface_ijk,free_surface_ispec, &
num_free_surface_faces,ispec_is_acoustic)
+ potential_acoustic_adj_coupling(:) = potential_acoustic(:) &
+ + deltat * potential_dot_acoustic(:) &
+ + deltatsqover2 * potential_dot_dot_acoustic(:)
+
! adjoint simulations
if (SIMULATION_TYPE == 3) &
call acoustic_enforce_free_surface(NSPEC_AB,NGLOB_ADJOINT, &
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_elastic.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_elastic.f90 2011-09-21 16:21:40 UTC (rev 18958)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/compute_forces_elastic.f90 2011-09-21 19:02:02 UTC (rev 18959)
@@ -445,13 +445,25 @@
! acoustic coupling
if( ACOUSTIC_SIMULATION ) then
- call compute_coupling_elastic_ac(NSPEC_AB,NGLOB_AB, &
- ibool,accel,potential_dot_dot_acoustic, &
- num_coupling_ac_el_faces, &
- coupling_ac_el_ispec,coupling_ac_el_ijk, &
- coupling_ac_el_normal, &
- coupling_ac_el_jacobian2Dw, &
- ispec_is_inner,phase_is_inner)
+ if( SIMULATION_TYPE == 1 ) then
+ ! forward definition: pressure=-potential_dot_dot
+ call compute_coupling_elastic_ac(NSPEC_AB,NGLOB_AB, &
+ ibool,accel,potential_dot_dot_acoustic, &
+ num_coupling_ac_el_faces, &
+ coupling_ac_el_ispec,coupling_ac_el_ijk, &
+ coupling_ac_el_normal, &
+ coupling_ac_el_jacobian2Dw, &
+ ispec_is_inner,phase_is_inner)
+ else
+ ! adoint definition: pressure^\dagger=potential^\dagger
+ call compute_coupling_elastic_ac(NSPEC_AB,NGLOB_AB, &
+ ibool,accel,-potential_acoustic_adj_coupling, &
+ num_coupling_ac_el_faces, &
+ coupling_ac_el_ispec,coupling_ac_el_ijk, &
+ coupling_ac_el_normal, &
+ coupling_ac_el_jacobian2Dw, &
+ ispec_is_inner,phase_is_inner)
+ endif
! adjoint simulations
if( SIMULATION_TYPE == 3 ) &
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/iterate_time.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/iterate_time.f90 2011-09-21 16:21:40 UTC (rev 18958)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/iterate_time.f90 2011-09-21 19:02:02 UTC (rev 18959)
@@ -324,6 +324,7 @@
if( ELASTIC_SIMULATION ) then
displ(:,:) = displ(:,:) + deltat*veloc(:,:) + deltatsqover2*accel(:,:)
veloc(:,:) = veloc(:,:) + deltatover2*accel(:,:)
+ accel_adj_coupling(:,:) = accel(:,:)
accel(:,:) = 0._CUSTOM_REAL
endif
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/read_mesh_databases.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/read_mesh_databases.f90 2011-09-21 16:21:40 UTC (rev 18958)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/read_mesh_databases.f90 2011-09-21 19:02:02 UTC (rev 18959)
@@ -86,6 +86,8 @@
if( ier /= 0 ) stop 'error allocating array potential_dot_acoustic'
allocate(potential_dot_dot_acoustic(NGLOB_AB),stat=ier)
if( ier /= 0 ) stop 'error allocating array potential_dot_dot_acoustic'
+ allocate(potential_acoustic_adj_coupling(NGLOB_AB),stat=ier)
+ if( ier /= 0 ) stop 'error allocating array potential_acoustic_adj_coupling'
! mass matrix, density
allocate(rmass_acoustic(NGLOB_AB),stat=ier)
@@ -107,6 +109,8 @@
if( ier /= 0 ) stop 'error allocating array veloc'
allocate(accel(NDIM,NGLOB_AB),stat=ier)
if( ier /= 0 ) stop 'error allocating array accel'
+ allocate(accel_adj_coupling(NDIM,NGLOB_AB),stat=ier)
+ if( ier /= 0 ) stop 'error allocating array accel_adj_coupling'
allocate(rmass(NGLOB_AB),stat=ier)
if( ier /= 0 ) stop 'error allocating array rmass'
Modified: seismo/3D/SPECFEM3D/trunk/src/specfem3D/specfem3D_par.f90
===================================================================
--- seismo/3D/SPECFEM3D/trunk/src/specfem3D/specfem3D_par.f90 2011-09-21 16:21:40 UTC (rev 18958)
+++ seismo/3D/SPECFEM3D/trunk/src/specfem3D/specfem3D_par.f90 2011-09-21 19:02:02 UTC (rev 18959)
@@ -273,6 +273,7 @@
! displacement, velocity, acceleration
real(kind=CUSTOM_REAL), dimension(:,:), allocatable :: displ,veloc,accel
+ real(kind=CUSTOM_REAL), dimension(:,:), allocatable :: accel_adj_coupling
! mass matrix
real(kind=CUSTOM_REAL), dimension(:), allocatable :: rmass
@@ -351,6 +352,7 @@
! potential
real(kind=CUSTOM_REAL), dimension(:), allocatable :: potential_acoustic, &
potential_dot_acoustic,potential_dot_dot_acoustic
+ real(kind=CUSTOM_REAL), dimension(:), allocatable :: potential_acoustic_adj_coupling
! density
real(kind=CUSTOM_REAL), dimension(:,:,:,:), allocatable :: rhostore
More information about the CIG-COMMITS
mailing list