[cig-commits] r12596 - in seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta: . DATA OUTPUT_FILES src
dkomati1 at geodynamics.org
dkomati1 at geodynamics.org
Fri Aug 8 13:42:50 PDT 2008
Author: dkomati1
Date: 2008-08-08 13:42:50 -0700 (Fri, 08 Aug 2008)
New Revision: 12596
Removed:
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_jacobian_discontinuities.f90
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/Makefile
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/declarations_mesher.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90
Log:
suppressed get_jacobian_discontinuities(), which was used only for sensitivity kernel calculations
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/DATA/Par_file 2008-08-08 20:42:50 UTC (rev 12596)
@@ -31,7 +31,8 @@
# fully 3D models:
# transversely_isotropic_prem_plus_3D_crust_2.0, 3D_anisotropic, 3D_attenuation,
# s20rts, s362ani, s362iso, s362wmani, s362ani_prem, s29ea, s29ea,sea99_jp3d1994,sea99,jp3d1994
-MODEL = 1D_transversely_isotropic_prem_onecrust
+#MODEL = 1D_transversely_isotropic_prem_onecrust
+MODEL = 1D_isotropic_prem_onecrust
# parameters describing the Earth model
OCEANS = .false.
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/Makefile
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/Makefile 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/Makefile 2008-08-08 20:42:50 UTC (rev 12596)
@@ -44,7 +44,7 @@
#
#FC = gfortran
#MPIFC = /opt/mpich2_gfortran/bin/mpif90
-#FLAGS_NO_CHECK = -std=f2003 -fimplicit-none -frange-check -O3 -fmax-errors=10 -pedantic -pedantic-errors -Waliasing -Wampersand -Wcharacter-truncation -Wline-truncation -Wsurprising -Wno-tabs -Wunderflow -fno-trapping-math # -fbounds-check
+#FLAGS_NO_CHECK = -std=f2003 -fimplicit-none -frange-check -O3 -fmax-errors=10 -pedantic -pedantic-errors -Waliasing -Wampersand -Wcharacter-truncation -Wline-truncation -Wsurprising -Wno-tabs -Wunderflow -fno-trapping-math -fbounds-check
#
# Portland pgf90
@@ -112,7 +112,6 @@
$O/get_event_info.o \
$O/get_global.o \
$O/get_jacobian_boundaries.o \
- $O/get_jacobian_discontinuities.o \
$O/get_model.o \
$O/get_shape2D.o \
$O/get_shape3D.o \
@@ -308,9 +307,6 @@
$O/get_jacobian_boundaries.o: $(SPECINC)/constants.h $S/get_jacobian_boundaries.f90
${FCCOMPILE_CHECK} -c -o $O/get_jacobian_boundaries.o ${FCFLAGS_f90} $S/get_jacobian_boundaries.f90
-$O/get_jacobian_discontinuities.o: $(SPECINC)/constants.h $S/get_jacobian_discontinuities.f90
- ${FCCOMPILE_CHECK} -c -o $O/get_jacobian_discontinuities.o ${FCFLAGS_f90} $S/get_jacobian_discontinuities.f90
-
$O/get_MPI_cutplanes_xi.o: $(SPECINC)/constants.h $S/get_MPI_cutplanes_xi.f90
${FCCOMPILE_CHECK} -c -o $O/get_MPI_cutplanes_xi.o ${FCFLAGS_f90} $S/get_MPI_cutplanes_xi.f90
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/OUTPUT_FILES/values_from_mesher.h 2008-08-08 20:42:50 UTC (rev 12596)
@@ -75,7 +75,7 @@
! approximate static memory needed by the solver:
! ----------------------------------------------
!
- ! size of static arrays per slice = 7.561429217457771E-002 GB
+ ! size of static arrays per slice = 7.239703834056854E-002 GB
!
! (should be below and typically equal to 80% of 1.5 GB = 1.2 GB on pangu
! at Caltech, and below and typically equal to 85% of 2 GB = 1.7 GB
@@ -83,8 +83,8 @@
! (if significantly more, the job will not run by lack of memory)
! (if significantly less, you waste a significant amount of memory)
!
- ! size of static arrays for all slices = 0.302457168698311 GB
- ! = 2.953683288069442E-004 TB
+ ! size of static arrays for all slices = 0.289588153362274 GB
+ ! = 2.828009310178459E-004 TB
!
integer, parameter :: NEX_XI_VAL = 64
@@ -101,7 +101,7 @@
integer, parameter :: NSPECMAX_ANISO_IC = 1
integer, parameter :: NSPECMAX_ISO_MANTLE = 7616
- integer, parameter :: NSPECMAX_TISO_MANTLE = 2304
+ integer, parameter :: NSPECMAX_TISO_MANTLE = 1
integer, parameter :: NSPECMAX_ANISO_MANTLE = 1
integer, parameter :: NSPEC_CRUST_MANTLE_ATTENUAT = 1
@@ -129,7 +129,7 @@
integer, parameter :: NGLOB_CRUST_MANTLE_OCEANS = 1
- logical, parameter :: TRANSVERSE_ISOTROPY_VAL = .true.
+ logical, parameter :: TRANSVERSE_ISOTROPY_VAL = .false.
logical, parameter :: ANISOTROPIC_3D_MANTLE_VAL = .false.
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/compute_element_properties.f90 2008-08-08 20:42:50 UTC (rev 12596)
@@ -329,6 +329,7 @@
! for model density and anisotropy
integer nspec_ani
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,nspec) :: kappavstore,kappahstore,muvstore,muhstore,eta_anisostore
+
!! DK DK added this for the merged version
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ) :: kappavstore_local
!! DK DK changed this for merged version
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/create_regions_mesh.f90 2008-08-08 20:42:50 UTC (rev 12596)
@@ -67,7 +67,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
rho_vp,rho_vs, Qmu_store, tau_e_store, ibelm_moho_top, ibelm_moho_bot, ibelm_400_top, ibelm_400_bot, &
- ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670, jacobian2D_moho, jacobian2D_400, jacobian2D_670)
+ ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670)
! create the different regions of the mesh
@@ -596,12 +596,7 @@
real(kind=CUSTOM_REAL) normal_moho(NDIM,NGLLX,NGLLY,NSPEC2D_MOHO)
real(kind=CUSTOM_REAL) normal_400(NDIM,NGLLX,NGLLY,NSPEC2D_400)
real(kind=CUSTOM_REAL) normal_670(NDIM,NGLLX,NGLLY,NSPEC2D_670)
- real(kind=CUSTOM_REAL) jacobian2D_moho(NGLLX,NGLLY,NSPEC2D_MOHO)
- real(kind=CUSTOM_REAL) jacobian2D_400(NGLLX,NGLLY,NSPEC2D_400)
- real(kind=CUSTOM_REAL) jacobian2D_670(NGLLX,NGLLY,NSPEC2D_670)
- logical :: is_superbrick
-
! the height at which the central cube is cut
integer :: nz_inf_limit
@@ -918,26 +913,15 @@
include "comp_mass_matrix_one_element.f90"
include "store_xelm_yelm_zelm.f90"
-! boundary mesh
- if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
- is_superbrick=.false.
- ispec_superbrick=0
- call get_jacobian_discontinuities(myrank,ispec,ix_elem,iy_elem,rmin,rmax,r1,r2,r3,r4,r5,r6,r7,r8, &
- xstore(:,:,:,ispec),ystore(:,:,:,ispec),zstore(:,:,:,ispec),dershape2D_bottom, &
- ibelm_moho_top,ibelm_moho_bot,ibelm_400_top,ibelm_400_bot,ibelm_670_top,ibelm_670_bot, &
- normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
- ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot, &
- NSPEC2D_MOHO,NSPEC2D_400,NSPEC2D_670,r_moho,r_400,r_670, &
- is_superbrick,USE_ONE_LAYER_SB,ispec_superbrick,nex_eta_moho,HONOR_1D_SPHERICAL_MOHO)
- endif
-
! end of loop on all the regular elements
enddo
enddo
enddo
+
!----
!---- mesh doubling elements
!----
+
! If there is a doubling at the top of this region, let us add these elements.
! The superbrick implements a symmetric four-to-two doubling and therefore replaces
! a basic regular block of 2 x 2 = 4 elements.
@@ -1125,18 +1109,6 @@
include "comp_mass_matrix_one_element.f90"
include "store_xelm_yelm_zelm.f90"
-! boundary mesh
- if (ipass == 2 .and. SAVE_BOUNDARY_MESH .and. iregion_code == IREGION_CRUST_MANTLE) then
- is_superbrick=.true.
- call get_jacobian_discontinuities(myrank,ispec,ix_elem,iy_elem,rmin,rmax,r1,r2,r3,r4,r5,r6,r7,r8, &
- xstore(:,:,:,ispec),ystore(:,:,:,ispec),zstore(:,:,:,ispec),dershape2D_bottom, &
- ibelm_moho_top,ibelm_moho_bot,ibelm_400_top,ibelm_400_bot,ibelm_670_top,ibelm_670_bot, &
- normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
- ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot, &
- NSPEC2D_MOHO,NSPEC2D_400,NSPEC2D_670,r_moho,r_400,r_670, &
- is_superbrick,USE_ONE_LAYER_SB,ispec_superbrick,nex_eta_moho,HONOR_1D_SPHERICAL_MOHO)
- endif
-
! end of loops on the mesh doubling elements
enddo
enddo
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/declarations_mesher.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/declarations_mesher.f90 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/declarations_mesher.f90 2008-08-08 20:42:50 UTC (rev 12596)
@@ -57,9 +57,6 @@
real(kind=CUSTOM_REAL) normal_moho(NDIM,NGLLX,NGLLY,NSPEC2D_MOHO)
real(kind=CUSTOM_REAL) normal_400(NDIM,NGLLX,NGLLY,NSPEC2D_400)
real(kind=CUSTOM_REAL) normal_670(NDIM,NGLLX,NGLLY,NSPEC2D_670)
- real(kind=CUSTOM_REAL) jacobian2D_moho(NGLLX,NGLLY,NSPEC2D_MOHO)
- real(kind=CUSTOM_REAL) jacobian2D_400(NGLLX,NGLLY,NSPEC2D_400)
- real(kind=CUSTOM_REAL) jacobian2D_670(NGLLX,NGLLY,NSPEC2D_670)
!!!!!!!!!!!!!!!! DK DK for merged version, all the arrays below are allocated statically instead
!!!!!!!!!!!!!!!! DK DK for merged version, all the arrays below are allocated statically instead
Deleted: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_jacobian_discontinuities.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_jacobian_discontinuities.f90 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/get_jacobian_discontinuities.f90 2008-08-08 20:42:50 UTC (rev 12596)
@@ -1,203 +0,0 @@
-!=====================================================================
-!
-! S p e c f e m 3 D G l o b e V e r s i o n 4 . 1
-! --------------------------------------------------
-!
-! Main authors: Dimitri Komatitsch and Jeroen Tromp
-! Seismological Laboratory, California Institute of Technology, USA
-! and University of Pau / CNRS / INRIA, France
-! (c) California Institute of Technology and University of Pau / CNRS / INRIA
-! August 2008
-!
-! This program is free software; you can redistribute it and/or modify
-! it under the terms of the GNU General Public License as published by
-! the Free Software Foundation; either version 2 of the License, or
-! (at your option) any later version.
-!
-! This program is distributed in the hope that it will be useful,
-! but WITHOUT ANY WARRANTY; without even the implied warranty of
-! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-! GNU General Public License for more details.
-!
-! You should have received a copy of the GNU General Public License along
-! with this program; if not, write to the Free Software Foundation, Inc.,
-! 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
-!
-!=====================================================================
-
-subroutine get_jacobian_discontinuities(myrank,ispec,ix_elem,iy_elem,rmin,rmax,r1,r2,r3,r4,r5,r6,r7,r8, &
- xstore,ystore,zstore,dershape2D_bottom, &
- ibelm_moho_top,ibelm_moho_bot,ibelm_400_top,ibelm_400_bot,ibelm_670_top,ibelm_670_bot, &
- normal_moho,normal_400,normal_670,jacobian2D_moho,jacobian2D_400,jacobian2D_670, &
- ispec2D_moho_top,ispec2D_moho_bot,ispec2D_400_top,ispec2D_400_bot,ispec2D_670_top,ispec2D_670_bot, &
- NSPEC2D_MOHO,NSPEC2D_400,NSPEC2D_670,r_moho,r_400,r_670, &
- is_superbrick,USE_ONE_LAYER_SB,ispec_superbrick,nex_eta_moho,HONOR_1D_SPHERICAL_MOHO)
-
- implicit none
-
- include 'constants.h'
-
- ! input
- integer myrank, ispec, ix_elem, iy_elem
- double precision rmin,rmax
- double precision xstore(NGLLX,NGLLY,NGLLZ)
- double precision ystore(NGLLX,NGLLY,NGLLZ)
- double precision zstore(NGLLX,NGLLY,NGLLZ)
- double precision dershape2D_bottom(NDIM2D,NGNOD2D,NGLLX,NGLLY)
- integer NSPEC2D_MOHO, NSPEC2D_400, NSPEC2D_670, nex_eta_moho, ispec_superbrick
- double precision r_moho, r_400, r_670
- logical :: is_superbrick, USE_ONE_LAYER_SB,HONOR_1D_SPHERICAL_MOHO
-
- ! output
- integer ispec2D_moho_top, ispec2D_moho_bot, ispec2D_400_top, ispec2D_400_bot, ispec2D_670_top, ispec2D_670_bot
- integer,dimension(NSPEC2D_MOHO) :: ibelm_moho_top, ibelm_moho_bot
- integer,dimension(NSPEC2D_400) :: ibelm_400_top, ibelm_400_bot
- integer,dimension(NSPEC2D_670) :: ibelm_670_top, ibelm_670_bot
- real(kind=CUSTOM_REAL) :: normal_moho(NDIM,NGLLX,NGLLY,NSPEC2D_MOHO), jacobian2D_moho(NGLLX,NGLLY,NSPEC2D_MOHO)
- real(kind=CUSTOM_REAL) :: normal_400(NDIM,NGLLX,NGLLY,NSPEC2D_400), jacobian2D_400(NGLLX,NGLLY,NSPEC2D_400)
- real(kind=CUSTOM_REAL) :: normal_670(NDIM,NGLLX,NGLLY,NSPEC2D_670), jacobian2D_670(NGLLX,NGLLY,NSPEC2D_670)
-
- ! local variables
- double precision, dimension(NGNOD2D) :: xelm2, yelm2, zelm2
- double precision :: r1, r2, r3, r4, r5, r6, r7, r8
- double precision :: target_moho_high, target_moho_low, target_400_high, target_400_low, target_670_high, target_670_low
- integer :: nele_sub_block, ispec_list(16), map_irem_ix_12(8), map_irem_ix_34(8), map_irem_iy_odd(8), map_irem_iy_even(8)
- integer :: map_isub_ix(4), map_isub_iy(4), map_ix(NSPEC_DOUBLING_SUPERBRICK), map_iy(NSPEC_DOUBLING_SUPERBRICK)
- integer :: i, ispec_superbrick_current, isub_block, irem_block, irem_ix, irem_iy, ix,iy,ix_top,iy_top, ispec2D_moho_bot_map
-
- ! ======================
-
-
- ! find the coordinates of 9 nodes for the bottom surface element to compute the jacobian if needed
- xelm2(1)=xstore(1,1,1)
- yelm2(1)=ystore(1,1,1)
- zelm2(1)=zstore(1,1,1)
- xelm2(2)=xstore(NGLLX,1,1)
- yelm2(2)=ystore(NGLLX,1,1)
- zelm2(2)=zstore(NGLLX,1,1)
- xelm2(3)=xstore(NGLLX,NGLLY,1)
- yelm2(3)=ystore(NGLLX,NGLLY,1)
- zelm2(3)=zstore(NGLLX,NGLLY,1)
- xelm2(4)=xstore(1,NGLLY,1)
- yelm2(4)=ystore(1,NGLLY,1)
- zelm2(4)=zstore(1,NGLLY,1)
- xelm2(5)=xstore((NGLLX+1)/2,1,1)
- yelm2(5)=ystore((NGLLX+1)/2,1,1)
- zelm2(5)=zstore((NGLLX+1)/2,1,1)
- xelm2(6)=xstore(NGLLX,(NGLLY+1)/2,1)
- yelm2(6)=ystore(NGLLX,(NGLLY+1)/2,1)
- zelm2(6)=zstore(NGLLX,(NGLLY+1)/2,1)
- xelm2(7)=xstore((NGLLX+1)/2,NGLLY,1)
- yelm2(7)=ystore((NGLLX+1)/2,NGLLY,1)
- zelm2(7)=zstore((NGLLX+1)/2,NGLLY,1)
- xelm2(8)=xstore(1,(NGLLY+1)/2,1)
- yelm2(8)=ystore(1,(NGLLY+1)/2,1)
- zelm2(8)=zstore(1,(NGLLY+1)/2,1)
- xelm2(9)=xstore((NGLLX+1)/2,(NGLLY+1)/2,1)
- yelm2(9)=ystore((NGLLX+1)/2,(NGLLY+1)/2,1)
- zelm2(9)=zstore((NGLLX+1)/2,(NGLLY+1)/2,1)
-
-! radii to determine if an element is on the discontinuity or not
- target_moho_high = r_moho * (ONE + SMALLVAL)
- target_moho_low = r_moho * (ONE - SMALLVAL)
- target_400_high = r_400 * (ONE + SMALLVAL)
- target_400_low = r_400 * (ONE - SMALLVAL)
- target_670_high = r_670 * (ONE + SMALLVAL)
- target_670_low = r_670 * (ONE - SMALLVAL)
-
-! setup the mapping array for superbrick case (only invoked for Moho bottom)
- if (is_superbrick) then
- map_irem_ix_12=(/2,2,0,1,0,1,0,0/)
- map_irem_ix_34=(/1,1,0,2,0,2,0,0/)
- map_irem_iy_odd=(/1,2,0,1,0,2,0,0/)
- map_irem_iy_even=(/2,1,0,2,0,1,0,0/)
- if (USE_ONE_LAYER_SB) then
- nele_sub_block = 7
- ispec_list=(/1,2,4,6,8,9,11,13,15,16,18,20,22,23,25,27/)
- else
- nele_sub_block = 8
- ispec_list=(/1,2,4,6,9,10,12,14,17,18,20,22,25,26,28,30/)
- endif
- map_isub_ix=(/2,2,1,1/)
- map_isub_iy=(/2,1,2,1/)
-
- map_ix(1:NSPEC_DOUBLING_SUPERBRICK) = 0
- map_iy(1:NSPEC_DOUBLING_SUPERBRICK) = 0
-
- do i = 1, 16
- ispec_superbrick_current=ispec_list(i)
- isub_block = ispec_superbrick_current/nele_sub_block + 1
- irem_block = mod(ispec_superbrick_current,nele_sub_block)
-
- if (isub_block > 2) then
- irem_ix = map_irem_ix_34(irem_block)
- else
- irem_ix = map_irem_ix_12(irem_block)
- endif
- if (mod(isub_block,2) == 0) then
- irem_iy = map_irem_iy_even(irem_block)
- else
- irem_iy = map_irem_iy_odd(irem_block)
- endif
- map_ix(ispec_list(i)) = (map_isub_ix(isub_block) - 1) * 2 + irem_ix
- map_iy(ispec_list(i)) = (map_isub_iy(isub_block) - 1) * 2 + irem_iy
- enddo
- endif
-
-! determine if the elements are on the discontinuity, and calculate the boundary jaocobian if needed
- if (.not. is_superbrick) then
-
-! Moho top
- if (.not. SUPPRESS_CRUSTAL_MESH .and. HONOR_1D_SPHERICAL_MOHO .and. &
- abs(rmin-r_moho)/r_moho < SMALLVAL .and. r1 < target_moho_high .and. r2 < target_moho_high &
- .and. r3 < target_moho_high .and. r4 < target_moho_high) then
- ispec2D_moho_top = ispec2D_moho_top + 1
- ibelm_moho_top(ispec2D_moho_top) = ispec
- call compute_jacobian_2D(myrank,ispec2D_moho_top,xelm2,yelm2,zelm2,dershape2D_bottom, &
- jacobian2D_moho,normal_moho,NGLLX,NGLLY,NSPEC2D_MOHO)
-! 400 top
- else if (abs(rmin-r_400)/r_400 < SMALLVAL .and. r1 < target_400_high .and. r2 < target_400_high &
- .and. r3 < target_400_high .and. r4 < target_400_high) then
- ispec2D_400_top = ispec2D_400_top + 1
- ibelm_400_top(ispec2D_400_top) = ispec
- call compute_jacobian_2D(myrank,ispec2D_400_top,xelm2,yelm2,zelm2,dershape2D_bottom, &
- jacobian2D_400,normal_400,NGLLX,NGLLY,NSPEC2D_400)
-
-! 400 bot
- else if (abs(rmax-r_400)/r_400 < SMALLVAL .and. r5 > target_400_low .and. r6 > target_400_low &
- .and. r7 > target_400_low .and. r8 > target_400_low) then
- ispec2D_400_bot = ispec2D_400_bot + 1
- ibelm_400_bot(ispec2D_400_bot) = ispec
-
-! 670 top
- else if (abs(rmin-r_670)/r_670 < SMALLVAL .and. r1 < target_670_high .and. r2 < target_670_high &
- .and. r3 < target_670_high .and. r4 < target_670_high) then
- ispec2D_670_top = ispec2D_670_top + 1
- ibelm_670_top(ispec2D_670_top) = ispec
- call compute_jacobian_2D(myrank,ispec2D_670_top,xelm2,yelm2,zelm2,dershape2D_bottom, &
- jacobian2D_670,normal_670,NGLLX,NGLLY,NSPEC2D_670)
-! 670 bot
- else if (abs(rmax-r_670)/r_670 < SMALLVAL .and. r5 > target_670_low .and. r6 > target_670_low &
- .and. r7 > target_670_low .and. r8 > target_670_low) then
- ispec2D_670_bot = ispec2D_670_bot + 1
- ibelm_670_bot(ispec2D_670_bot) = ispec
- endif
-
- else ! superbrick case
- ! Moho bot (special care should be taken to deal with mapping 2D element indices)
- if (.not. SUPPRESS_CRUSTAL_MESH .and. HONOR_1D_SPHERICAL_MOHO .and. &
- abs(rmax-r_moho)/r_moho < SMALLVAL .and. r5 > target_moho_low .and. r6 > target_moho_low &
- .and. r7 > target_moho_low .and. r8 > target_moho_low) then
- ispec2D_moho_bot = ispec2D_moho_bot + 1
- ix=map_ix(ispec_superbrick)
- iy=map_iy(ispec_superbrick)
- if (ix == 0 .or. iy == 0) call exit_mpi(myrank, 'Check (ix,iy) on the Moho bot is 0')
- ix_top = (ix_elem - 1) + ix
- iy_top = (iy_elem - 1) + iy
- ispec2D_moho_bot_map = (ix_top - 1) * nex_eta_moho + iy_top
- ibelm_moho_bot(ispec2D_moho_bot_map) = ispec
- endif
- endif
-
-end subroutine get_jacobian_discontinuities
-
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/meshfem3D.f90 2008-08-08 20:42:50 UTC (rev 12596)
@@ -1225,7 +1225,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
rho_vp,rho_vs, Qmu_store, tau_e_store, ibelm_moho_top, ibelm_moho_bot, ibelm_400_top, ibelm_400_bot, &
- ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670, jacobian2D_moho, jacobian2D_400, jacobian2D_670)
+ ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670)
else if(iregion_code == IREGION_OUTER_CORE) then
! outer_core
@@ -1271,7 +1271,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
rho_vp,rho_vs, Qmu_store, tau_e_store, ibelm_moho_top, ibelm_moho_bot, ibelm_400_top, ibelm_400_bot, &
- ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670, jacobian2D_moho, jacobian2D_400, jacobian2D_670)
+ ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670)
else if(iregion_code == IREGION_INNER_CORE) then
! inner_core
@@ -1316,7 +1316,7 @@
c36store,c44store,c45store,c46store,c55store,c56store,c66store, &
iboun, locval, ifseg, xp,yp,zp, rmass_ocean_load, mask_ibool, copy_ibool_ori, iMPIcut_xi,iMPIcut_eta, &
rho_vp,rho_vs, Qmu_store, tau_e_store, ibelm_moho_top, ibelm_moho_bot, ibelm_400_top, ibelm_400_bot, &
- ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670, jacobian2D_moho, jacobian2D_400, jacobian2D_670)
+ ibelm_670_top, ibelm_670_bot, normal_moho, normal_400, normal_670)
else
stop 'DK DK incorrect region in merged code'
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90 2008-08-08 19:31:31 UTC (rev 12595)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/version41_beta/src/specfem3D.f90 2008-08-08 20:42:50 UTC (rev 12596)
@@ -181,7 +181,6 @@
kappahstore_crust_mantle,muhstore_crust_mantle,eta_anisostore_crust_mantle
! arrays for full anisotropy only when needed
- integer nspec_iso,nspec_tiso,nspec_ani
real(kind=CUSTOM_REAL), dimension(NGLLX,NGLLY,NGLLZ,NSPECMAX_ANISO_MANTLE) :: &
c11store_crust_mantle,c12store_crust_mantle,c13store_crust_mantle, &
c22store_crust_mantle,c23store_crust_mantle,c33store_crust_mantle, &
@@ -426,9 +425,6 @@
double precision radius,rho,drhodr,vp,vs,Qkappa,Qmu
double precision, dimension(NRAD_GRAVITY) :: d_ln_density_dr_table
-! flags to read kappa and mu and anisotropy arrays in regions where needed
- logical READ_KAPPA_MU,READ_TISO
-
! names of the data files for all the processors in MPI
character(len=150) outputname
@@ -879,50 +875,6 @@
! total number of messages corresponding to these common faces
NUMMSGS_FACES = NPROC_ONE_DIRECTION*NUM_FACES*NUM_MSG_TYPES
-! start reading the databases
-
-! read arrays created by the mesher
-
-! crust and mantle
-
- if(ANISOTROPIC_3D_MANTLE_VAL) then
- READ_KAPPA_MU = .false.
- READ_TISO = .false.
- nspec_iso = 1
- nspec_tiso = 1
- nspec_ani = NSPEC_CRUST_MANTLE
- else
- nspec_iso = NSPEC_CRUST_MANTLE
- if(TRANSVERSE_ISOTROPY_VAL) then
- nspec_tiso = NSPECMAX_TISO_MANTLE
- else
- nspec_tiso = 1
- endif
- nspec_ani = 1
- READ_KAPPA_MU = .true.
- READ_TISO = .true.
- endif
-
-! outer core (no anisotropy nor S velocity)
-! rmass_ocean_load is not used in this routine because it is meaningless in the outer core
- READ_KAPPA_MU = .false.
- READ_TISO = .false.
- nspec_iso = NSPEC_OUTER_CORE
- nspec_tiso = 1
- nspec_ani = 1
-
-! inner core (no anisotropy)
-! rmass_ocean_load is not used in this routine because it is meaningless in the inner core
- READ_KAPPA_MU = .true.
- READ_TISO = .false.
- nspec_iso = NSPEC_INNER_CORE
- nspec_tiso = 1
- if(ANISOTROPIC_INNER_CORE_VAL) then
- nspec_ani = NSPEC_INNER_CORE
- else
- nspec_ani = 1
- endif
-
! check that the number of points in this slice is correct
if(minval(ibool_crust_mantle) /= 1 .or. maxval(ibool_crust_mantle) /= NGLOB_CRUST_MANTLE) &
call exit_MPI(myrank,'incorrect global numbering: iboolmax does not equal nglob in crust and mantle')
More information about the cig-commits
mailing list