[cig-commits] r18119 - in seismo/3D/SPECFEM3D_GLOBE/trunk: setup src/meshfem3D src/specfem3D
danielpeter at geodynamics.org
danielpeter at geodynamics.org
Wed Mar 16 14:39:41 PDT 2011
Author: danielpeter
Date: 2011-03-16 14:39:41 -0700 (Wed, 16 Mar 2011)
New Revision: 18119
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_gll.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
Log:
fixes bug in array allocations; puts GLL model constants to constants.h
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in 2011-03-16 21:18:25 UTC (rev 18118)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in 2011-03-16 21:39:41 UTC (rev 18119)
@@ -82,8 +82,14 @@
! average density in the full Earth to normalize equation
double precision, parameter :: RHOAV = 5514.3d0
-! for topography/bathymetry model
+!!-----------------------------------------------------------
+!!
+!! for topography/bathymetry model
+!!
+!!-----------------------------------------------------------
+!! (uncomment desired resolution)
+
!!--- ETOPO5 5-minute model, smoothed Harvard version
!! size of topography and bathymetry file
! integer, parameter :: NX_BATHY = 4320,NY_BATHY = 2160
@@ -129,7 +135,15 @@
! minimum thickness in meters to include the effect of the oceans and topo
double precision, parameter :: MINIMUM_THICKNESS_3D_OCEANS = 100.d0
-!-- crustal models
+
+!!-----------------------------------------------------------
+!!
+!! for crustal model
+!!
+!!-----------------------------------------------------------
+!! (uncomment desired model)
+
+! crustal model constants
integer, parameter :: ICRUST_CRUST2 = 1
integer, parameter :: ICRUST_CRUSTMAPS = 2
integer, parameter :: ICRUST_EPCRUST = 3
@@ -141,18 +155,32 @@
logical, parameter :: INCLUDE_SEDIMENTS_CRUST = .true.
double precision, parameter :: MINIMUM_SEDIMENT_THICKNESS = 2.d0 ! minimim thickness in km
-!-- uncomment for using Crust2.0 (used when CRUSTAL flag is set for simulation)
+!!-- uncomment for using Crust2.0 (used as default when CRUSTAL flag is set for simulation)
integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUST2
!!-- uncomment for using General Crustmaps instead
-!! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUSTMAPS
-!!-- uncomment for using EPcrust instead
-!! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_EPCRUST
+! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUSTMAPS
+!!-- uncomment for using EPcrust instead (European regional model)
+! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_EPCRUST
+
+!!-----------------------------------------------------------
+!!
+!! Gauss-Lobatto-Legrendre resolution
+!!
+!!-----------------------------------------------------------
+
! number of GLL points in each direction of an element (degree plus one)
integer, parameter :: NGLLX = 5
integer, parameter :: NGLLY = NGLLX
integer, parameter :: NGLLZ = NGLLX
+
+!!-----------------------------------------------------------
+!!
+!! source/receiver setup
+!!
+!!-----------------------------------------------------------
+
! flag to exclude elements that are too far from target in source detection
logical, parameter :: USE_DISTANCE_CRITERION = .true.
@@ -188,6 +216,13 @@
! is located outside the mesh and therefore excluded from the station list
double precision, parameter :: THRESHOLD_EXCLUDE_STATION = 50.d0
+
+!!-----------------------------------------------------------
+!!
+!! mesh optimization
+!!
+!!-----------------------------------------------------------
+
! the first doubling is implemented right below the Moho
! it seems optimal to implement the three other doublings at these depths
! in the mantle
@@ -201,6 +236,13 @@
! the mesher) and use them for the computation of boundary kernel (in the solver)
logical, parameter :: SAVE_BOUNDARY_MESH = .false.
+
+!!-----------------------------------------------------------
+!!
+!! adjoint kernel outputs
+!!
+!!-----------------------------------------------------------
+
! this parameter must be set to .true. to compute anisotropic kernels
! in crust and mantle (related to the 21 Cij in geographical coordinates)
! default is .false. to compute isotropic kernels (related to alpha and beta)
@@ -216,6 +258,13 @@
! output kernel mask to zero out source region
logical,parameter :: SAVE_SOURCE_MASK = .false.
+
+!!-----------------------------------------------------------
+!!
+!! time stamp information
+!!
+!!-----------------------------------------------------------
+
! print date and time estimate of end of run in another country,
! in addition to local time.
! For instance: the code runs at Caltech in California but the person
@@ -225,10 +274,13 @@
integer, parameter :: HOURS_TIME_DIFFERENCE = +9
integer, parameter :: MINUTES_TIME_DIFFERENCE = +0
-!
-!--- debugging flags
-!
+!!-----------------------------------------------------------
+!!
+!! debugging flags
+!!
+!!-----------------------------------------------------------
+
! flags to actually assemble with MPI or not
! and to actually match fluid and solid regions of the Earth or not
! should always be set to true except when debugging code
@@ -242,6 +294,7 @@
! can be useful for benchmarks of a spherical Earth with fictitious sources and stations
logical, parameter :: ASSUME_PERFECT_SPHERE = .false.
+
!------------------------------------------------------
!----------- do not modify anything below -------------
!------------------------------------------------------
@@ -480,6 +533,11 @@
integer, parameter :: NRAD_GRAVITY = 70000
!!!!!!!!!!!!!! parameters added for the thread-safe version of the code
+
+!!
+!! model constants
+!!
+
! number of layers in DATA/1066a/1066a.dat
integer, parameter :: NR_1066A = 160
@@ -502,7 +560,7 @@
integer, parameter :: MPA=42,MRA=48,MHA=21,MPB=42,MRB=48,MHB=18
integer, parameter :: MKA=2101,MKB=2101
-!QRFSI12 constants
+! QRFSI12 constants
integer,parameter :: NKQ=8,MAXL_Q=12
integer,parameter :: NSQ=(MAXL_Q+1)**2,NDEPTHS_REFQ=913
@@ -541,8 +599,19 @@
integer, parameter :: NTHETA_EP = 4, NPHI_EP = 20
double precision, parameter :: cap_degree_EP = 0.2d0
+! parameters for GLL model (used for iterative model inversions)
+ character(len=*), parameter :: PATHNAME_GLL_modeldir = 'DATA/GLL/'
+ integer, parameter :: GLL_REFERENCE_MODEL = THREE_D_MODEL_S29EA
+!!-- uncomment for using S362ANI as reference model (Hejun Zhu)
+! integer, parameter :: GLL_REFERENCE_MODEL = THREE_D_MODEL_S362ANI
+
+
!!!!!!!!!!!!!! end of parameters added for the thread-safe version of the code
+!!
+!! crustal stretching
+!!
+
! for the stretching of crustal elements in the case of 3D models
! (values are chosen for 3D models to have RMOHO_FICTICIOUS at 35 km
! and RMIDDLE_CRUST to become 15 km with stretching function stretch_tab)
@@ -556,11 +625,14 @@
logical, parameter :: REGIONAL_MOHO_MESH_EUROPE = .false. ! used only for fixing time step
logical, parameter :: REGIONAL_MOHO_MESH_ASIA = .false. ! used only for fixing time step
logical, parameter :: HONOR_DEEP_MOHO = .false.
-! uncomment for e.g. Europe case, where deep moho is rare
-!! double precision, parameter :: RMOHO_STRETCH_ADJUSTEMENT = -15000.d0 ! moho mesh boundary down to 55km
-! uncomment for deep moho cases, e.g. Asia case (Himalayan moho)
-!! double precision, parameter :: RMOHO_STRETCH_ADJUSTEMENT = -20000.d0 ! moho mesh boundary down to 60km
+!!-- uncomment for e.g. Europe case, where deep moho is rare
+! double precision, parameter :: RMOHO_STRETCH_ADJUSTEMENT = -15000.d0 ! moho mesh boundary down to 55km
+!!-- uncomment for deep moho cases, e.g. Asia case (Himalayan moho)
+! double precision, parameter :: RMOHO_STRETCH_ADJUSTEMENT = -20000.d0 ! moho mesh boundary down to 60km
+!!
+!! mesh tweaking
+!!
! to suppress the crustal layers
! (replaced by an extension of the mantle: R_EARTH is not modified, but no more crustal doubling)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90 2011-03-16 21:18:25 UTC (rev 18118)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90 2011-03-16 21:39:41 UTC (rev 18119)
@@ -471,9 +471,7 @@
if( THREE_D_MODEL == THREE_D_MODEL_GLL ) then
MGLL_V%MODEL_GLL = .true.
! sets to initial reference model from which iterations started
- THREE_D_MODEL = THREE_D_MODEL_S29EA
- ! Hejun Zhu
- !THREE_D_MODEL = THREE_D_MODEL_S362ANI
+ THREE_D_MODEL = GLL_REFERENCE_MODEL
else
MGLL_V%MODEL_GLL = .false.
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_gll.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_gll.f90 2011-03-16 21:18:25 UTC (rev 18118)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_gll.f90 2011-03-16 21:39:41 UTC (rev 18119)
@@ -238,23 +238,17 @@
integer, dimension(MAX_NUM_REGIONS) :: NSPEC
integer :: myrank
- !--------------------------------------------------------------------
- ! USER PARAMETER
-
- character(len=150),parameter:: MGLL_path = 'DATA/GLL/'
- !--------------------------------------------------------------------
-
! local parameters
integer :: ier
character(len=150) :: prname
if( myrank == 0) then
write(IMAIN,*)
- write(IMAIN,*)'reading in model from ',trim(MGLL_path)
+ write(IMAIN,*)'reading in model from ',trim(PATHNAME_GLL_modeldir)
endif
! only crust and mantle
- write(prname,'(a,i6.6,a)') MGLL_path(1:len_trim(MGLL_path))//'proc',myrank,'_reg1_'
+ write(prname,'(a,i6.6,a)') PATHNAME_GLL_modeldir(1:len_trim(PATHNAME_GLL_modeldir))//'proc',myrank,'_reg1_'
! reads in model for each partition
if( .not. TRANSVERSE_ISOTROPY ) then
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90 2011-03-16 21:18:25 UTC (rev 18118)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90 2011-03-16 21:39:41 UTC (rev 18119)
@@ -208,7 +208,7 @@
real(kind=CUSTOM_REAL) :: adj_src(NDIM,NSTEP_BLOCK)
double precision, dimension(NDIM,NSTEP_BLOCK) :: adj_src_u
- integer icomp, itime, i, j, k, ios
+ integer icomp, itime, ios
integer it_start,it_end,index_i
real(kind=CUSTOM_REAL) :: junk
character(len=3),dimension(NDIM) :: comp
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90 2011-03-16 21:18:25 UTC (rev 18118)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90 2011-03-16 21:39:41 UTC (rev 18119)
@@ -1327,6 +1327,13 @@
! initialize seismograms
seismograms(:,:,:) = 0._CUSTOM_REAL
nit_written = 0
+ else
+ ! allocate dummy array since we need it to pass as argument e.g. in write_seismograms() routine
+ ! note: nrec_local is zero, fortran 90/95 should allow zero-sized array allocation...
+ allocate(seismograms(NDIM,nrec_local,NTSTEP_BETWEEN_OUTPUT_SEISMOS),stat=ier)
+ if( ier /= 0) stop 'error while allocating zero seismograms'
+ allocate(number_receiver_global(nrec_local),stat=ier)
+ if( ier /= 0) stop 'error while allocating zero number_receiver_global'
endif
! get information about event name and location for SAC seismograms
More information about the CIG-COMMITS
mailing list