[cig-commits] r18199 - in seismo/3D/SPECFEM3D_GLOBE/trunk: EXAMPLES EXAMPLES/global_s362ani EXAMPLES/regional_MiddleEast UTILS/Cluster/lsf doc/USER_MANUAL setup src/auxiliaries src/meshfem3D src/shared src/specfem3D
danielpeter at geodynamics.org
danielpeter at geodynamics.org
Thu Apr 7 20:43:59 PDT 2011
Author: danielpeter
Date: 2011-04-07 20:43:58 -0700 (Thu, 07 Apr 2011)
New Revision: 18199
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/README
seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/global_s362ani/README
seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/regional_MiddleEast/README
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/Cluster/lsf/go_mesher_solver_lsf_globe.kernel
seismo/3D/SPECFEM3D_GLOBE/trunk/doc/USER_MANUAL/manual_SPECFEM3D_GLOBE.pdf
seismo/3D/SPECFEM3D_GLOBE/trunk/doc/USER_MANUAL/manual_SPECFEM3D_GLOBE.tex
seismo/3D/SPECFEM3D_GLOBE/trunk/setup/config.h.in
seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in
seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_AVS_DX.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_GMT_global.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s362ani.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_receivers.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_surface.f90
Log:
removes trailing spaces from source code files; adds transverse isotropy for STW 1D reference models to have tiso between moho and 670km depth
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/README
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/README 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/README 2011-04-08 03:43:58 UTC (rev 18199)
@@ -10,6 +10,12 @@
- basic usage:
+ * regional_Greece_small/
+
+ contains an example for a small regional simulation and an event located in southern Greece;
+ the example can be run as a small test on a single desktop machine
+ (4 CPUs, forward simulation lasts ~15min, kernel simulation lasts ~40min)
+
* regional_MiddleEast/
contains an example for a regional simulation for the Middle East region
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/global_s362ani/README
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/global_s362ani/README 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/global_s362ani/README 2011-04-08 03:43:58 UTC (rev 18199)
@@ -15,7 +15,7 @@
PBS-job script "go_mesher_solver_pbs.bash".
The simulation is accurate down to a shortest period of ~ 18 s.
-Reference seismograms for comparions are provided in the directory: REF_SEIS/
+Reference seismograms for comparisons are provided in the directory: REF_SEIS/
2. kernel simulation:
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/regional_MiddleEast/README
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/regional_MiddleEast/README 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/EXAMPLES/regional_MiddleEast/README 2011-04-08 03:43:58 UTC (rev 18199)
@@ -15,7 +15,7 @@
PBS-job script "go_mesher_solver_pbs.bash".
The simulation is accurate down to a shortest period of ~ 17 s.
-Reference seismograms for comparions are provided in the directory: REF_SEIS/
+Reference seismograms for comparisons are provided in the directory: REF_SEIS/
2. kernel simulation:
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/Cluster/lsf/go_mesher_solver_lsf_globe.kernel
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/Cluster/lsf/go_mesher_solver_lsf_globe.kernel 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/Cluster/lsf/go_mesher_solver_lsf_globe.kernel 2011-04-08 03:43:58 UTC (rev 18199)
@@ -32,7 +32,9 @@
mv ../OUTPUT_FILES/*.ascii .
rename .semd.ascii .for *.semd.ascii
# P waves
-/home/lqy/bin/utils_bin/xcut_velocity 460 525 3 T60*.for
+#/home/lqy/bin/utils_bin/xcut_velocity 460 525 3 T60*.for
+~/SPECFEM3D_GLOBE/UTILS/adjoint_sources/traveltime/xcreate_adjsrc_traveltime 460. 535. 3 T60*.for
+
rename .for.ad .adj *.for.ad
cd $current_pwd
##########################################################################
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/doc/USER_MANUAL/manual_SPECFEM3D_GLOBE.pdf
===================================================================
(Binary files differ)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/doc/USER_MANUAL/manual_SPECFEM3D_GLOBE.tex
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/doc/USER_MANUAL/manual_SPECFEM3D_GLOBE.tex 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/doc/USER_MANUAL/manual_SPECFEM3D_GLOBE.tex 2011-04-08 03:43:58 UTC (rev 18199)
@@ -1771,7 +1771,6 @@
in \texttt{UTILS/adjoint\_sources/traveltime} to cut a certain portion of the original
displacement seismograms and convert them into the proper adjoint
source to compute the finite-frequency kernel.
-
\begin{lyxcode}
xcreate\_adjsrc\_traveltime~t1~t2~ifile{[}0-5]~E/N/Z-ascii-files~{[}baz]
\end{lyxcode}
@@ -1784,8 +1783,8 @@
back-azimuth of the station from the event location. Note that \texttt{baz}
is only supplied when \texttt{ifile} = 4 or 5.
-\item Similarly, a sample program to compute adjoint sources for amplitude finite-frequency kernels may be found in \texttt{UTILS/adjoint\_sources/amplitude} and used in the same way as described for traveltime measurement
-
+\item Similarly, a sample program to compute adjoint sources for amplitude finite-frequency kernels may be found in \texttt{UTILS/adjoint\_sources/amplitude} and used in the same way as described
+for traveltime measurements
\begin{lyxcode}
xcreate\_adjsrc\_amplitude~t1~t2~ifile{[}0-5]~E/N/Z-ascii-files~{[}baz].
\end{lyxcode}
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/setup/config.h.in
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/setup/config.h.in 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/setup/config.h.in 2011-04-08 03:43:58 UTC (rev 18199)
@@ -33,4 +33,4 @@
#undef PACKAGE_VERSION
/* Uncomment and define to select optimized file i/o for regional simulations */
-#undef USE_MAP_FUNCTION
+#undef USE_MAP_FUNCTION
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/setup/constants.h.in 2011-04-08 03:43:58 UTC (rev 18199)
@@ -147,7 +147,7 @@
integer, parameter :: ICRUST_CRUST2 = 1
integer, parameter :: ICRUST_CRUSTMAPS = 2
integer, parameter :: ICRUST_EPCRUST = 3
-
+
! increase smoothing for critical regions (increases mesh stability)
logical, parameter :: SMOOTH_CRUST = .true.
@@ -258,13 +258,13 @@
! output kernel mask to zero out source region
logical, parameter :: SAVE_SOURCE_MASK = .false.
-! forces transverse isotropy for all mantle elements
+! forces transverse isotropy for all mantle elements
! (default is to use transverse isotropy only between MOHO and 220 )
logical, parameter :: USE_FULL_TISO_MANTLE = .false.
!!-----------------------------------------------------------
!!
-!! time stamp information
+!! time stamp information
!!
!!-----------------------------------------------------------
@@ -539,7 +539,7 @@
!!
!! model constants
-!!
+!!
! number of layers in DATA/1066a/1066a.dat
integer, parameter :: NR_1066A = 160
@@ -588,7 +588,7 @@
integer, parameter :: CRUSTMAP_RESOLUTION = 4 !means 1/4 degrees
integer, parameter :: NLAYERS_CRUSTMAP = 5
-! parameters for EPCRUST , from Molinari & Morelli model(2011)
+! parameters for EPCRUST , from Molinari & Morelli model(2011)
! latitude : 9.0N - 89.5N
! longitude: 56.0W - 70.0E
character(len=*), parameter :: PATHNAME_EPCRUST = 'DATA/epcrust/EPcrust_0_5.txt'
@@ -596,17 +596,17 @@
double precision, parameter :: EPCRUST_LON_MIN = -56.0d0
double precision, parameter :: EPCRUST_LON_MAX = 70.0d0
double precision, parameter :: EPCRUST_LAT_MIN = 9.0d0
- double precision, parameter :: EPCRUST_LAT_MAX = 89.5d0
- double precision, parameter :: EPCRUST_SAMPLE = 0.5d0
+ double precision, parameter :: EPCRUST_LAT_MAX = 89.5d0
+ double precision, parameter :: EPCRUST_SAMPLE = 0.5d0
logical, parameter :: flag_smooth_epcrust = .true.
integer, parameter :: NTHETA_EP = 4, NPHI_EP = 20
- double precision, parameter :: cap_degree_EP = 0.2d0
+ double precision, parameter :: cap_degree_EP = 0.2d0
-! parameters for GLL model (used for iterative model inversions)
+! 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
+! integer, parameter :: GLL_REFERENCE_MODEL = THREE_D_MODEL_S362ANI
!!!!!!!!!!!!!! end of parameters added for the thread-safe version of the code
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_AVS_DX.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_AVS_DX.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_AVS_DX.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -103,7 +103,7 @@
! flag to apply non linear scaling to normalized norm of displacement
logical, parameter :: NONLINEAR_SCALING = .false.
! uses fixed max_value to normalize instead of max of current wavefield
- logical, parameter :: FIX_SCALING = .true.
+ logical, parameter :: FIX_SCALING = .true.
real,parameter:: MAX_VALUE = 1.0
! coefficient of power law used for non linear scaling
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_GMT_global.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_GMT_global.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/auxiliaries/create_movie_GMT_global.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -58,7 +58,7 @@
! muting source region
logical, parameter :: MUTE_SOURCE = .true.
real(kind=CUSTOM_REAL) :: RADIUS_TO_MUTE = 0.5 ! start radius in degrees
- real(kind=CUSTOM_REAL) :: STARTTIME_TO_MUTE = 0.5 ! adds seconds to shift starttime
+ real(kind=CUSTOM_REAL) :: STARTTIME_TO_MUTE = 0.5 ! adds seconds to shift starttime
!---------------------
@@ -158,7 +158,7 @@
real(kind=CUSTOM_REAL) :: cmt_hdur,cmt_t_shift,t0,hdur
real(kind=CUSTOM_REAL) :: xmesh,ymesh,zmesh
integer :: istamp1,istamp2
-
+
! ************** PROGRAM STARTS HERE **************
print *
@@ -209,7 +209,7 @@
print *
if(MOVIE_COARSE) then
! note:
- ! nex_per_proc_xi*nex_per_proc_eta = nex_xi*nex_eta/nproc = nspec2d_top(iregion_crust_mantle)
+ ! nex_per_proc_xi*nex_per_proc_eta = nex_xi*nex_eta/nproc = nspec2d_top(iregion_crust_mantle)
! used in specfem3D.f90
! and ilocnum = nmovie_points = 2 * 2 * NEX_XI * NEX_ETA / NPROC
ilocnum = 2 * 2 * NEX_PER_PROC_XI*NEX_PER_PROC_ETA
@@ -317,7 +317,7 @@
if(ierror /= 0) stop 'error while allocating field_display'
- ! initializes maxima history
+ ! initializes maxima history
if( USE_AVERAGED_MAXIMUM ) then
! determines length of history
nmax_history = AVERAGE_MINIMUM + int( HDUR_MOVIE / (DT*NTSTEP_BETWEEN_FRAMES) * 1.5 )
@@ -334,7 +334,7 @@
print *,'Normalization by averaged maxima over ',nmax_history,'snapshots'
print *
endif
-
+
if( MUTE_SOURCE ) then
! initializes
LAT_SOURCE = -1000.0
@@ -342,7 +342,7 @@
DEP_SOURCE = 0.0
cmt_t_shift = 0.0
cmt_hdur = 0.0
-
+
! reads in source lat/lon
open(22,file="DATA/CMTSOLUTION",status='old',action='read',iostat=ierror )
if( ierror == 0 ) then
@@ -350,10 +350,10 @@
read(22,*,iostat=ierror ) line ! PDE line
read(22,*,iostat=ierror ) line ! event name
! timeshift
- read(22,'(a256)',iostat=ierror ) line
+ read(22,'(a256)',iostat=ierror ) line
if( ierror == 0 ) read(line(12:len_trim(line)),*) cmt_t_shift
! halfduration
- read(22,'(a256)',iostat=ierror ) line
+ read(22,'(a256)',iostat=ierror ) line
if( ierror == 0 ) read(line(15:len_trim(line)),*) cmt_hdur
! latitude
read(22,'(a256)',iostat=ierror ) line
@@ -379,9 +379,9 @@
! note: this starttime is supposed to be the time when displacements at the surface
! can be observed;
! it helps to mute out numerical noise before the source effects actually start showing up
- STARTTIME_TO_MUTE = STARTTIME_TO_MUTE + DEP_SOURCE
+ STARTTIME_TO_MUTE = STARTTIME_TO_MUTE + DEP_SOURCE
if( STARTTIME_TO_MUTE < 0.0 ) STARTTIME_TO_MUTE = 0.0
-
+
print *,'mutes source area'
print *
print *,'source lat/lon/dep: ',LAT_SOURCE,LON_SOURCE,DEP_SOURCE
@@ -447,11 +447,11 @@
mute_factor = 1.0
print*,'simulation time: ',(it-1)*DT - t0,'(s)'
-
+
! muting radius grows/shrinks with time
if( (it-1)*DT - t0 > STARTTIME_TO_MUTE ) then
- ! approximate wavefront travel distance in degrees
+ ! approximate wavefront travel distance in degrees
! (~3.5 km/s wave speed for surface waves)
distance = 3.5 * ((it-1)*DT-t0) / 6371.0 * 180./PI
@@ -465,7 +465,7 @@
! shrinks when waves reached antipode
!if( distance > 180. ) distance = 360. - distance
! shrinks when waves reached half-way to antipode
- if( distance > 90.0 ) distance = 90.0 - distance
+ if( distance > 90.0 ) distance = 90.0 - distance
! limit size around source (in degrees)
if( distance < 0.0 ) distance = 0.0
@@ -481,7 +481,7 @@
! mute factor 1: no masking out
! 0: masks out values (within mute radius)
! linear scaling between [0,1]:
- ! from 0 (simulation time equal to zero )
+ ! from 0 (simulation time equal to zero )
! to 1 (simulation time equals starttime_to_mute)
mute_factor = 1.0 - ( STARTTIME_TO_MUTE - ((it-1)*DT-t0) ) / (STARTTIME_TO_MUTE+t0)
! threshold value for mute_factor
@@ -560,7 +560,7 @@
! distance in colatitude (in rad)
! note: this mixes geocentric (point location) and geographic (source location) coordinates;
- ! since we only need approximate distances here,
+ ! since we only need approximate distances here,
! this should be fine for the muting region
dist_lat = thetaval - LAT_SOURCE
@@ -583,8 +583,8 @@
dist_lon = phival - LON_SOURCE
endif
! distance of point to source (in rad)
- distance = sqrt(dist_lat**2 + dist_lon**2)
-
+ distance = sqrt(dist_lat**2 + dist_lon**2)
+
! mutes source region values
if ( distance < RADIUS_TO_MUTE ) then
! muting takes account of the event time
@@ -683,7 +683,7 @@
if(zmesh < SMALL_VAL_ANGLE .and. zmesh >= ZERO) zmesh = SMALL_VAL_ANGLE
thetaval = atan2(sqrt(xmesh*xmesh+ymesh*ymesh),zmesh)
! thetaval between 0 and PI / 2
- !print*,'thetaval:',thetaval * 180. / PI
+ !print*,'thetaval:',thetaval * 180. / PI
! close to north pole
if( thetaval >= 0.495 * PI ) istamp1 = ieoff
! close to south pole
@@ -744,33 +744,33 @@
if( istamp1 == 0 ) istamp1 = ieoff
if( istamp2 == 0 ) istamp2 = ieoff-1
!print *, 'stamp: ',istamp1,istamp2
-
+
if( max_absol < max_average ) then
! distance (in degree) of surface waves travelled
distance = 3.5 * ((it-1)*DT-t0) / 6371.0 * 180./PI
if( distance > 10.0 .and. distance <= 20.0 ) then
! smooth transition between 10 and 20 degrees
- ! sets positive and negative maximum
+ ! sets positive and negative maximum
field_display(istamp1) = + max_absol + (max_average-max_absol) * (distance - 10.0)/10.0
field_display(istamp2) = - ( max_absol + (max_average-max_absol) * (distance - 10.0)/10.0 )
else if( distance > 20.0 ) then
- ! sets positive and negative maximum
+ ! sets positive and negative maximum
field_display(istamp1) = + max_average
- field_display(istamp2) = - max_average
+ field_display(istamp2) = - max_average
endif
else
! thresholds positive & negative maximum values
where( field_display(:) > max_average ) field_display = max_average
where( field_display(:) < - max_average ) field_display = -max_average
- ! sets positive and negative maximum
+ ! sets positive and negative maximum
field_display(istamp1) = + max_average
- field_display(istamp2) = - max_average
+ field_display(istamp2) = - max_average
endif
! updates current wavefield maxima
min_field_current = minval(field_display(:))
max_field_current = maxval(field_display(:))
- max_absol = (abs(min_field_current)+abs(max_field_current))/2.0
- endif
+ max_absol = (abs(min_field_current)+abs(max_field_current))/2.0
+ endif
! scales field values up to match average
if( abs(max_absol) > TINYVAL) &
@@ -787,33 +787,33 @@
if( (it-1)*DT - t0 > STARTTIME_TO_MUTE ) then
! wavefield should be visible at surface now
! normalizes wavefield
- if( abs(max_average) > TINYVAL ) field_display = field_display / max_average
+ if( abs(max_average) > TINYVAL ) field_display = field_display / max_average
else
! no wavefield yet assumed
! we set two single field values (last in array)
- ! to: +/- 100 * max_average
+ ! to: +/- 100 * max_average
! to avoid further amplifying when
- ! a normalization routine is used for rendering images
+ ! a normalization routine is used for rendering images
! (which for example is the case for shakemovies)
- if( STARTTIME_TO_MUTE > TINYVAL ) then
- ! with additional scale factor:
+ if( STARTTIME_TO_MUTE > TINYVAL ) then
+ ! with additional scale factor:
! linear scaling between [0,1]:
- ! from 0 (simulation time equal to -t0 )
+ ! from 0 (simulation time equal to -t0 )
! to 1 (simulation time equals starttime_to_mute)
mute_factor = 1.0 - ( STARTTIME_TO_MUTE - ((it-1)*DT-t0) ) / (STARTTIME_TO_MUTE+t0)
! takes complement and shifts scale to (1,100)
! thus, mute factor is 100 at simulation start and 1.0 at starttime_to_mute
mute_factor = abs(1.0 - mute_factor) * 99.0 + 1.0
! positive and negative maximum reach average when wavefield appears
- val = mute_factor * max_average
+ val = mute_factor * max_average
else
- ! uses a constant factor
+ ! uses a constant factor
val = 100.0 * max_average
- endif
- ! positive and negative maximum
+ endif
+ ! positive and negative maximum
field_display(istamp1) = + val
- field_display(istamp2) = - val
+ field_display(istamp2) = - val
if( abs(max_average) > TINYVAL ) field_display = field_display / val
endif
else
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/compute_element_properties.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -107,13 +107,13 @@
double precision:: yigll(NGLLY)
double precision:: zigll(NGLLZ)
- logical, dimension(nspec) :: ispec_is_tiso
+ logical, dimension(nspec) :: ispec_is_tiso
! Parameter used to decide whether this element is in the crust or not
logical:: elem_in_crust,elem_in_mantle
! flag for transverse isotropic elements
logical:: elem_is_tiso
-
+
! add topography of the Moho *before* adding the 3D crustal velocity model so that the streched
! mesh gets assigned the right model values
elem_in_crust = .false.
@@ -136,7 +136,7 @@
call moho_stretching_honor_crust(myrank, &
xelm,yelm,zelm,RMOHO_FICTITIOUS_IN_MESHER,&
R220,RMIDDLE_CRUST,elem_in_crust,elem_in_mantle)
- endif
+ endif
else
! element below 220km
! sets element flag for mantle
@@ -148,37 +148,46 @@
if( idoubling(ispec) == IFLAG_CRUST ) then
elem_in_crust = .true.
else
- elem_in_mantle = .true.
- endif
+ elem_in_mantle = .true.
+ endif
endif
-
+
! sets transverse isotropic flag for elements in mantle
if( TRANSVERSE_ISOTROPY ) then
- ! modifies tiso to have it for all mantle elements
+ ! modifies tiso to have it for all mantle elements
! preferred for example, when using 1Dref (STW model)
- ! note: this will increase the computation time by ~ 45 %
if( USE_FULL_TISO_MANTLE ) then
- ! fully transverse isotropic mantle
+ ! all elements below the actual moho will be used for transverse isotropy
+ ! note: this will increase the computation time by ~ 45 %
+ if( elem_in_mantle ) then
+ elem_is_tiso = .true.
+ endif
+ else if( REFERENCE_1D_MODEL == REFERENCE_MODEL_1DREF ) then
+ ! transverse isotropic mantle between fictitious moho to 670km depth
! preferred for harvard (kustowski's) models using STW 1D reference, i.e.
- ! THREE_D_MODEL_S362ANI
- ! THREE_D_MODEL_S362WMANI
+ ! THREE_D_MODEL_S362ANI
+ ! THREE_D_MODEL_S362WMANI
! THREE_D_MODEL_S29EA
! THREE_D_MODEL_GLL
- if( elem_in_mantle ) elem_is_tiso = .true.
+ ! which show significant transverse isotropy also below 220km depth
+ if( idoubling(ispec)==IFLAG_220_80 .or. idoubling(ispec)==IFLAG_80_MOHO &
+ .or. idoubling(ispec)==IFLAG_670_220 ) then
+ elem_is_tiso = .true.
+ endif
else if( idoubling(ispec)==IFLAG_220_80 .or. idoubling(ispec)==IFLAG_80_MOHO ) then
- ! default case:
+ ! default case for PREM reference models:
! models use only transverse isotropy between moho and 220 km depth
elem_is_tiso = .true.
! checks mantle flag to be sure
- if( elem_in_mantle .eqv. .false. ) stop 'error mantle flag confused between moho and 220'
+ !if( elem_in_mantle .eqv. .false. ) stop 'error mantle flag confused between moho and 220'
endif
endif
-
+
endif ! IREGION_CRUST_MANTLE
-
+
! sets element tiso flag
ispec_is_tiso(ispec) = elem_is_tiso
-
+
! interpolates and stores GLL point locations
call compute_element_GLL_locations(xelm,yelm,zelm,ispec,nspec, &
xstore,ystore,zstore,shape3D)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_central_cube.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -121,7 +121,7 @@
logical :: ACTUALLY_STORE_ARRAYS,ABSORBING_CONDITIONS
- logical, dimension(nspec) :: ispec_is_tiso
+ logical, dimension(nspec) :: ispec_is_tiso
!local parameters
double precision, dimension(NGNOD) :: xelm,yelm,zelm
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_doubling_elements.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -155,7 +155,7 @@
integer :: offset_proc_xi,offset_proc_eta
logical :: CUT_SUPERBRICK_XI,CUT_SUPERBRICK_ETA
- logical, dimension(nspec) :: ispec_is_tiso
+ logical, dimension(nspec) :: ispec_is_tiso
! local parameters
double precision, dimension(NGLOB_DOUBLING_SUPERBRICK) :: x_superbrick,y_superbrick,z_superbrick
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regions_mesh.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -227,7 +227,7 @@
double precision r_moho,r_400,r_670
! flags for transverse isotropic elements
- logical, dimension(:), allocatable :: ispec_is_tiso
+ logical, dimension(:), allocatable :: ispec_is_tiso
! create the name for the database of the current slide and region
call create_name_database(prname,myrank,iregion_code,LOCAL_PATH)
@@ -244,7 +244,7 @@
allocate(Qmu_store(NGLLX,NGLLY,NGLLZ,nspec_att), &
tau_e_store(N_SLS,NGLLX,NGLLY,NGLLZ,nspec_att),stat=ier)
if(ier /= 0) stop 'error in allocate 1'
-
+
! Gauss-Lobatto-Legendre points of integration
allocate(xigll(NGLLX), &
yigll(NGLLY), &
@@ -286,7 +286,7 @@
eta_anisostore(NGLLX,NGLLY,NGLLZ,nspec), &
ispec_is_tiso(nspec),stat=ier)
if(ier /= 0) stop 'error in allocate 7'
-
+
! initializes flags for transverse isotropic elements
ispec_is_tiso(:) = .false.
@@ -393,7 +393,7 @@
gammaystore(NGLLX,NGLLY,NGLLZ,nspec_actually), &
gammazstore(NGLLX,NGLLY,NGLLZ,nspec_actually),stat=ier)
if(ier /= 0) stop 'error in allocate 16'
-
+
! boundary mesh
if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
NSPEC2D_MOHO = NSPEC2D_TOP
@@ -636,7 +636,7 @@
yp(npointot), &
zp(npointot),stat=ier)
if(ier /= 0) stop 'error in allocate 20'
-
+
locval = 0
ifseg = .false.
xp = 0.d0
@@ -746,7 +746,7 @@
! (used only for checks in meshfem3D() routine)
!nspec_tiso = count(idoubling(1:nspec) == IFLAG_220_80) + count(idoubling(1:nspec) == IFLAG_80_MOHO)
nspec_tiso = count(ispec_is_tiso(:))
-
+
! precomputes jacobian for 2d absorbing boundary surfaces
call get_jacobian_boundaries(myrank,iboun,nspec,xstore,ystore,zstore, &
dershape2D_x,dershape2D_y,dershape2D_bottom,dershape2D_top, &
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/create_regular_elements.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -161,7 +161,7 @@
integer ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top, &
ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot
- logical, dimension(nspec) :: ispec_is_tiso
+ logical, dimension(nspec) :: ispec_is_tiso
! local parameters
double precision, dimension(NGNOD) :: offset_x,offset_y,offset_z
@@ -266,7 +266,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
nspec_ani,nspec_stacey,nspec_att,Qmu_store,tau_e_store,tau_s,T_c_source,&
rho_vp,rho_vs,ACTUALLY_STORE_ARRAYS,&
- xigll,yigll,zigll,ispec_is_tiso)
+ xigll,yigll,zigll,ispec_is_tiso)
! boundary mesh
if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -686,7 +686,7 @@
if(iregion_code == IREGION_CRUST_MANTLE .and. nspec_tiso == 0) &
call exit_MPI(myrank,'found no anisotropic elements in the mantle')
endif
-
+
! computes total area and volume
call meshfem3D_compute_area(myrank,NCHUNKS,iregion_code, &
area_local_bottom,area_local_top,&
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/meshfem3D_models.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -286,17 +286,17 @@
end type model_eucrust_variables
type (model_eucrust_variables) EUCM_V
-! type for EPCRUST 1.0
- type model_epcrust_variables
+! type for EPCRUST 1.0
+ type model_epcrust_variables
sequence
- double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
+ double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: thickness_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vp_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vs_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: rho_ep
- end type model_epcrust_variables
+ end type model_epcrust_variables
type (model_epcrust_variables) EPCRUST
-
+
! model_crustmaps_variables combined crustal maps
type model_crustmaps_variables
sequence
@@ -620,7 +620,7 @@
case (ICRUST_EPCRUST)
! EPcrust
- call model_epcrust_broadcast(myrank,EPCRUST)
+ call model_epcrust_broadcast(myrank,EPCRUST)
case default
stop 'crustal model type not defined'
@@ -1188,17 +1188,17 @@
call model_crustmaps(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,GC_V,elem_in_crust)
case (ICRUST_EPCRUST)
-! call model_crust(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,CM_V,elem_in_crust)
- ! within EPCRUST region
+! call model_crust(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,CM_V,elem_in_crust)
+ ! within EPCRUST region
! if (lat >= EPCRUST_LAT_MIN .and. lat <= EPCRUST_LAT_MAX &
! .and. lon >= EPCRUST_LON_MIN .and. lon<=EPCRUST_LON_MAX ) then
! vpc=0.0d0
! vsc=0.0d0
! rhoc=0.0d0
! moho=0.0d0
-! found_crust = .false.
- call model_epcrust(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,EPCRUST,elem_in_crust)
-! end if
+! found_crust = .false.
+ call model_epcrust(lat,lon,r,vpc,vsc,rhoc,moho,found_crust,EPCRUST,elem_in_crust)
+! end if
case default
stop 'crustal model type not defined'
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_epcrust.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -27,34 +27,34 @@
!--------------------------------------------------------------------------------------------------
! EPCRUST 1.0
-! I. Molinari and A. Morelli, 2011.
+! I. Molinari and A. Morelli, 2011.
! EPcrust: A reference crustal model for the European plate
-! GJI, 185 (1), pages 352-364
+! GJI, 185 (1), pages 352-364
!--------------------------------------------------------------------------------------------------
subroutine model_epcrust_broadcast(myrank,EPCRUST)
-
- implicit none
+ implicit none
+
include "constants.h"
include 'mpif.h'
- type model_epcrust_variables
+ type model_epcrust_variables
sequence
- double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
+ double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: thickness_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vp_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vs_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: rho_ep
- end type model_epcrust_variables
+ end type model_epcrust_variables
type (model_epcrust_variables) EPCRUST
- integer:: myrank,ierr
+ integer:: myrank,ierr
- ! read EPCRUST model on master
- if(myrank == 0) call read_epcrust_model(EPCRUST)
+ ! read EPCRUST model on master
+ if(myrank == 0) call read_epcrust_model(EPCRUST)
- ! broadcast EPCRUST model
+ ! broadcast EPCRUST model
call MPI_BCAST(EPCRUST%lon_ep,EPCRUST_NLON*EPCRUST_NLAT, &
MPI_DOUBLE_PRECISION,0,MPI_COMM_WORLD,ierr)
call MPI_BCAST(EPCRUST%lat_ep,EPCRUST_NLON*EPCRUST_NLAT, &
@@ -76,19 +76,19 @@
!-------------------------------------------------------------------------------------------------
!
- subroutine read_epcrust_model(EPCRUST)
-
- implicit none
- include "constants.h"
+ subroutine read_epcrust_model(EPCRUST)
- type model_epcrust_variables
+ implicit none
+ include "constants.h"
+
+ type model_epcrust_variables
sequence
- double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
+ double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: thickness_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vp_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vs_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: rho_ep
- end type model_epcrust_variables
+ end type model_epcrust_variables
type (model_epcrust_variables) EPCRUST
character(len=150) EPCRUST_FNM
@@ -98,11 +98,11 @@
call get_value_string(EPCRUST_FNM,'model.EPCRUST_FNM',PATHNAME_EPCRUST)
open(unit=1001,file=EPCRUST_FNM,status='old',action='read')
- read(1001,*) header
+ read(1001,*) header
do jlat = 1,EPCRUST_NLAT
do ilon=1,EPCRUST_NLON
- read(1001,*) tmp
+ read(1001,*) tmp
EPCRUST%lon_ep(ilon,jlat) = tmp(1)
EPCRUST%lat_ep(ilon,jlat) = tmp(2)
@@ -110,51 +110,51 @@
EPCRUST%thickness_ep(ilon,jlat,1:3) = tmp(4:6)
EPCRUST%vp_ep(ilon,jlat,1:3) = tmp(7:9)
EPCRUST%vs_ep(ilon,jlat,1:3) = tmp(10:12)
- EPCRUST%rho_ep(ilon,jlat,1:3) = tmp(13:15)
- end do
- end do
- close(1001)
+ EPCRUST%rho_ep(ilon,jlat,1:3) = tmp(13:15)
+ end do
+ end do
+ close(1001)
- end subroutine read_epcrust_model
+ end subroutine read_epcrust_model
!
!-------------------------------------------------------------------------------------------------
!
- subroutine model_epcrust(lat,lon,dep,vp,vs,rho,moho,found_crust,EPCRUST,elem_in_crust)
- implicit none
- include "constants.h"
+ subroutine model_epcrust(lat,lon,dep,vp,vs,rho,moho,found_crust,EPCRUST,elem_in_crust)
+ implicit none
+ include "constants.h"
- ! INPUT & OUTPUT
- type model_epcrust_variables
+ ! INPUT & OUTPUT
+ type model_epcrust_variables
sequence
- double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
+ double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT):: lon_ep,lat_ep,topo_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: thickness_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vp_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: vs_ep
double precision,dimension(EPCRUST_NLON,EPCRUST_NLAT,EPCRUST_NLAYER):: rho_ep
- end type model_epcrust_variables
+ end type model_epcrust_variables
type (model_epcrust_variables) EPCRUST
double precision:: lat, lon, dep, vp, vs, rho, moho
- logical :: found_crust, elem_in_crust
+ logical :: found_crust, elem_in_crust
- ! INTERIOR
+ ! INTERIOR
integer:: ilon, jlat, k
- double precision :: z0 , topo, basement, conrad, moho_top, scaleval
- double precision,dimension(3):: zsmooth, vpsmooth, vssmooth, rhosmooth
- double precision,dimension(NTHETA_EP*NPHI_EP) :: x1,y1,weight
- double precision:: weightl
+ double precision :: z0 , topo, basement, conrad, moho_top, scaleval
+ double precision,dimension(3):: zsmooth, vpsmooth, vssmooth, rhosmooth
+ double precision,dimension(NTHETA_EP*NPHI_EP) :: x1,y1,weight
+ double precision:: weightl
double precision:: cut, min_sed
!if ( lat < EPCRUST_LAT_MIN .or. lat > EPCRUST_LAT_MAX &
- ! .or. lon < EPCRUST_LON_MIN .or. lon > EPCRUST_LON_MAX ) then
+ ! .or. lon < EPCRUST_LON_MIN .or. lon > EPCRUST_LON_MAX ) then
! stop 'incorrect enter EPCRUST model, check lat and lon'
- !end if
+ !end if
vp=0.0d0
vs=0.0d0
- rho=0.0d0
+ rho=0.0d0
if ( .not. flag_smooth_epcrust) then
call ilon_jlat(lon,lat,ilon,jlat)
@@ -163,14 +163,14 @@
vpsmooth(:)=EPCRUST%vp_ep(ilon,jlat,:)
vssmooth(:)=EPCRUST%vs_ep(ilon,jlat,:)
rhosmooth(:)=EPCRUST%rho_ep(ilon,jlat,:)
- else
+ else
call epcrust_smooth_base(lon,lat,x1,y1,weight)
z0=0.d0
zsmooth(:)=0.0d0
vpsmooth(:)=0.0d0
vssmooth(:)=0.0d0
rhosmooth(:)=0.0d0
-
+
do k = 1,NTHETA_EP*NPHI_EP
call ilon_jlat(x1(k),y1(k),ilon,jlat)
weightl=weight(k)
@@ -180,7 +180,7 @@
vpsmooth(:)=vpsmooth(:)+weightl*EPCRUST%vp_ep(ilon,jlat,:)
vssmooth(:)=vssmooth(:)+weightl*EPCRUST%vs_ep(ilon,jlat,:)
rhosmooth(:)=rhosmooth(:)+weightl*EPCRUST%rho_ep(ilon,jlat,:)
- end do
+ end do
end if
!topo=(R_EARTH_KM+z0)/R_EARTH_KM
@@ -197,7 +197,7 @@
found_crust=.true.
if ( dep > basement .and. INCLUDE_SEDIMENTS_CRUST &
- .and. zsmooth(1) >= MINIMUM_SEDIMENT_THICKNESS ) then ! Hejun Zhu add minimum sediment thickness
+ .and. zsmooth(1) >= MINIMUM_SEDIMENT_THICKNESS ) then ! Hejun Zhu add minimum sediment thickness
vp=vpsmooth(1)
vs=vssmooth(1)
rho=rhosmooth(1)
@@ -205,29 +205,29 @@
vp=vpsmooth(2)
vs=vssmooth(2)
rho=rhosmooth(2)
- else if ( dep > moho_top .or. elem_in_crust ) then
+ else if ( dep > moho_top .or. elem_in_crust ) then
vp=vpsmooth(3)
vs=vssmooth(3)
rho=rhosmooth(3)
else
- found_crust=.false.
- end if
+ found_crust=.false.
+ end if
- if (found_crust ) then
+ if (found_crust ) then
scaleval=dsqrt(PI*GRAV*RHOAV)
vp=vp*1000.d0/(R_EARTH*scaleval)
vs=vs*1000.d0/(R_EARTH*scaleval)
rho=rho*1000.d0/RHOAV
end if
- !moho = -(z0-zsmooth(1)-zsmooth(2)-zsmooth(3))/R_EARTH_KM
- moho = (zsmooth(1)+zsmooth(2)+zsmooth(3))/R_EARTH_KM
+ !moho = -(z0-zsmooth(1)-zsmooth(2)-zsmooth(3))/R_EARTH_KM
+ moho = (zsmooth(1)+zsmooth(2)+zsmooth(3))/R_EARTH_KM
- ! Hejun Zhu, delete moho thickness less than 7 km
+ ! Hejun Zhu, delete moho thickness less than 7 km
cut=7.0/R_EARTH_KM
- if ( moho < cut ) then
- moho = cut
- end if
+ if ( moho < cut ) then
+ moho = cut
+ end if
end subroutine model_epcrust
@@ -237,17 +237,17 @@
subroutine epcrust_smooth_base(x,y,x1,y1,weight)
-
+
implicit none
include "constants.h"
- ! INPUT & OUTPUT
+ ! INPUT & OUTPUT
double precision:: x, y
double precision,dimension(NTHETA_EP*NPHI_EP):: x1,y1,weight
! INTERIOR
double precision:: CAP,dtheta,dphi,cap_area,dweight,pi_over_nphi,total,wght
- double precision:: theta,phi,sint,cost,sinp,cosp
+ double precision:: theta,phi,sint,cost,sinp,cosp
double precision:: r_rot,theta_rot,phi_rot
double precision,dimension(3,3):: rotation_matrix
double precision,dimension(3):: xx,xc
@@ -259,18 +259,18 @@
y1(:)=0.0d0
weight(:)=0.0d0
- if (cap_degree_EP < TINYVAL ) then
+ if (cap_degree_EP < TINYVAL ) then
print*, 'error cap:', cap_degree_EP
print*, 'lat/lon:', x,y
stop 'error cap_degree too small'
- end if
+ end if
CAP=cap_degree_EP*PI/180.0d0
dtheta=0.5d0*CAP/dble(NTHETA_EP)
dphi=TWO_PI/dble(NPHI_EP)
cap_area=TWO_PI*(1.0d0-dcos(CAP))
dweight=CAP/dble(NTHETA_EP)*dphi/cap_area
- pi_over_nphi=PI/dble(NPHI_EP)
+ pi_over_nphi=PI/dble(NPHI_EP)
phi=x*DEGREES_TO_RADIANS
theta=(90.0d0-y)*DEGREES_TO_RADIANS
@@ -313,30 +313,30 @@
xx(j)=0.0d0
do k = 1,3
xx(j)=xx(j)+rotation_matrix(j,k)*xc(k)
- end do
- end do
+ end do
+ end do
call xyz_2_rthetaphi_dble(xx(1),xx(2),xx(3),r_rot,theta_rot,phi_rot)
call reduce(theta_rot,phi_rot)
x1(i)=phi_rot*RADIANS_TO_DEGREES
y1(i)=(PI_OVER_TWO-theta_rot)*RADIANS_TO_DEGREES
if (x1(i) > 180.d0) x1(i)=x1(i)-360.d0
- end do
- end do
+ end do
+ end do
- if (abs(total-1.0d0) > 0.001d0) then
+ if (abs(total-1.0d0) > 0.001d0) then
print*,'error cap:',total,cap_degree_EP
stop
- end if
+ end if
- end subroutine epcrust_smooth_base
+ end subroutine epcrust_smooth_base
!
!-------------------------------------------------------------------------------------------------
!
subroutine ilon_jlat(lon,lat,ilon,jlat)
-
- implicit none
+
+ implicit none
include "constants.h"
double precision:: lon,lat
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s362ani.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s362ani.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/model_s362ani.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -1116,7 +1116,7 @@
real(kind=4) verlon(nver)
real(kind=4) verrad(nver)
- integer, parameter :: maxver=1000
+ integer, parameter :: maxver=1000
integer icon(maxver)
real(kind=4) con(maxver)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/meshfem3D/save_arrays_solver.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -152,7 +152,7 @@
logical ABSORBING_CONDITIONS,SAVE_MESH_FILES
- logical, dimension(nspec) :: ispec_is_tiso
+ logical, dimension(nspec) :: ispec_is_tiso
! local parameters
integer i,j,k,ispec,iglob,nspec1, nglob1
@@ -341,7 +341,7 @@
write(27) is_on_a_slice_edge
write(27) ispec_is_tiso
-
+
close(27)
! absorbing boundary parameters
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/shared/memory_eval.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -243,11 +243,11 @@
! idoubling_crust_mantle (not needed anymore..)
! static_memory_size = static_memory_size + NSPEC(IREGION_CRUST_MANTLE)*dble(SIZE_INTEGER)
-! idoubling_outer_core
+! idoubling_outer_core
static_memory_size = static_memory_size + NSPEC(IREGION_OUTER_CORE)*dble(SIZE_INTEGER)
! idoubling_inner_core
static_memory_size = static_memory_size + NSPEC(IREGION_INNER_CORE)*dble(SIZE_INTEGER)
-
+
! ispec_is_tiso_crust_mantle
static_memory_size = static_memory_size + NSPEC(IREGION_CRUST_MANTLE)*dble(SIZE_LOGICAL)
! ispec_is_tiso_outer_core
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_arrays_source.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -204,8 +204,8 @@
double precision :: hxir(NGLLX), hpxir(NGLLX), hetar(NGLLY), hpetar(NGLLY), &
hgammar(NGLLZ), hpgammar(NGLLZ)
double precision, dimension(NDIM,NGLLX,NGLLY,NGLLZ) :: sourcearrayd
-
- real(kind=CUSTOM_REAL) :: adj_src(NDIM,NSTEP_BLOCK)
+
+ real(kind=CUSTOM_REAL) :: adj_src(NDIM,NSTEP_BLOCK)
double precision, dimension(NDIM,NSTEP_BLOCK) :: adj_src_u
integer icomp, itime, ios
@@ -313,16 +313,16 @@
! adds interpolated source contribution to all GLL points within this element
do itime = 1, NSTEP_BLOCK
- ! multiply with interpolators
- call multiply_arrays_adjoint(sourcearrayd,hxir,hetar,hgammar,adj_src_u(:,itime))
-
+ ! multiply with interpolators
+ call multiply_arrays_adjoint(sourcearrayd,hxir,hetar,hgammar,adj_src_u(:,itime))
+
! distinguish between single and double precision for reals
if(CUSTOM_REAL == SIZE_REAL) then
adj_sourcearray(:,:,:,:,itime) = sngl(sourcearrayd(:,:,:,:))
else
adj_sourcearray(:,:,:,:,itime) = sourcearrayd(:,:,:,:)
endif
-
+
enddo
! do k = 1, NGLLZ
! do j = 1, NGLLY
@@ -352,9 +352,9 @@
double precision, dimension(NGLLY) :: hetar
double precision, dimension(NGLLZ) :: hgammar
double precision, dimension(NDIM) :: adj_src_ud
-
+
integer :: i,j,k
-
+
! adds interpolated source contribution to all GLL points within this element
do k = 1, NGLLZ
do j = 1, NGLLY
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_boundary_kernel.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -258,7 +258,7 @@
integer, dimension(NGLLX,NGLLY,NGLLZ,*) :: ibool
! integer, dimension(*) :: idoubling
logical, dimension(*) :: ispec_is_tiso
-
+
! --- local variables ---
real(kind=CUSTOM_REAL) :: duxdxl_plus_duydyl,duxdxl_plus_duzdzl,duydyl_plus_duzdzl
real(kind=CUSTOM_REAL) :: duxdyl_plus_duydxl,duzdxl_plus_duxdzl,duzdyl_plus_duydzl
@@ -339,7 +339,7 @@
!else if(.not. (TRANSVERSE_ISOTROPY_VAL .and. (idoubling(ispec) == IFLAG_80_MOHO .or. idoubling(ispec) == IFLAG_220_80))) then
! isotropic elements
-
+
kappal = kappavstore(i,j,k,ispec)
mul = muvstore(i,j,k,ispec)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/compute_forces_crust_mantle_Dev.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -548,7 +548,7 @@
! if(.not. (TRANSVERSE_ISOTROPY_VAL .and. (idoubling(ispec)==IFLAG_220_80 &
! .or. idoubling(ispec)==IFLAG_80_MOHO))) then
if( .not. ispec_is_tiso(ispec) ) then
-
+
! layer with no transverse isotropy, use kappav and muv
kappal = kappavstore(i,j,k,ispec)
mul = muvstore(i,j,k,ispec)
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_receivers.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_receivers.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/locate_receivers.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -612,7 +612,7 @@
write(IMAIN,*) '*****************************************************************'
else
nrec_found = nrec_found + 1
-
+
islice_selected_rec_found(nrec_found) = islice_selected_rec(irec)
ispec_selected_rec_found(nrec_found) = ispec_selected_rec(irec)
xi_receiver_found(nrec_found) = xi_receiver(irec)
@@ -626,9 +626,9 @@
stbur_found(nrec_found) = stbur(irec)
nu_found(:,:,nrec_found) = nu(:,:,irec)
epidist_found(nrec_found) = epidist(irec)
-
+
! writes out actual receiver location to vtk file
- write(IOVTK,*) sngl(x_found(irec)), sngl(y_found(irec)), sngl(z_found(irec))
+ write(IOVTK,*) sngl(x_found(irec)), sngl(y_found(irec)), sngl(z_found(irec))
endif
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/noise_tomography.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -124,7 +124,7 @@
open(unit=IOUT_NOISE,file='OUTPUT_FILES/irec_master_noise', &
status='unknown',action='write',iostat=ios)
if( ios /= 0 ) call exit_MPI(myrank,'error opening output file irec_master_noise')
-
+
WRITE(IOUT_NOISE,*) 'The master receiver is: (RECEIVER ID)', irec_master_noise
close(IOUT_NOISE)
endif
@@ -204,14 +204,14 @@
open(unit=IOUT_NOISE,file='OUTPUT_FILES/mask_noise', &
status='unknown',form='unformatted',action='write',iostat=ios)
if( ios /= 0 ) call exit_MPI(myrank,'error opening output file mask_noise')
-
+
write(IOUT_NOISE) store_val_x_all
write(IOUT_NOISE) store_val_y_all
write(IOUT_NOISE) store_val_z_all
write(IOUT_NOISE) store_val_ux_all
write(IOUT_NOISE) store_val_uy_all
write(IOUT_NOISE) store_val_uz_all
-
+
close(IOUT_NOISE)
endif
!!!END!!! save mask_noise for check, a file called "mask_noise" is saved in "./OUTPUT_FIELS/"
@@ -246,7 +246,7 @@
open(unit=IOUT_NOISE,file='OUTPUT_FILES/NOISE_SIMULATION', &
status='unknown',action='write',iostat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error opening output file NOISE_SIMULATION')
-
+
WRITE(IOUT_NOISE,*) '*******************************************************************************'
WRITE(IOUT_NOISE,*) '*******************************************************************************'
WRITE(IOUT_NOISE,*) 'WARNING!!!!!!!!!!!!'
@@ -465,7 +465,7 @@
real(kind=CUSTOM_REAL), dimension(NDIM,NGLOB_CRUST_MANTLE) :: displ_crust_mantle
! output parameters
! local parameters
- integer :: ispec2D,ispec,i,j,k,iglob
+ integer :: ispec2D,ispec,i,j,k,iglob
real(kind=CUSTOM_REAL),dimension(NDIM,NGLLX,NGLLY,nspec_top) :: noise_surface_movie
! get coordinates of surface mesh and surface displacement
@@ -729,7 +729,7 @@
integer :: i,j,k,ispec,iglob,ipoin,ispec2D
real(kind=CUSTOM_REAL) :: eta
real(kind=CUSTOM_REAL), dimension(NDIM,NGLLX,NGLLY,nspec_top) :: noise_surface_movie
-
+
! read surface movie, needed for Sigma_kl_crust_mantle
call read_abs(9,noise_surface_movie,CUSTOM_REAL*NDIM*NGLLX*NGLLY*nspec_top,it)
@@ -853,10 +853,10 @@
open(unit=IOUT_NOISE,file=trim(prname)//'Sigma_kernel.bin', &
status='unknown',form='unformatted',action='write',iostat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error opening file Sigma_kernel.bin')
-
+
write(IOUT_NOISE) Sigma_kl_crust_mantle ! need to put dimensions back (not done yet)
close(IOUT_NOISE)
-
+
end subroutine save_kernels_strength_noise
! =============================================================================================================
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/prepare_timerun.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -503,7 +503,7 @@
c33store_crust_mantle,c44store_crust_mantle, &
c55store_crust_mantle,c66store_crust_mantle, &
muvstore_crust_mantle,muhstore_crust_mantle,ispec_is_tiso_crust_mantle, &
- !----- idoubling_crust_mantle, &
+ !----- idoubling_crust_mantle, &
muvstore_inner_core, &
SIMULATION_TYPE,MOVIE_VOLUME,muvstore_crust_mantle_3dmovie, &
c11store_inner_core,c12store_inner_core,c13store_inner_core, &
@@ -646,8 +646,8 @@
muvstore_crust_mantle_3dmovie(i,j,k,ispec)=muvstore_crust_mantle(i,j,k,ispec)
endif
muvstore_crust_mantle(i,j,k,ispec) = muvstore_crust_mantle(i,j,k,ispec) * scale_factor
-
- ! scales transverse isotropic values for mu_h
+
+ ! scales transverse isotropic values for mu_h
!if(TRANSVERSE_ISOTROPY_VAL .and. (idoubling_crust_mantle(ispec) == IFLAG_220_80 &
! .or. idoubling_crust_mantle(ispec) == IFLAG_80_MOHO)) &
if( ispec_is_tiso_crust_mantle(ispec) ) then
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_arrays_solver.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -93,7 +93,7 @@
! this for non blocking MPI
logical, dimension(nspec) :: is_on_a_slice_edge
- logical, dimension(nspec) :: ispec_is_tiso
+ logical, dimension(nspec) :: ispec_is_tiso
! processor identification
character(len=150) prname
@@ -195,7 +195,7 @@
read(IIN) is_on_a_slice_edge
read(IIN) ispec_is_tiso
-
+
close(IIN)
end subroutine read_arrays_solver
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/read_mesh_databases.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -41,7 +41,7 @@
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
- ! -- idoubling_crust_mantle
+ ! -- idoubling_crust_mantle
is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
vp_outer_core,xstore_outer_core,ystore_outer_core,zstore_outer_core, &
xix_outer_core,xiy_outer_core,xiz_outer_core, &
@@ -97,10 +97,10 @@
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_CRUST_MANTLE) :: ibool_crust_mantle
-
+
! integer, dimension(NSPEC_CRUST_MANTLE) :: idoubling_crust_mantle
logical, dimension(NSPEC_CRUST_MANTLE) :: ispec_is_tiso_crust_mantle
-
+
! mass matrix
real(kind=CUSTOM_REAL), dimension(NGLOB_CRUST_MANTLE) :: rmass_crust_mantle
! additional mass matrix for ocean load
@@ -120,7 +120,7 @@
integer, dimension(NGLLX,NGLLY,NGLLZ,NSPEC_OUTER_CORE) :: ibool_outer_core
integer, dimension(NSPEC_OUTER_CORE) :: idoubling_outer_core
logical, dimension(NSPEC_OUTER_CORE) :: ispec_is_tiso_outer_core
-
+
real(kind=CUSTOM_REAL), dimension(NGLOB_OUTER_CORE) :: rmass_outer_core
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPEC_INNER_CORE) :: &
@@ -146,7 +146,7 @@
logical READ_KAPPA_MU,READ_TISO
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,1) :: dummy_array
integer, dimension(NSPEC_CRUST_MANTLE) :: dummy_i
-
+
! this for non blocking MPI
logical, dimension(NSPEC_CRUST_MANTLE) :: is_on_a_slice_edge_crust_mantle
logical, dimension(NSPEC_OUTER_CORE) :: is_on_a_slice_edge_outer_core
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/setup_sources_receivers.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -336,7 +336,7 @@
! creates source/receiver location file
filename = trim(OUTPUT_FILES)//'/sr_tmp.vtk'
filename_new = trim(OUTPUT_FILES)//'/sr.vtk'
- write(system_command, &
+ write(system_command, &
"('sed -e ',a1,'s/POINTS.*/POINTS',i6,' float/',a1,' < ',a,' > ',a)")&
"'",NSOURCES + nrec,"'",trim(filename),trim(filename_new)
call system(system_command)
@@ -347,7 +347,7 @@
"('awk ',a1,'{if(NR<5) print $0;if(NR==6)print ',a1,'POINTS',i6,' float',a1,';if(NR>5+',i6,')print $0}',a1,' < ',a,' > ',a)")&
"'",'"',nrec,'"',NSOURCES,"'",trim(filename),trim(filename_new)
call system(system_command)
-
+
! only extract source locations and remove temporary file
filename_new = trim(OUTPUT_FILES)//'/source.vtk'
write(system_command, &
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/specfem3D.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -871,7 +871,7 @@
real(kind=CUSTOM_REAL), dimension(:,:,:,:,:), allocatable :: noise_sourcearray
real(kind=CUSTOM_REAL), dimension(:), allocatable :: &
normal_x_noise,normal_y_noise,normal_z_noise, mask_noise
- real(kind=CUSTOM_REAL), dimension(:,:,:,:), allocatable :: noise_surface_movie
+ real(kind=CUSTOM_REAL), dimension(:,:,:,:), allocatable :: noise_surface_movie
integer :: irec_master_noise
! ************** PROGRAM STARTS HERE **************
@@ -997,7 +997,7 @@
c44store_crust_mantle,c45store_crust_mantle,c46store_crust_mantle, &
c55store_crust_mantle,c56store_crust_mantle,c66store_crust_mantle, &
ibool_crust_mantle,ispec_is_tiso_crust_mantle, &
- ! -- idoubling_crust_mantle,
+ ! -- idoubling_crust_mantle,
is_on_a_slice_edge_crust_mantle,rmass_crust_mantle,rmass_ocean_load, &
vp_outer_core,xstore_outer_core,ystore_outer_core,zstore_outer_core, &
xix_outer_core,xiy_outer_core,xiz_outer_core, &
@@ -1283,7 +1283,7 @@
NSTEP_SUB_ADJ = ceiling( dble(NSTEP)/dble(NTSTEP_BETWEEN_READ_ADJSRC) )
allocate(iadj_vec(NSTEP),stat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error allocating iadj_vec')
-
+
! initializes iadj_vec
do it=1,NSTEP
iadj_vec(it) = NSTEP-it+1 ! default is for reversing entire record
@@ -1349,7 +1349,7 @@
allocate(moment_der(NDIM,NDIM,nrec_local),sloc_der(NDIM,nrec_local), &
stshift_der(nrec_local),shdur_der(nrec_local),stat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error allocating frechet derivatives arrays')
-
+
moment_der = 0._CUSTOM_REAL
sloc_der = 0._CUSTOM_REAL
stshift_der = 0._CUSTOM_REAL
@@ -1363,9 +1363,9 @@
! 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'
+ 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'
+ if( ier /= 0) stop 'error while allocating zero number_receiver_global'
endif
! get information about event name and location for SAC seismograms
@@ -1598,7 +1598,7 @@
store_val_uy(nmovie_points), &
store_val_uz(nmovie_points),stat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error allocating movie surface arrays')
-
+
if (MOVIE_SURFACE) then ! those arrays are not neccessary for noise tomography, so only allocate them in MOVIE_SURFACE case
allocate(store_val_x_all(nmovie_points,0:NPROCTOT_VAL-1), &
store_val_y_all(nmovie_points,0:NPROCTOT_VAL-1), &
@@ -1646,7 +1646,7 @@
if(myrank == 0) then
write(IMAIN,*)
- write(IMAIN,*) 'Movie volume:'
+ write(IMAIN,*) 'Movie volume:'
write(IMAIN,*) ' Writing to movie3D*** files on local disk databases directory'
if(MOVIE_VOLUME_TYPE == 1) then
write(IMAIN,*) ' movie output: strain'
@@ -1655,7 +1655,7 @@
else if(MOVIE_VOLUME_TYPE == 3) then
write(IMAIN,*) ' movie output: potency or integral of strain'
else if(MOVIE_VOLUME_TYPE == 4) then
- write(IMAIN,*) ' movie output: divergence and curl'
+ write(IMAIN,*) ' movie output: divergence and curl'
else if(MOVIE_VOLUME_TYPE == 5) then
write(IMAIN,*) ' movie output: displacement'
else if(MOVIE_VOLUME_TYPE == 6) then
@@ -1847,14 +1847,14 @@
mask_noise(nmovie_points), &
noise_surface_movie(NDIM,NGLLX,NGLLY,NSPEC2D_TOP(IREGION_CRUST_MANTLE)),stat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error allocating noise arrays')
-
+
noise_sourcearray(:,:,:,:,:) = 0._CUSTOM_REAL
normal_x_noise(:) = 0._CUSTOM_REAL
normal_y_noise(:) = 0._CUSTOM_REAL
normal_z_noise(:) = 0._CUSTOM_REAL
mask_noise(:) = 0._CUSTOM_REAL
noise_surface_movie(:,:,:,:) = 0._CUSTOM_REAL
-
+
call read_parameters_noise(myrank,nrec,NSTEP,nmovie_points, &
islice_selected_rec,xi_receiver,eta_receiver,gamma_receiver,nu, &
noise_sourcearray,xigll,yigll,zigll,NSPEC2D_TOP(IREGION_CRUST_MANTLE), &
@@ -4393,7 +4393,7 @@
! save/read the surface movie using the same c routine as we do for absorbing boundaries (file ID is 9)
if (NOISE_TOMOGRAPHY/=0) then
- call close_file_abs(9)
+ call close_file_abs(9)
deallocate(noise_surface_movie)
endif
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_surface.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_surface.f90 2011-04-07 22:44:10 UTC (rev 18198)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/src/specfem3D/write_movie_surface.f90 2011-04-08 03:43:58 UTC (rev 18199)
@@ -104,7 +104,7 @@
store_val_uy(ipoin) = veloc_crust_mantle(2,iglob)*scale_veloc
store_val_uz(ipoin) = veloc_crust_mantle(3,iglob)*scale_veloc
endif
-
+
enddo
enddo
More information about the CIG-COMMITS
mailing list