[cig-commits] r15499 - in seismo/2D/SPECFEM2D/branches/BIOT: . DATA
cmorency at geodynamics.org
cmorency at geodynamics.org
Sat Aug 1 10:51:09 PDT 2009
Author: cmorency
Date: 2009-08-01 10:51:08 -0700 (Sat, 01 Aug 2009)
New Revision: 15499
Modified:
seismo/2D/SPECFEM2D/branches/BIOT/DATA/CMTSOLUTION
seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file
seismo/2D/SPECFEM2D/branches/BIOT/DATA/STATIONS
seismo/2D/SPECFEM2D/branches/BIOT/assemble_MPI.F90
seismo/2D/SPECFEM2D/branches/BIOT/checkgrid.F90
seismo/2D/SPECFEM2D/branches/BIOT/define_external_model.f90
seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90
Log:
MPI ok for forward (coupled or not) calculations + define_external_model.f90 updated (myrank)
Modified: seismo/2D/SPECFEM2D/branches/BIOT/DATA/CMTSOLUTION
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/DATA/CMTSOLUTION 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/DATA/CMTSOLUTION 2009-08-01 17:51:08 UTC (rev 15499)
@@ -1,11 +1,11 @@
# source 1
source_surf = .false. # source inside the medium or at the surface
-xs = 2500. # source location x in meters
-zs = 2500. # source location z in meters
-source_type = 1 # elastic force or acoustic pressure = 1 or moment tensor = 2
+xs = 1600. # source location x in meters
+zs = 2900. # source location z in meters
+source_type = 2 # 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 = 3 # dominant source frequency (Hz) if not Dirac or Heaviside
-t0 = 0. # offset of the source, irrelevant if NSOURCE=1
+f0 = 15 # dominant source frequency (Hz) if not Dirac or Heaviside
+t0 = 0.
angleforce = 00. # 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)
Modified: seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/DATA/Par_file 2009-08-01 17:51:08 UTC (rev 15499)
@@ -1,7 +1,7 @@
# title of job, and file that contains interface data
title = Test for M2 UPPA
-interfacesfile = interfaces_1layer.dat
+interfacesfile = interfaces_reg2layers.dat
# data concerning mesh, when generated using third-party app (more info in README)
read_external_mesh = .false.
@@ -19,8 +19,8 @@
# geometry of the model (origin lower-left corner = 0,0) and mesh description
xmin = 0.d0 # abscissa of left side of the model
-xmax = 10.d3 # abscissa of right side of the model
-nx = 100 # number of elements along X
+xmax = 4800.d0 # abscissa of right side of the model
+nx = 200 # number of elements along X
ngnod = 9 # number of control nodes per element (4 or 9)
initialfield = .false. # use a plane wave as source or not
add_Bielak_conditions = .false. # add Bielak conditions or not if initial plane wave
@@ -39,9 +39,9 @@
absorbleft = .true.
# time step parameters
-nt = 6000 # total number of time steps
-deltat = 1d-3 # duration of a time step
-isolver = 2 # type of simulation 1=forward 2=adjoint + kernels
+nt = 2600 # total number of time steps
+deltat = 0.5d-3 # duration of a time step
+isolver = 1 # type of simulation 1=forward 2=adjoint + kernels
# source parameters
NSOURCE = 1 # number of sources [source info read in CMTSOLUTION file]
@@ -56,24 +56,32 @@
f0_attenuation = 5.196152422706633 # (Hz) relevant only if source is a Dirac or a Heaviside, else it is f0
# receiver line parameters for seismograms
-seismotype = 1 # record 1=displ 2=veloc 3=accel 4=pressure 5=curl 6=potential
+seismotype = 2 # record 1=displ 2=veloc 3=accel 4=pressure 5=curl 6=potential
save_forward = .false. # save the last frame, needed for adjoint simulation
generate_STATIONS = .true. # creates a STATION file in ./DATA
-nreceiverlines = 1 # number of receiver lines
+nreceiverlines = 2 # number of receiver lines
anglerec = 0.d0 # angle to rotate components at receivers
rec_normal_to_surface = .false. # base anglerec normal to surface (external mesh and curve file needed)
# first receiver line
nrec = 1 # number of receivers
-xdeb = 7500. # first receiver x in meters
-zdeb = 2500. # first receiver z in meters
+xdeb = 2000. # first receiver x in meters
+zdeb = 2933.33 # first receiver z in meters
xfin = 3777. # last receiver x in meters (ignored if onlyone receiver)
zfin = 1866.67 # last receiver z in meters (ignored if onlyone receiver)
enreg_surf = .false. # receivers inside the medium or at the surface
+# second receiver line
+nrec = 1 # number of receivers
+xdeb = 2000. # first receiver x in meters
+zdeb = 1866.67 # first receiver z in meters
+xfin = 3777. # last receiver x in meters (ignored if onlyone receiver)
+zfin = 1866.67 # last receiver z in meters (ignored if onlyone receiver)
+enreg_surf = .false. # receivers inside the medium or at the surface
+
# display parameters
NTSTEP_BETWEEN_OUTPUT_INFO = 200 # display frequency in time steps
-output_postscript_snapshot = .false. # output Postscript snapshot of the results
+output_postscript_snapshot = .true. # output Postscript snapshot of the results
output_color_image = .true. # output color image of the results
imagetype = 1 # display 1=displ 2=veloc 3=accel 4=pressure
cutsnaps = 1. # minimum amplitude in % for snapshots
@@ -88,11 +96,13 @@
outputgrid = .false. # save the grid in a text file or not
# velocity and density models
-nbmodels = 1 # nb of different models
+nbmodels = 2 # nb of different models
# define models as I: (model_number,1,rho,Vp,Vs,0,0,Qp,Qs) or II: (model_number,2,rho,c11,c13,c33,c44,Qp,Qs) or III: (model_number,3,rhos,rhof,phi,c,kxx,kxz,kzz,Ks,Kf,Kfr,etaf,mufr,Qs).
# For istropic elastic/acoustic material use I and set Vs to zero to make a given model acoustic, for anisotropic elastic use II,
# and for isotropic poroelastic material use III. The mesh can contain acoustic, elastic, & poroelastic models simultaneously
-1 1 2500.d0 3000.d0 1800.d0 0 0 10.d0 10.d0 0 0 0 0 0 0 #1558.89d0 0 0 136.d0 136.d0 0 0 0 0 0 0
+1 3 2200.d0 1040.d0 0.4d0 2.0 1d-11 0.d0 1d-11 6.9d9 4.0d9 6.7d9 0.0d-3 3.d9 10.d0 # 2700.d0 3000.d0 1732.051d0 0 0 10.d0 10.d0 0 0 0 0 0 0
+2 1 2500.d0 3000.d0 2000.d0 0 0 10.d0 10.d0 0 0 0 0 0 0 #1558.89d0 0 0 136.d0 136.d0 0 0 0 0 0 0
# define the different regions of the model in the (nx,nz) spectral element mesh
-nbregions = 1 # nb of regions and model number for each
-1 100 1 80 1
+nbregions = 2 # nb of regions and model number for each
+1 200 1 70 2
+1 200 71 140 1
Modified: seismo/2D/SPECFEM2D/branches/BIOT/DATA/STATIONS
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/DATA/STATIONS 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/DATA/STATIONS 2009-08-01 17:51:08 UTC (rev 15499)
@@ -1 +1,2 @@
-S0001 AA 7500.0000000 2500.0000000 0.0 0.0
+S0001 AA 2000.0000000 2933.3300000 0.0 0.0
+S0002 AA 2000.0000000 1866.6700000 0.0 0.0
Modified: seismo/2D/SPECFEM2D/branches/BIOT/assemble_MPI.F90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/assemble_MPI.F90 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/assemble_MPI.F90 2009-08-01 17:51:08 UTC (rev 15499)
@@ -340,11 +340,11 @@
integer, intent(in) :: max_ibool_interfaces_size_ac,max_ibool_interfaces_size_el,max_ibool_interfaces_size_po
integer, dimension(NGLLX*max_interface_size,ninterface), intent(in) :: &
ibool_interfaces_acoustic,ibool_interfaces_elastic,ibool_interfaces_poroelastic
- integer, dimension(ninterface), intent(in) :: nibool_interfaces_acoustic,nibool_interfaces_elastic &
+ integer, dimension(ninterface), intent(in) :: nibool_interfaces_acoustic,nibool_interfaces_elastic, &
nibool_interfaces_poroelastic
integer, dimension(ninterface), intent(in) :: my_neighbours
- double precision, dimension(max_ibool_interfaces_size_ac+max_ibool_interfaces_size_el, ninterface) :: &
+ double precision, dimension(max_ibool_interfaces_size_ac+max_ibool_interfaces_size_el+2*max_ibool_interfaces_size_po, ninterface) :: &
buffer_send_faces_scalar, &
buffer_recv_faces_scalar
integer :: msg_status(MPI_STATUS_SIZE)
@@ -687,7 +687,7 @@
integer, intent(in) :: max_ibool_interfaces_size_po
integer, dimension(NGLLX*max_interface_size,ninterface), intent(in) :: ibool_interfaces_poroelastic
integer, dimension(ninterface), intent(in) :: nibool_interfaces_poroelastic
- integer, dimension(ninterface_elastic*2), intent(inout) :: tab_requests_send_recv_poroelastic
+ integer, dimension(ninterface_poroelastic*4), intent(inout) :: tab_requests_send_recv_poroelastic
real(CUSTOM_REAL), dimension(max_ibool_interfaces_size_po,ninterface_poroelastic), intent(inout) :: &
buffer_send_faces_vector_pos,buffer_send_faces_vector_pow
real(CUSTOM_REAL), dimension(max_ibool_interfaces_size_po,ninterface_poroelastic), intent(inout) :: &
@@ -765,7 +765,7 @@
do inum_interface = 1, ninterface_poroelastic*4
- call MPI_Wait (tab_requests_send_recv_elastic(inum_interface), status_poroelastic, ier)
+ call MPI_Wait (tab_requests_send_recv_poroelastic(inum_interface), status_poroelastic, ier)
enddo
Modified: seismo/2D/SPECFEM2D/branches/BIOT/checkgrid.F90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/checkgrid.F90 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/checkgrid.F90 2009-08-01 17:51:08 UTC (rev 15499)
@@ -107,7 +107,6 @@
#ifdef USE_MPI
double precision :: vpImin_glob,vpImax_glob,vsmin_glob,vsmax_glob,densmin_glob,densmax_glob
double precision :: vpIImin_glob,vpIImax_glob
- double precision :: densmin_glob,densmax_glob
double precision :: distance_min_glob,distance_max_glob
double precision :: courant_stability_max_glob,lambdaPImin_glob,lambdaPImax_glob,&
lambdaPIImin_glob,lambdaPIImax_glob,lambdaSmin_glob,lambdaSmax_glob
Modified: seismo/2D/SPECFEM2D/branches/BIOT/define_external_model.f90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/define_external_model.f90 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/define_external_model.f90 2009-08-01 17:51:08 UTC (rev 15499)
@@ -40,7 +40,7 @@
!
!========================================================================
- subroutine define_external_model(x,y,iflag_element,rho,vp,vs)
+ subroutine define_external_model(x,y,iflag_element,myrank,rho,vp,vs)
implicit none
@@ -49,14 +49,14 @@
! user can modify this routine to assign any different external Earth model (rho, vp, vs)
! based on the x and y coordinates of that grid point and the flag of the region it belongs to
- integer, intent(in) :: iflag_element
+ integer, intent(in) :: iflag_element,myrank
double precision, intent(in) :: x,y
double precision, intent(out) :: rho,vp,vs
! dummy routine here, just to demonstrate how the model can be assigned
- if(iflag_element == 1 .or. x < 1700.d0 .or. y >= 2300.d0) then
+ if(myrank == 0 .and. iflag_element == 1 .or. x < 1700.d0 .or. y >= 2300.d0) then
rho = 2000.d0
vp = 3000.d0
vs = vp / sqrt(3.d0)
Modified: seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90
===================================================================
--- seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90 2009-08-01 11:06:19 UTC (rev 15498)
+++ seismo/2D/SPECFEM2D/branches/BIOT/specfem2D.F90 2009-08-01 17:51:08 UTC (rev 15499)
@@ -1529,7 +1529,7 @@
do j = 1,NGLLZ
do i = 1,NGLLX
iglob = ibool(i,j,ispec)
- call define_external_model(coord(1,iglob),coord(2,iglob),kmato(ispec), &
+ call define_external_model(coord(1,iglob),coord(2,iglob),kmato(ispec),myrank, &
rhoext(i,j,ispec),vpext(i,j,ispec),vsext(i,j,ispec),myrank)
! stop if the same element is assigned both acoustic and elastic points in external model
if(.not. (i == 1 .and. j == 1) .and. &
@@ -1654,7 +1654,7 @@
enddo ! do i_source=1,NSOURCE
! compute source array for adjoint source
- if(isolver == 2 .and. ipass == 1) then ! adjoint calculation
+ if(isolver == 2) then ! adjoint calculation
nadj_rec_local = 0
do irec = 1,nrec
if(myrank == which_proc_receiver(irec))then
@@ -1664,8 +1664,8 @@
nadj_rec_local = nadj_rec_local + 1
endif
enddo
- allocate(adj_sourcearray(NSTEP,NDIM,NGLLX,NGLLZ))
- if (nadj_rec_local > 0) allocate(adj_sourcearrays(nadj_rec_local,NSTEP,NDIM,NGLLX,NGLLZ))
+ if(ipass == 1) allocate(adj_sourcearray(NSTEP,NDIM,NGLLX,NGLLZ))
+ if (nadj_rec_local > 0 .and. ipass == 1) allocate(adj_sourcearrays(nadj_rec_local,NSTEP,NDIM,NGLLX,NGLLZ))
irec_local = 0
do irec = 1, nrec
! compute only adjoint source arrays in the local proc
@@ -4005,6 +4005,24 @@
endif
+! initiation
+ if(any_poroelastic .and. anyabs) then
+! loop on all the absorbing elements
+ do ispecabs = 1,nelemabs
+ jbegin_left_poro(ispecabs) = 1
+ jbegin_right_poro(ispecabs) = 1
+
+ jend_left_poro(ispecabs) = NGLLZ
+ jend_right_poro(ispecabs) = NGLLZ
+
+ ibegin_bottom_poro(ispecabs) = 1
+ ibegin_top_poro(ispecabs) = 1
+
+ iend_bottom_poro(ispecabs) = NGLLX
+ iend_top_poro(ispecabs) = NGLLX
+ enddo
+ endif
+
! exclude common points between poroelastic absorbing edges and elastic/poroelastic matching interfaces
if(coupled_elastic_poroelastic .and. anyabs) then
@@ -4759,6 +4777,15 @@
tab_requests_send_recv_acoustic,buffer_send_faces_vector_ac, &
buffer_recv_faces_vector_ac, my_neighbours)
endif
+
+ if ( nproc > 1 .and. any_acoustic .and. ninterface_acoustic > 0 .and. isolver == 2) then
+ call assemble_MPI_vector_ac(b_potential_dot_dot_acoustic,npoin, &
+ ninterface, ninterface_acoustic,inum_interfaces_acoustic, &
+ max_interface_size, max_ibool_interfaces_size_ac,&
+ ibool_interfaces_acoustic, nibool_interfaces_acoustic, &
+ tab_requests_send_recv_acoustic,buffer_send_faces_vector_ac, &
+ buffer_recv_faces_vector_ac, my_neighbours)
+ endif
#endif
@@ -5289,6 +5316,15 @@
tab_requests_send_recv_elastic,buffer_send_faces_vector_el, &
buffer_recv_faces_vector_el, my_neighbours)
endif
+
+ if (nproc > 1 .and. any_elastic .and. ninterface_elastic > 0 .and. isolver == 2) then
+ call assemble_MPI_vector_el(b_accel_elastic,npoin, &
+ ninterface, ninterface_elastic,inum_interfaces_elastic, &
+ max_interface_size, max_ibool_interfaces_size_el,&
+ ibool_interfaces_elastic, nibool_interfaces_elastic, &
+ tab_requests_send_recv_elastic,buffer_send_faces_vector_el, &
+ buffer_recv_faces_vector_el, my_neighbours)
+ endif
#endif
@@ -5764,6 +5800,16 @@
buffer_recv_faces_vector_pos,buffer_recv_faces_vector_pow, &
my_neighbours)
endif
+
+ if (nproc > 1 .and. any_poroelastic .and. ninterface_poroelastic > 0 .and. isolver == 2) then
+ call assemble_MPI_vector_po(b_accels_poroelastic,b_accelw_poroelastic,npoin, &
+ ninterface, ninterface_poroelastic,inum_interfaces_poroelastic, &
+ max_interface_size, max_ibool_interfaces_size_po,&
+ ibool_interfaces_poroelastic, nibool_interfaces_poroelastic, &
+ tab_requests_send_recv_poroelastic,buffer_send_faces_vector_pos,buffer_send_faces_vector_pow, &
+ buffer_recv_faces_vector_pos,buffer_recv_faces_vector_pow, &
+ my_neighbours)
+ endif
#endif
More information about the CIG-COMMITS
mailing list