[cig-commits] r13356 - in seismo/2D/SPECFEM2D/branches/BIOT: . DATA
cmorency at geodynamics.org
cmorency at geodynamics.org
Thu Nov 20 08:35:46 PST 2008
Author: cmorency
Date: 2008-11-20 08:35:46 -0800 (Thu, 20 Nov 2008)
New Revision: 13356
Modified:
seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file
seismo/2D/SPECFEM2D/branches/BIOT/DATA/interfaces_poro_flat.dat
seismo/2D/SPECFEM2D/branches/BIOT/meshfem2D.F90
seismo/2D/SPECFEM2D/branches/BIOT/part_unstruct.F90
seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90
Log:
Solved some compatibility bugs between multiple sources and irregular meshes.
Modified: seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file 2008-11-20 14:33:01 UTC (rev 13355)
+++ seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file 2008-11-20 16:35:46 UTC (rev 13356)
@@ -1,7 +1,7 @@
# title of job, and file that contains interface data
title = Test for 2 layers: acoustic/poroelastic
-interfacesfile = interfaces_acoustic.dat
+interfacesfile = interfaces_poro_flat.dat
# data concerning mesh, when generated using third-party app (more info in README)
read_external_mesh = .false.
@@ -28,7 +28,6 @@
TURN_ANISOTROPY_ON = .false. # turn anisotropy on or off for solid medium
TURN_ATTENUATION_ON = .false. # turn attenuation on or off for solid medium
-
# absorbing boundaries parameters
absorbing_conditions = .true. # absorbing boundary active or not
absorbbottom = .true.
@@ -41,13 +40,17 @@
deltat = 3d-4 # duration of a time step
isolver = 1 # type of simulation 1=forward 2=adjoint + kernels
-# source parameters
+#source parameters
+NSOURCE = 1 # number of sources
+
+#source 1
source_surf = .false. # source inside the medium or at the surface
xs = 1600. # source location x in meters
zs = 2900. # source location z in meters
source_type = 1 # elastic force or acoustic pressure = 1 or moment tensor = 2
time_function_type = 1 # Ricker = 1, first derivative = 2, Gaussian = 3, Dirac = 4, Heaviside = 5
f0 = 15.0 # dominant source frequency (Hz) if not Dirac or Heaviside
+t0 = 0.0 # time shift when multi sources (if one source, must be zero)
angleforce = 0. # angle of the source (for a force only)
Mxx = 1. # Mxx component (for a moment tensor source only)
Mzz = 1. # Mzz component (for a moment tensor source only)
@@ -100,13 +103,11 @@
outputgrid = .false. # save the grid in a text file or not
# velocity and density models
-nbmodels = 2 # nb of different models
+nbmodels = 1 # nb of different models
# define models as (model_number,1,rho_s,rho_f,phi,tort,permx,permz,kappa_s,kappa_f,kappa_fr,mu_s,eta_f,mu_fr) or (Anisotropic: to be defined)
# set the porosity phi to 1 to make a given model acoustic, and to 0 to make it elastic
# the mesh can contain acoustic, elastic and poroelastic models simultaneously
-1 1 2500.d0 1020.d0 0.4d0 2.0 1d-11 0.d0 1d-11 1.60554d10 2.295d9 1.0d10 9.63342d9 0.0d-4 9.63342d9
-2 1 2650.d0 1020.d0 1.10d0 2.0 1d-11 0.d0 1d-11 1.60554d10 2.295d9 0.0d9 0.0d9 0.0d-3 0.0d9
+1 1 2500.d0 1020.d0 0.0d0 2.0 1d-11 0.d0 1d-11 1.60554d10 2.295d9 1.0d10 9.63342d9 0.0d-4 9.63342d9
# define the different regions of the model in the (nx,nz) spectral element mesh
-nbregions = 2 # nb of regions and model number for each
-1 260 1 110 1
-1 260 111 220 2
+nbregions = 1 # nb of regions and model number for each
+1 260 1 160 1
Modified: seismo/2D/SPECFEM2D/branches/BIOT/DATA/interfaces_poro_flat.dat
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/DATA/interfaces_poro_flat.dat 2008-11-20 14:33:01 UTC (rev 13355)
+++ seismo/2D/SPECFEM2D/branches/BIOT/DATA/interfaces_poro_flat.dat 2008-11-20 16:35:46 UTC (rev 13356)
@@ -10,13 +10,13 @@
#
2
0 0
- 1000 0
+ 4800 0
#
# interface number 4 (topography, top of the mesh)
#
2
- 0 1000
- 1000 1000
+ 0 4800
+ 4800 4800
#
# for each layer, we give the number of spectral elements in the vertical direction
#
Modified: seismo/2D/SPECFEM2D/branches/BIOT/meshfem2D.F90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/meshfem2D.F90 2008-11-20 14:33:01 UTC (rev 13355)
+++ seismo/2D/SPECFEM2D/branches/BIOT/meshfem2D.F90 2008-11-20 16:35:46 UTC (rev 13356)
@@ -430,7 +430,6 @@
call read_value_integer(IIN,IGNORE_JUNK,isolver)
! read source parameters
-!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
call read_value_integer(IIN,IGNORE_JUNK,NSOURCE) !yang
allocate(source_surf(NSOURCE))
allocate(xs(NSOURCE))
@@ -444,6 +443,7 @@
allocate(Mxz(NSOURCE))
allocate(Mzz(NSOURCE))
allocate(factor(NSOURCE))
+
do i_source=1,NSOURCE
call read_value_logical(IIN,IGNORE_JUNK,source_surf(i_source))
call read_value_double_precision(IIN,IGNORE_JUNK,xs(i_source))
@@ -480,6 +480,7 @@
print *,'Mxz of the source if moment tensor = ',Mxz(i_source)
print *,'Multiplying factor = ',factor(i_source)
enddo
+
! read constants for attenuation
call read_value_integer(IIN,IGNORE_JUNK,N_SLS)
call read_value_double_precision(IIN,IGNORE_JUNK,Qp_attenuation)
@@ -1328,7 +1329,9 @@
enddo
!--- compute position of the receivers and write the STATIONS file
if (read_external_mesh) then
- call read_receivers(receivers_file,xs,zs)
+
+ call read_receivers(receivers_file,xs,zs,NSOURCE)
+
else
if (generate_STATIONS) then
print *
Modified: seismo/2D/SPECFEM2D/branches/BIOT/part_unstruct.F90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/part_unstruct.F90 2008-11-20 14:33:01 UTC (rev 13355)
+++ seismo/2D/SPECFEM2D/branches/BIOT/part_unstruct.F90 2008-11-20 16:35:46 UTC (rev 13356)
@@ -249,12 +249,13 @@
! first line: number of receivers
! following lines: xrec zrec
!-----------------------------------------------
- subroutine read_receivers(filename,xs,zs)
+ subroutine read_receivers(filename,xs,zs,NSOURCE)
character(len=256), intent(in) :: filename
- integer :: nrec, irec_global_number
- double precision :: xrec,zrec,xs,zs
- integer :: i
+ integer :: nrec, irec_global_number,NSOURCE
+ double precision :: xrec,zrec
+ double precision :: xs(NSOURCE),zs(NSOURCE)
+ integer :: i,i_source
open(unit=996,file='DATA/STATIONS',status='unknown')
open(unit=997,file='OUTPUT_FILES/receivers_file',status='unknown')
@@ -270,8 +271,9 @@
print *,'Position (x,z) of the ',nrec,' receivers'
print *
-
- write(997,"(f20.7,1x,f20.7)") xs,zs
+ do i_source = 1,NSOURCE
+ write(997,"(f20.7,1x,f20.7)") xs(i_source),zs(i_source)
+ enddo
do i = 1, nrec
read(995,*) xrec,zrec
write(996,"('S',i4.4,' AA ',f20.7,1x,f20.7,' 0.0 0.0')") i,xrec,zrec
Modified: seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90 2008-11-20 14:33:01 UTC (rev 13355)
+++ seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90 2008-11-20 16:35:46 UTC (rev 13356)
@@ -195,7 +195,6 @@
real(kind=CUSTOM_REAL), dimension(:,:), allocatable :: accel_elastic,veloc_elastic,displ_elastic
double precision, dimension(:,:), allocatable :: coord, flagrange,xinterp,zinterp,Uxinterp,Uzinterp,vector_field_display
- double precision, dimension(:,:), allocatable :: vector_field_display1
! material properties of the poroelastic medium (solid phase:s and fluid phase [defined as w=phi(u_f-u_s)]: w)
real(kind=CUSTOM_REAL), dimension(:,:), allocatable :: accels_poroelastic,velocs_poroelastic,displs_poroelastic
@@ -359,7 +358,7 @@
double precision :: xmin_color_image,xmax_color_image, &
zmin_color_image,zmax_color_image,size_pixel_horizontal,size_pixel_vertical
integer, dimension(:,:), allocatable :: iglob_image_color,copy_iglob_image_color
- double precision, dimension(:,:), allocatable :: image_color_data, image_color_data1
+ double precision, dimension(:,:), allocatable :: image_color_data
double precision, dimension(:,:), allocatable :: image_color_vp_display
double precision :: xmin_color_image_loc, xmax_color_image_loc, zmin_color_image_loc, &
@@ -1127,7 +1126,6 @@
! to display acoustic elements
allocate(vector_field_display(NDIM,npoin))
- allocate(vector_field_display1(NDIM,npoin))
if(assign_external_model) then
allocate(vpext(NGLLX,NGLLZ,nspec))
@@ -1856,7 +1854,6 @@
! allocate an array for image data
allocate(image_color_data(NX_IMAGE_color,NZ_IMAGE_color))
- allocate(image_color_data1(NX_IMAGE_color,NZ_IMAGE_color))
allocate(image_color_vp_display(NX_IMAGE_color,NZ_IMAGE_color))
! allocate an array for the grid point that corresponds to a given image data point
@@ -2659,7 +2656,7 @@
write(IOUT,*)
open(unit=55,file='OUTPUT_FILES/source.txt',status='unknown')
endif
-do i_source=1,NSOURCE !yang
+ do i_source=1,NSOURCE !yang
! loop on all the time steps
do it = 1,NSTEP
@@ -2705,8 +2702,7 @@
write(55,*) sngl(time),real(source_time_function(i_source,it),4),sngl(time-t0(i_source))
endif
enddo
- write(*,*) '!!!!saving source time function in a text file has been omitted!!!!!',i_source
-enddo ! i_source=1,NSOURCE !yang
+ enddo ! i_source=1,NSOURCE !yang
if ( myrank == 0 ) then
close(55)
endif
@@ -5381,13 +5377,10 @@
if ( myrank == 0 ) then
write(IOUT,*) 'drawing image of vertical component of displacement vector...'
endif
-!!!!!!!!!!!!!!!!!!!!!!!!yang!!!!!!!!!!!!!!!!!!!!!!!
+
call compute_vector_whole_medium(potential_acoustic,displ_elastic,displs_poroelastic,&
elastic,poroelastic,vector_field_display, &
xix,xiz,gammax,gammaz,ibool,hprime_xx,hprime_zz,nspec,npoin)
- call compute_vector_whole_medium(potential_acoustic,displ_elastic,displw_poroelastic,&
- elastic,poroelastic,vector_field_display1, &
- xix,xiz,gammax,gammaz,ibool,hprime_xx,hprime_zz,nspec,npoin)
else if(imagetype == 2) then
@@ -5431,7 +5424,6 @@
j = ceiling(real(num_pixel_loc(k)) / real(NX_IMAGE_color))
i = num_pixel_loc(k) - (j-1)*NX_IMAGE_color
image_color_data(i,j) = vector_field_display(2,iglob_image_color(i,j))
- image_color_data1(i,j) = vector_field_display1(2,iglob_image_color(i,j))
end do
! assembling array image_color_data on process zero for color output
@@ -5470,7 +5462,6 @@
if ( myrank == 0 ) then
call create_color_image(image_color_data,iglob_image_color,NX_IMAGE_color,NZ_IMAGE_color,it,cutsnaps,image_color_vp_display)
- call create_color_image1(image_color_data1,iglob_image_color,NX_IMAGE_color,NZ_IMAGE_color,it,cutsnaps,image_color_vp_display)
write(IOUT,*) 'Color image created'
endif
More information about the CIG-COMMITS
mailing list