[cig-commits] [commit] master: computing time derivative of the stf and dump it to netcf (0661c3b)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Thu Aug 21 05:16:54 PDT 2014
Repository : https://github.com/geodynamics/axisem
On branch : master
Link : https://github.com/geodynamics/axisem/compare/0732e253c0a4b9db6b5482bab2f2d7cfb9171a65...0661c3b200d42c88acd94351f610be9bf52b934c
>---------------------------------------------------------------
commit 0661c3b200d42c88acd94351f610be9bf52b934c
Author: martinvandriel <martin at vandriel.de>
Date: Thu Aug 21 14:15:21 2014 +0200
computing time derivative of the stf and dump it to netcf
this is needed to deconvolve it in kerner/rdbm, because for a heaviside like stf,
computing the fft is useless. using sliprate circumvents this
>---------------------------------------------------------------
0661c3b200d42c88acd94351f610be9bf52b934c
SOLVER/UTILS/field_transform.F90 | 51 +++++++++++++++++++---------------
SOLVER/nc_routines.F90 | 60 ++++++++++++++++++++++++++++++++++++----
2 files changed, 83 insertions(+), 28 deletions(-)
diff --git a/SOLVER/UTILS/field_transform.F90 b/SOLVER/UTILS/field_transform.F90
index 3fd7ff1..d540ca0 100644
--- a/SOLVER/UTILS/field_transform.F90
+++ b/SOLVER/UTILS/field_transform.F90
@@ -44,7 +44,7 @@ program field_transformation
integer :: ncout_mesh_mpz_varid, ncin_mesh_mpz_varid
integer :: ncout_surf_grpid, ncout_surfstrain_dimid, ncout_surf_dimid
integer :: ncout_comp_dimid
- integer :: ncout_surf_varids(6), ncin_surf_grpid, ncin_surf_varids(6)
+ integer :: ncout_surf_varids(7), ncin_surf_grpid, ncin_surf_varids(7)
integer :: nc_mesh_npol_dimid, nc_mesh_cntrlpts_dimid, &
nc_mesh_elem_dimid
@@ -341,9 +341,10 @@ program field_transformation
len = nsurfelem, &
dimid = ncout_surf_dimid) )
- allocate(varname_surf(6))
+ allocate(varname_surf(7))
varname_surf = ['elem_theta ', 'displacement ', 'velocity ', &
- 'disp_src ', 'strain ', 'stf_dump ']
+ 'disp_src ', 'strain ', 'stf_dump ', &
+ 'stf_d_dump ']
do ivar = 1, size(varname_surf)
call check( nf90_inq_varid(ncid = ncin_surf_grpid, &
varid = ncin_surf_varids(ivar), &
@@ -385,16 +386,18 @@ program field_transformation
shuffle = 1, deflate = 1, &
deflate_level = deflate_level) )
- call check( nf90_def_var(ncid = ncout_surf_grpid, &
- name = varname_surf(6), &
- xtype = NF90_FLOAT, &
- dimids = [ncout_snap_dimid], &
- chunksizes = [nsnap], &
- varid = ncout_surf_varids(6)) )
- call check( nf90_def_var_deflate( ncid = ncout_surf_grpid, &
- varid = ncout_surf_varids(6), &
- shuffle = 1, deflate = 1, &
- deflate_level = deflate_level) )
+ do ivar = 6, 7
+ call check( nf90_def_var(ncid = ncout_surf_grpid, &
+ name = varname_surf(ivar), &
+ xtype = NF90_FLOAT, &
+ dimids = [ncout_snap_dimid], &
+ chunksizes = [nsnap], &
+ varid = ncout_surf_varids(ivar)) )
+ call check( nf90_def_var_deflate( ncid = ncout_surf_grpid, &
+ varid = ncout_surf_varids(ivar), &
+ shuffle = 1, deflate = 1, &
+ deflate_level = deflate_level) )
+ enddo
print *, 'Surface variables defined'
@@ -608,16 +611,18 @@ program field_transformation
deallocate(data_surf_3d)
allocate(data_surf_1d(nsnap))
- call check( nf90_get_var( ncid = ncin_surf_grpid, &
- varid = ncin_surf_varids(6), &
- start = [1], &
- count = [nsnap], &
- values = data_surf_1d) )
- call check( nf90_put_var( ncid = ncout_surf_grpid, &
- varid = ncout_surf_varids(6), &
- start = [1], &
- count = [nsnap], &
- values = data_surf_1d))
+ do ivar = 6, 7
+ call check( nf90_get_var( ncid = ncin_surf_grpid, &
+ varid = ncin_surf_varids(ivar), &
+ start = [1], &
+ count = [nsnap], &
+ values = data_surf_1d) )
+ call check( nf90_put_var( ncid = ncout_surf_grpid, &
+ varid = ncout_surf_varids(ivar), &
+ start = [1], &
+ count = [nsnap], &
+ values = data_surf_1d))
+ end do
deallocate(data_surf_1d)
diff --git a/SOLVER/nc_routines.F90 b/SOLVER/nc_routines.F90
index 4daddc2..d793774 100644
--- a/SOLVER/nc_routines.F90
+++ b/SOLVER/nc_routines.F90
@@ -87,13 +87,13 @@ module nc_routines
integer :: ncid_out, ncid_recout, ncid_snapout, ncid_surfout, ncid_meshout
integer :: nc_snap_dimid, nc_proc_dimid, nc_rec_dimid, nc_recproc_dimid
- integer :: nc_times_dimid, nc_comp_dimid, nc_disp_varid, nc_stf_seis_varid
- integer :: nc_time_varid, nc_iter_dimid, nc_stf_iter_varid
+ integer :: nc_times_dimid, nc_comp_dimid, nc_disp_varid, nc_stf_seis_varid, nc_stf_d_seis_varid
+ integer :: nc_time_varid, nc_iter_dimid, nc_stf_iter_varid, nc_stf_d_iter_varid
integer :: nc_strcomp_dimid
integer :: nc_surfelem_disp_varid, nc_surfelem_velo_varid
integer :: nc_surfelem_strain_varid, nc_surfelem_disp_src_varid
- integer :: nc_mesh_sol_varid, nc_mesh_flu_varid, nc_stf_dump_varid
+ integer :: nc_mesh_sol_varid, nc_mesh_flu_varid, nc_stf_dump_varid, nc_stf_d_dump_varid
integer :: nc_point_dimid, nc_pt_sol_dimid, nc_pt_flu_dimid
integer :: nc_szcoord_dimid
integer :: nc_snaptime_varid, nc_elem_dom_varid, nc_surfelem_theta_varid
@@ -122,6 +122,11 @@ module nc_routines
real(kind=sp), allocatable :: stf_seis_dumpvar(:)
real(kind=sp), allocatable :: stf_dumpvar(:)
+ !! Buffer variables for time derivative of the STF.
+ real(kind=sp), allocatable :: stf_d_dump_dumpvar(:)
+ real(kind=sp), allocatable :: stf_d_seis_dumpvar(:)
+ real(kind=sp), allocatable :: stf_d_dumpvar(:)
+
!> How many snaps should be buffered in RAM?
integer :: nc_dumpbuffersize
@@ -494,7 +499,7 @@ end subroutine nc_dump_strain_to_disk
!-----------------------------------------------------------------------------------------
subroutine nc_dump_stf(stf)
use data_io, only : nseismo, nstrain, dump_wavefields
- use data_time, only : seis_it, strain_it, niter
+ use data_time, only : seis_it, strain_it, niter, deltat
real(kind=sp), intent(in), dimension(:) :: stf
#ifdef unc
integer :: it_s, it_d, i
@@ -502,7 +507,27 @@ subroutine nc_dump_stf(stf)
allocate(stf_dumpvar(niter))
allocate(stf_seis_dumpvar(nseismo))
allocate(stf_dump_dumpvar(nstrain))
- stf_seis_dumpvar = 0.0
+
+ allocate(stf_d_dumpvar(niter))
+ allocate(stf_d_seis_dumpvar(nseismo))
+ allocate(stf_d_dump_dumpvar(nstrain))
+
+ stf_d_dumpvar = 0
+
+ stf_seis_dumpvar = 0
+ stf_dump_dumpvar = 0
+
+ stf_d_seis_dumpvar = 0
+ stf_d_dump_dumpvar = 0
+
+ ! compute derivative using central differences
+ do i = 2, niter-1
+ stf_d_dumpvar(i) = (stf(i+1) - stf(i-1)) / (2 * deltat)
+ enddo
+ ! compute derivative using fwd/bwd difference at the first and last sample
+ stf_d_dumpvar(1) = (stf(2) - stf(1)) / deltat
+ stf_d_dumpvar(niter) = (stf(niter) - stf(niter-1)) / deltat
+
it_s = 1
it_d = 1
@@ -510,6 +535,7 @@ subroutine nc_dump_stf(stf)
! Dumping the STF in the fine time stepping of the seismogram output
if ( mod(i,seis_it) == 0) then
stf_seis_dumpvar(it_s) = stf(i)
+ stf_d_seis_dumpvar(it_s) = stf_d_dumpvar(i)
it_s = it_s + 1
end if
@@ -517,6 +543,7 @@ subroutine nc_dump_stf(stf)
! Dumping the STF in the coarse time stepping of the strain (KERNER) output
if ( mod(i,strain_it) == 0) then
stf_dump_dumpvar(it_d) = stf(i)
+ stf_d_dump_dumpvar(it_s) = stf_d_dumpvar(i)
it_d = it_d + 1
end if
end if
@@ -1072,11 +1099,21 @@ subroutine nc_define_outputfile(nrec, rec_names, rec_th, rec_th_req, rec_ph, rec
deflate_level = deflate_level, &
varid=nc_stf_seis_varid) )
+ call check( nf90_def_var(ncid=ncid_recout, name="stf_d_seis", xtype=NF90_FLOAT,&
+ dimids=[nc_times_dimid], &
+ deflate_level = deflate_level, &
+ varid=nc_stf_d_seis_varid) )
+
call check( nf90_def_var(ncid=ncid_recout, name="stf_iter", xtype=NF90_FLOAT,&
dimids=[nc_iter_dimid], &
deflate_level = deflate_level, &
varid=nc_stf_iter_varid) )
+ call check( nf90_def_var(ncid=ncid_recout, name="stf_d_iter", xtype=NF90_FLOAT,&
+ dimids=[nc_iter_dimid], &
+ deflate_level = deflate_level, &
+ varid=nc_stf_d_iter_varid) )
+
call check( nf90_def_var(ncid=ncid_recout, name="time", xtype=NF90_DOUBLE,&
dimids=[nc_times_dimid], &
deflate_level = deflate_level, &
@@ -1311,6 +1348,10 @@ subroutine nc_define_outputfile(nrec, rec_names, rec_th, rec_th_req, rec_ph, rec
call check( nf90_def_var( ncid_surfout, "stf_dump", NF90_FLOAT, &
[nc_snap_dimid], &
nc_stf_dump_varid) )
+
+ call check( nf90_def_var( ncid_surfout, "stf_d_dump", NF90_FLOAT, &
+ [nc_snap_dimid], &
+ nc_stf_d_dump_varid) )
end if
@@ -1340,8 +1381,14 @@ subroutine nc_define_outputfile(nrec, rec_names, rec_th, rec_th_req, rec_ph, rec
varid = nc_stf_iter_varid, &
values = stf_dumpvar) )
call check( nf90_put_var(ncid = ncid_recout, &
+ varid = nc_stf_d_iter_varid, &
+ values = stf_d_dumpvar) )
+ call check( nf90_put_var(ncid = ncid_recout, &
varid = nc_stf_seis_varid, &
values = stf_seis_dumpvar) )
+ call check( nf90_put_var(ncid = ncid_recout, &
+ varid = nc_stf_d_seis_varid, &
+ values = stf_d_seis_dumpvar) )
if (dump_wavefields) then
! Write out strain dump times
@@ -1362,6 +1409,9 @@ subroutine nc_define_outputfile(nrec, rec_names, rec_th, rec_th_req, rec_ph, rec
call check( nf90_put_var(ncid = ncid_surfout, &
varid = nc_stf_dump_varid, &
values = stf_dump_dumpvar) )
+ call check( nf90_put_var(ncid = ncid_surfout, &
+ varid = nc_stf_d_dump_varid, &
+ values = stf_d_dump_dumpvar) )
if (verbose > 1) write(6,*) '...done'
end if
More information about the CIG-COMMITS
mailing list