[cig-commits] [commit] devel, master: updates interpolation tool and rules.mk in src/tomography (72d3705)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Thu Nov 6 08:32:44 PST 2014
Repository : https://github.com/geodynamics/specfem3d_globe
On branches: devel,master
Link : https://github.com/geodynamics/specfem3d_globe/compare/bc58e579b3b0838a0968725a076f5904845437ca...be63f20cbb6f462104e949894dbe205d2398cd7f
>---------------------------------------------------------------
commit 72d3705bffdf469a7ddbfffae4290ad36db89064
Author: daniel peter <peterda at ethz.ch>
Date: Wed Oct 22 14:53:14 2014 +0200
updates interpolation tool and rules.mk in src/tomography
>---------------------------------------------------------------
72d3705bffdf469a7ddbfffae4290ad36db89064
src/meshfem3D/compute_element_properties.f90 | 6 +--
src/meshfem3D/moho_stretching.f90 | 4 +-
src/tomography/interpolate_model.F90 | 73 +++++++++++++++++-----------
src/tomography/interpolate_model_kdtree.f90 | 12 +++--
src/tomography/rules.mk | 6 +++
5 files changed, 64 insertions(+), 37 deletions(-)
diff --git a/src/meshfem3D/compute_element_properties.f90 b/src/meshfem3D/compute_element_properties.f90
index c0683e0..05f2fab 100644
--- a/src/meshfem3D/compute_element_properties.f90
+++ b/src/meshfem3D/compute_element_properties.f90
@@ -217,8 +217,8 @@
! adds surface topography
if (TOPOGRAPHY) then
if (idoubling(ispec) == IFLAG_CRUST .or. &
- idoubling(ispec) == IFLAG_220_80 .or. &
- idoubling(ispec) == IFLAG_80_MOHO) then
+ idoubling(ispec) == IFLAG_220_80 .or. &
+ idoubling(ispec) == IFLAG_80_MOHO) then
! stretches mesh between surface and R220 accordingly
if (USE_GLL) then
! stretches every GLL point accordingly
@@ -235,7 +235,7 @@
.or. THREE_D_MODEL == THREE_D_MODEL_S362ANI_PREM .or. THREE_D_MODEL == THREE_D_MODEL_S29EA) then
! stretching between 220 and 770
if (idoubling(ispec) == IFLAG_670_220 .or. &
- idoubling(ispec) == IFLAG_MANTLE_NORMAL) then
+ idoubling(ispec) == IFLAG_MANTLE_NORMAL) then
if (USE_GLL) then
! stretches every GLL point accordingly
call add_topography_410_650_gll(myrank,xstore,ystore,zstore,ispec,nspec)
diff --git a/src/meshfem3D/moho_stretching.f90 b/src/meshfem3D/moho_stretching.f90
index 3c43256..7ba55f8 100644
--- a/src/meshfem3D/moho_stretching.f90
+++ b/src/meshfem3D/moho_stretching.f90
@@ -26,7 +26,7 @@
!=====================================================================
subroutine moho_stretching_honor_crust(myrank,xelm,yelm,zelm, &
- elem_in_crust,elem_in_mantle)
+ elem_in_crust,elem_in_mantle)
! stretching the moho according to the crust 2.0
! input: myrank, xelm, yelm, zelm
@@ -224,7 +224,7 @@
!
subroutine moho_stretching_honor_crust_reg(myrank,xelm,yelm,zelm, &
- elem_in_crust,elem_in_mantle)
+ elem_in_crust,elem_in_mantle)
! regional routine: for REGIONAL_MOHO_MESH adaptations
!
diff --git a/src/tomography/interpolate_model.F90 b/src/tomography/interpolate_model.F90
index 0bd1a36..2633b62 100644
--- a/src/tomography/interpolate_model.F90
+++ b/src/tomography/interpolate_model.F90
@@ -44,7 +44,13 @@
!
! usage: ./xinterpolate_model vsv MODEL_M10/
!
-
+! note on mid-point search option:
+! - mid-point-search == 1: looking for mid-points only is a good approach when changing number of processes (NPROC) only,
+! while keeping the mesh resolution (NEX) the same
+! (by default set to .true.)
+! - mid-point-search == 0: searching for each single gll point is a good approach when changing resolution (NEX) of meshes;
+! in general, interpolation suffers and might lead to differences at internal interfaces (e.g. 410)
+!
!------------------------------------------------------------------------------
@@ -132,12 +138,7 @@
character(len=MAX_STRING_LEN) :: dir_topo2
character(len=MAX_STRING_LEN) :: input_model_dir,output_model_dir
- ! kdtree search:
- ! note: looking for mid-points only is a good approach when changing number of processes (NPROC) only,
- ! while keeping the mesh resolution (NEX) the same;
- ! searching for each single gll point is a good approach when changing resolution (NEX) of meshes;
- ! in general, interpolation suffers and might lead to differences at internal interfaces (e.g. 410);
- ! by default set to .true.
+ ! kdtree search
integer :: want_midpoint
logical :: USE_MIDPOINT_SEARCH
@@ -745,17 +746,26 @@
enddo
enddo
if (inodes /= kdtree_nnodes_local ) stop 'Error index inodes does not match nnodes_local'
+ call synchronize_all()
! creates kd-tree for searching
- do i = 0,NPROCTOT_VAL-1
- if (myrank == i) then
- !print*,'kd-tree setup for process: ',myrank
- call kdtree_setup()
- endif
- call synchronize_all()
- enddo
+ ! serial way
+ !do i = 0,NPROCTOT_VAL-1
+ ! if (myrank == i) then
+ ! print*,'kd-tree setup for process: ',myrank
+ ! call kdtree_setup()
+ ! endif
+ ! call synchronize_all()
+ !enddo
+ ! parallel way
+ call kdtree_setup()
endif
+ call synchronize_all()
+
+ ! user output
+ if (myrank == 0) print*,'looping over elements:'
+ call synchronize_all()
! loop over all elements (mainly those in this layer)
do ispec = 1, nspec
@@ -786,6 +796,8 @@
USE_MIDPOINT_SEARCH)
endif
enddo ! ispec
+ if (myrank == 0) print*
+ call synchronize_all()
! frees tree memory
if (.not. DO_BRUTE_FORCE_SEARCH) then
@@ -805,9 +817,8 @@
enddo
! user output
- if (myrank == 0) then
- print*
- endif
+ if (myrank == 0) print*
+ call synchronize_all()
enddo ! ilayer
! frees memory
@@ -1083,8 +1094,9 @@
stop 'Error locate_single: specifying closest ispec element'
endif
! checks minimum distance within a typical element size
- if (dist_min > typical_size ) then
- print*,'Error rank ',myrank,' - invalid dist_min: ',dist_min * R_EARTH_KM,'(km)'
+ if (dist_min > 2 * typical_size ) then
+ print*,'Warning: rank ',myrank,' - large dist_min: ',dist_min * R_EARTH_KM,'(km)', &
+ 'element size:',typical_size * R_EARTH_KM
print*,'element',ispec,'selected ispec:',ispec_selected,'in rank:',rank_selected,'iglob_min:',iglob_min
print*,'typical element size:',typical_size * 0.5 * R_EARTH_KM
print*,'target location:',xyz_target(:)
@@ -1093,7 +1105,8 @@
print*,'found radius :',sqrt(kdtree_nodes_local(1,iglob_min)**2 &
+ kdtree_nodes_local(2,iglob_min)**2 &
+ kdtree_nodes_local(3,iglob_min)**2 ) * R_EARTH_KM,'(km)'
- stop 'Error dist_min too large'
+ !debug
+ !stop 'Error dist_min too large'
endif
! debug
!if (myrank == 0 .and. iglob < 100) &
@@ -1141,8 +1154,9 @@
stop 'Error locate_single: specifying closest ispec element'
endif
! checks minimum distance within a typical element size
- if (dist_min > typical_size ) then
- print*,'Error rank ',myrank,' - invalid dist_min: ',dist_min * R_EARTH_KM,'(km)'
+ if (dist_min > 2 * typical_size ) then
+ print*,'Warning: rank ',myrank,' - large dist_min: ',dist_min * R_EARTH_KM,'(km)', &
+ 'element size:',typical_size * R_EARTH_KM
print*,'element',ispec,'selected ispec:',ispec_selected,'in rank:',rank_selected,'iglob_min:',iglob_min
print*,'typical element size:',typical_size * 0.5 * R_EARTH_KM
print*,'target location:',xyz_target(:)
@@ -1151,7 +1165,8 @@
print*,'found radius :',sqrt(kdtree_nodes_local(1,iglob_min)**2 &
+ kdtree_nodes_local(2,iglob_min)**2 &
+ kdtree_nodes_local(3,iglob_min)**2 ) * R_EARTH_KM,'(km)'
- stop 'Error dist_min too large'
+ ! debug
+ !stop 'Error dist_min too large'
endif
! debug
!if (myrank == 0 .and. iglob < 100) &
@@ -1175,14 +1190,16 @@
! checks distance
dist_min = sqrt((x_found-x_target)**2 + (y_found-y_target)**2 + (z_found-z_target)**2)
- if (dist_min > typical_size ) then
- print*,'Error selected gll point: rank ',myrank,' - invalid dist_min: ',dist_min * R_EARTH_KM,'(km)'
+ if (dist_min > 2 * typical_size ) then
+ print*,'Warning: rank ',myrank,' - large dist_min: ',dist_min * R_EARTH_KM,'(km)', &
+ 'element size:',typical_size * R_EARTH_KM
print*,'target location:',xyz_target(:)
print*,'target radius :',sqrt(xyz_target(1)**2 + xyz_target(2)**2 + xyz_target(3)**2) * R_EARTH_KM,'(km)'
print*,'gll location :',x_found,y_found,z_found
print*,'gll radius :',sqrt(x_found**2 + y_found**2 + z_found**2) * R_EARTH_KM,'(km)'
print*,'distance gll:',dist_min * R_EARTH_KM,'(km)'
- stop 'Error model value invalid'
+ ! debug
+ !stop 'Error gll model value invalid'
endif
! debug
!if (myrank == 0 .and. iglob < 100) &
@@ -1615,11 +1632,11 @@
! add warning if estimate is poor
! (usually means receiver outside the mesh given by the user)
if (DO_WARNING) then
- if (final_distance > typical_size / NGLLX * R_EARTH_KM ) then
+ if (final_distance > typical_size * R_EARTH_KM ) then
print*, '*****************************************************************'
print*, '***** WARNING: location estimate is poor *****'
print*, '*****************************************************************'
- print*, 'closest estimate found: ',sngl(final_distance),'km away',' - not within:',typical_size /NGLLX * R_EARTH_KM
+ print*, 'closest estimate found: ',sngl(final_distance),'km away',' - not within:',typical_size * R_EARTH_KM
print*, ' in rank ',rank_selected,' in element ',ispec_selected,ix_initial_guess,iy_initial_guess,iz_initial_guess
print*, ' at xi,eta,gamma coordinates = ',xi,eta,gamma
print*, ' at radius ',sqrt(x**2 + y**2 + z**2) * R_EARTH_KM,'(km)'
diff --git a/src/tomography/interpolate_model_kdtree.f90 b/src/tomography/interpolate_model_kdtree.f90
index b51a330..613a246 100644
--- a/src/tomography/interpolate_model_kdtree.f90
+++ b/src/tomography/interpolate_model_kdtree.f90
@@ -172,7 +172,8 @@ contains
if (numnodes > 2147483646 - i ) stop 'Error number of nodes might exceed integer limit'
enddo
if (be_verbose) then
- print*,' theoretical number of nodes: ',numnodes
+ print*,' theoretical number of nodes: ',numnodes
+ print*,' tree memory size: ',( numnodes * 32 )/1024./1024.,'MB'
endif
! local ordering
@@ -197,9 +198,9 @@ contains
if (.not. associated(kdtree) ) stop 'Error creation of kd-tree failed'
if (be_verbose) then
- print*,' actual number of nodes: ',numnodes
+ print*,' actual number of nodes: ',numnodes
! tree node size: 4 (idim) + 8 (cut_value) + 4 (ipoint) + 2*4 (ibound_**) + 2*4 (left,right) = 32 bytes
- print*,' tree memory size: ', ( numnodes * 32 )/1024./1024.,'MB'
+ print*,' tree memory size: ',( numnodes * 32 )/1024./1024.,'MB'
print*,' maximum depth : ',maxdepth
! timing
@@ -333,7 +334,10 @@ contains
! creates new node
allocate(node,stat=ier)
- if (ier /= 0) stop 'Error allocating kd-tree node'
+ if (ier /= 0) then
+ print*,'Error creating node: ',numnodes
+ stop 'Error allocating kd-tree node'
+ endif
! initializes new node
node%idim = -1
diff --git a/src/tomography/rules.mk b/src/tomography/rules.mk
index 285f3ae..2dd608a 100644
--- a/src/tomography/rules.mk
+++ b/src/tomography/rules.mk
@@ -112,6 +112,7 @@ xadd_model_OBJECTS = \
$(EMPTY_MACRO)
xadd_model_SHARED_OBJECTS = \
+ $O/shared_par.shared_module.o \
$O/parallel.sharedmpi.o \
$O/exit_mpi.shared.o \
$O/gll_library.shared.o \
@@ -172,6 +173,7 @@ xconvert_model_file_adios_OBJECTS = \
$(EMPTY_MACRO)
xconvert_model_file_adios_SHARED_OBJECTS = \
+ $O/shared_par.shared_module.o \
$O/parallel.sharedmpi.o \
$O/adios_helpers_definitions.shared_adios_module.o \
$O/adios_helpers_writers.shared_adios_module.o \
@@ -192,6 +194,7 @@ xinterpolate_model_OBJECTS = \
$(EMPTY_MACRO)
xinterpolate_model_SHARED_OBJECTS = \
+ $O/shared_par.shared_module.o \
$O/parallel.sharedmpi.o \
$O/gll_library.shared.o \
$O/hex_nodes.shared.o \
@@ -213,6 +216,7 @@ xsmooth_sem_OBJECTS = \
$(EMPTY_MACRO)
xsmooth_sem_SHARED_OBJECTS = \
+ $O/shared_par.shared_module.o \
$O/parallel.sharedmpi.o \
$O/exit_mpi.shared.o \
$O/get_all_eight_slices.shared.o \
@@ -233,6 +237,7 @@ xsum_kernels_OBJECTS = \
$(EMPTY_MACRO)
xsum_kernels_SHARED_OBJECTS = \
+ $O/shared_par.shared_module.o \
$O/parallel.sharedmpi.o \
$(EMPTY_MACRO)
@@ -249,6 +254,7 @@ xsum_preconditioned_kernels_OBJECTS = \
$(EMPTY_MACRO)
xsum_preconditioned_kernels_SHARED_OBJECTS = \
+ $O/shared_par.shared_module.o \
$O/parallel.sharedmpi.o \
$(EMPTY_MACRO)
More information about the CIG-COMMITS
mailing list