Skip to content

Commit

Permalink
Merge branch 'master' into release_1.1
Browse files Browse the repository at this point in the history
Pull intermediate updates of master branch into release branch. In case that should ever be released, it might be useful to have the latest version.

Conflicts:
	MANUAL/manual_axisem1.1.tex
  • Loading branch information
sstaehler committed Aug 27, 2014
2 parents cbc40ea + c73403f commit edd15e0
Show file tree
Hide file tree
Showing 10 changed files with 298 additions and 422 deletions.
51 changes: 28 additions & 23 deletions SOLVER/UTILS/field_transform.F90
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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), &
Expand Down Expand Up @@ -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'

Expand Down Expand Up @@ -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)


Expand Down
2 changes: 1 addition & 1 deletion SOLVER/UTILS/post_processing.F90
Original file line number Diff line number Diff line change
Expand Up @@ -282,7 +282,7 @@ program post_processing_seis
! define time series
allocate(time(nt_seis), seis(nt_seis,3))
do iseis=1, nt_seis
time(iseis) = iseis * dt_seis
time(iseis) = (iseis - 1) * dt_seis
enddo
allocate(seis_sglcomp(nt_seis,3))

Expand Down
9 changes: 3 additions & 6 deletions SOLVER/data_io.f90
Original file line number Diff line number Diff line change
Expand Up @@ -34,19 +34,16 @@ module data_io
logical :: dump_snaps_glob
logical :: dump_vtk
logical :: dump_xdmf
logical :: dump_snaps_solflu
!> N.B. This is not wavefield snapshots, but kernel wavefields. Belongs to nstrain and
!! istrain
logical :: dump_wavefields
logical :: checkpointing
logical :: diagfiles !< Write diagnostic files (seismograms at antipodes,
!! list of surface elements, blabla), default: false

logical :: need_fluid_displ
real(kind=dp) :: strain_samp
integer :: iseismo !< current seismogram sample
integer :: istrain !< current kernel wavefield sample
integer :: isnap !< current wavefield sample for movies
!integer :: iseismo !< current seismogram sample
!integer :: istrain !< current kernel wavefield sample
!integer :: isnap !< current wavefield sample for movies
integer :: nseismo !< Number of seismogram samples
integer :: nstrain !< Number of wavefield dumps for kernels
integer :: nsnap !< Number of wavefield snapshots for movies
Expand Down
131 changes: 69 additions & 62 deletions SOLVER/meshes_io.F90
Original file line number Diff line number Diff line change
Expand Up @@ -40,8 +40,8 @@ module meshes_io
public :: dump_wavefields_mesh_1d
public :: dump_glob_grid_midpoint
public :: dump_xdmf_grid
public :: dump_solid_grid
public :: dump_fluid_grid
!public :: dump_solid_grid
!public :: dump_fluid_grid
public :: prepare_mesh_memoryvar_vtk
public :: build_kwf_grid
public :: dump_kwf_midpoint_xdmf
Expand Down Expand Up @@ -788,14 +788,21 @@ subroutine dump_kwf_grid()

points_mp = 0.

ct = 1
do iel=1, nel_solid
points_mp(iel,1) = scoord(npol/2,npol/2,ielsolid(iel))
points_mp(iel,2) = zcoord(npol/2,npol/2,ielsolid(iel))
if (kwf_mask(npol/2,jpol/2,iel)) then
points_mp(ct,1) = scoord(npol/2,npol/2,ielsolid(iel))
points_mp(ct,2) = zcoord(npol/2,npol/2,ielsolid(iel))
ct = ct + 1
endif
enddo

do iel=1, nel_fluid
points_mp(iel + nel_solid,1) = scoord(npol/2,npol/2,ielfluid(iel))
points_mp(iel + nel_solid,2) = zcoord(npol/2,npol/2,ielfluid(iel))
if (kwf_mask(npol/2,jpol/2,iel + nel_solid)) then
points_mp(ct,1) = scoord(npol/2,npol/2,ielfluid(iel))
points_mp(ct,2) = zcoord(npol/2,npol/2,ielfluid(iel))
ct = ct + 1
endif
enddo

if (use_netcdf) then
Expand Down Expand Up @@ -1012,25 +1019,25 @@ subroutine dump_kwf_sem_xdmf(filename, npoints, nelem)
!> Dumps the mesh (s,z) [m] in ASCII format as needed to visualize snapshots
!! in the solid region only.
!! Convention for order in the file: First the fluid, then the solid domain.
subroutine dump_solid_grid(ibeg,iend,jbeg,jend)


integer, intent(in) :: ibeg,iend,jbeg,jend
integer :: iel, ipol,jpol

open(unit=2500+mynum,file=datapath(1:lfdata)//'/solid_grid_'&
//appmynum//'.dat')
do iel=1,nel_solid
do jpol=jbeg,jend
do ipol=ibeg,iend
write(2500+mynum,*)scoord(ipol,jpol,ielsolid(iel)), &
zcoord(ipol,jpol,ielsolid(iel))
enddo
enddo
enddo
close(2500+mynum)

end subroutine dump_solid_grid
!subroutine dump_solid_grid(ibeg,iend,jbeg,jend)
!
!
! integer, intent(in) :: ibeg,iend,jbeg,jend
! integer :: iel, ipol,jpol
!
! open(unit=2500+mynum,file=datapath(1:lfdata)//'/solid_grid_'&
! //appmynum//'.dat')
! do iel=1,nel_solid
! do jpol=jbeg,jend
! do ipol=ibeg,iend
! write(2500+mynum,*)scoord(ipol,jpol,ielsolid(iel)), &
! zcoord(ipol,jpol,ielsolid(iel))
! enddo
! enddo
! enddo
! close(2500+mynum)
!
!end subroutine dump_solid_grid
!-----------------------------------------------------------------------------------------

!-----------------------------------------------------------------------------------------
Expand All @@ -1040,43 +1047,43 @@ end subroutine dump_solid_grid
!! When reading the fluid wavefield, one therefore needs to multiply all
!! components with inv_rho_fluid and the phi component with one/scoord!
!! Convention for order in the file: First the fluid, then the solid domain.
subroutine dump_fluid_grid(ibeg,iend,jbeg,jend)

use data_pointwise, only : inv_rho_fluid


integer, intent(in) :: ibeg,iend,jbeg,jend
integer :: iel, ipol,jpol

! When reading the fluid wavefield, one needs to multiply all components
! with inv_rho_fluid and the phi component with one/scoord!!

open(unit=2500+mynum,file=datapath(1:lfdata)//&
'/fluid_grid_'//appmynum//'.dat')
open(unit=2600+mynum,file=datapath(1:lfdata)//&
'/inv_rho_scoord_fluid_flusnaps_'&
//appmynum//'.dat', STATUS="REPLACE")
do iel=1,nel_fluid
do jpol=jbeg,jend
do ipol=ibeg,iend
write(2500+mynum,*)scoord(ipol,jpol,ielfluid(iel)), &
zcoord(ipol,jpol,ielfluid(iel))
if ( axis_fluid(iel) .and. ipol==0 ) then
! Axis s=0! write 1 instead of 1/s and then multiply
! with the correct factor dsdchi, obtained by L'Hospital's rule
! (see routine fluid_snapshot below).
write(2600+mynum,*)inv_rho_fluid(ipol,jpol,iel),one
else
write(2600+mynum,*)inv_rho_fluid(ipol,jpol,iel), &
one/scoord(ipol,jpol,ielfluid(iel))
endif
enddo
enddo
enddo
close(2500+mynum)
close(2600+mynum)

end subroutine dump_fluid_grid
!subroutine dump_fluid_grid(ibeg,iend,jbeg,jend)
!
! use data_pointwise, only : inv_rho_fluid
!
!
! integer, intent(in) :: ibeg,iend,jbeg,jend
! integer :: iel, ipol,jpol
!
! ! When reading the fluid wavefield, one needs to multiply all components
! ! with inv_rho_fluid and the phi component with one/scoord!!
!
! open(unit=2500+mynum,file=datapath(1:lfdata)//&
! '/fluid_grid_'//appmynum//'.dat')
! open(unit=2600+mynum,file=datapath(1:lfdata)//&
! '/inv_rho_scoord_fluid_flusnaps_'&
! //appmynum//'.dat', STATUS="REPLACE")
! do iel=1,nel_fluid
! do jpol=jbeg,jend
! do ipol=ibeg,iend
! write(2500+mynum,*)scoord(ipol,jpol,ielfluid(iel)), &
! zcoord(ipol,jpol,ielfluid(iel))
! if ( axis_fluid(iel) .and. ipol==0 ) then
! ! Axis s=0! write 1 instead of 1/s and then multiply
! ! with the correct factor dsdchi, obtained by L'Hospital's rule
! ! (see routine fluid_snapshot below).
! write(2600+mynum,*)inv_rho_fluid(ipol,jpol,iel),one
! else
! write(2600+mynum,*)inv_rho_fluid(ipol,jpol,iel), &
! one/scoord(ipol,jpol,ielfluid(iel))
! endif
! enddo
! enddo
! enddo
! close(2500+mynum)
! close(2600+mynum)
!
!end subroutine dump_fluid_grid
!-----------------------------------------------------------------------------------------

!-----------------------------------------------------------------------------------------
Expand Down
Loading

0 comments on commit edd15e0

Please sign in to comment.