[cig-commits] [commit] devel, master: renamed model_crust to model_crust_2_0 to be consistent with existing model_crust_1_0 (ffd536a)
cig_noreply at geodynamics.org
cig_noreply at geodynamics.org
Thu Nov 6 08:15:25 PST 2014
Repository : https://github.com/geodynamics/specfem3d_globe
On branches: devel,master
Link : https://github.com/geodynamics/specfem3d_globe/compare/bc58e579b3b0838a0968725a076f5904845437ca...be63f20cbb6f462104e949894dbe205d2398cd7f
>---------------------------------------------------------------
commit ffd536a8ab4e52fcf0becb7f741e027800a9d77d
Author: Dimitri Komatitsch <komatitsch at lma.cnrs-mrs.fr>
Date: Sun May 11 14:32:01 2014 +0200
renamed model_crust to model_crust_2_0 to be consistent with existing model_crust_1_0
>---------------------------------------------------------------
ffd536a8ab4e52fcf0becb7f741e027800a9d77d
setup/constants.h.in | 28 ++++++++++++-------
src/meshfem3D/meshfem3D_models.f90 | 2 +-
src/meshfem3D/model_crust_1_0.f90 | 31 ++++++----------------
.../{model_crust.f90 => model_crust_2_0.f90} | 27 +++++++++----------
src/meshfem3D/model_crustmaps.f90 | 5 +++-
src/meshfem3D/rules.mk | 4 +--
6 files changed, 46 insertions(+), 51 deletions(-)
diff --git a/setup/constants.h.in b/setup/constants.h.in
index b6dd72a..ac1b88f 100644
--- a/setup/constants.h.in
+++ b/setup/constants.h.in
@@ -163,7 +163,7 @@
! the code will print an error message and stop if it finds that the topography input file
! contains values outside this range.
! we take a safety margin just in case of a smoothed or modified model, which can locally create slightly different values
- integer, parameter :: TOPO_MINIMUM = - 11500 ! (max depth in m, Mariana trench)
+ integer, parameter :: TOPO_MINIMUM = - 11200 ! (max depth in m, Mariana trench)
integer, parameter :: TOPO_MAXIMUM = + 9000 ! (height in m, Mount Everest)
! maximum depth of the oceans in trenches and height of topo in mountains
@@ -190,8 +190,16 @@
integer, parameter :: ICRUST_CRUSTMAPS = 3
integer, parameter :: ICRUST_EPCRUST = 4
-! increase smoothing for critical regions (increases mesh stability)
- logical, parameter :: SMOOTH_CRUST = .true.
+!-------------------------------
+! crustal model cap smoothing - extension range in degree
+! note: using a smaller range, e.g. 0.5 degrees, leads to undefined Jacobian error at different places.
+! this is probably due to stretching elements below sharp gradients, especially with deep moho values.
+! so far, the only thing that works is to smooth out values and take special care of the Andes...
+! TODO: one could try to adapt this degree range to the simulation resolution in the future
+ double precision, parameter :: CAP_SMOOTHING_DEGREE_DEFAULT = 1.0d0
+
+! increase smoothing for critical regions (for instance high mountains in Europe and in the Andes) to increases mesh stability
+ logical, parameter :: SMOOTH_CRUST_EVEN_MORE = .true.
! use sedimentary layers in crustal model
logical, parameter :: INCLUDE_SEDIMENTS_CRUST = .true.
@@ -199,14 +207,14 @@
! default crustal model
! (used as default when CRUSTAL flag is set for simulation)
-!!-- uncomment for using Crust1.0
-! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUST1
-!!-- uncomment for using Crust2.0
+!-- uncomment for using Crust1.0
+! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUST1
+!-- uncomment for using Crust2.0
integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUST2
-!!-- uncomment for using General Crustmaps instead
-! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUSTMAPS
-!!-- uncomment for using EPcrust instead (European regional model)
-! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_EPCRUST
+!-- uncomment for using General Crustmaps instead
+! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_CRUSTMAPS
+!-- uncomment for using EPcrust instead (European regional model)
+! integer, parameter :: ITYPE_CRUSTAL_MODEL = ICRUST_EPCRUST
!!-----------------------------------------------------------
diff --git a/src/meshfem3D/meshfem3D_models.f90 b/src/meshfem3D/meshfem3D_models.f90
index 087e895..60cb170 100644
--- a/src/meshfem3D/meshfem3D_models.f90
+++ b/src/meshfem3D/meshfem3D_models.f90
@@ -200,7 +200,7 @@
case (ICRUST_CRUST2)
! default
! crust 2.0
- call model_crust_broadcast(myrank)
+ call model_crust_2_0_broadcast(myrank)
case (ICRUST_CRUSTMAPS)
! general crustmaps
diff --git a/src/meshfem3D/model_crust_1_0.f90 b/src/meshfem3D/model_crust_1_0.f90
index bd59928..a3e4c54 100644
--- a/src/meshfem3D/model_crust_1_0.f90
+++ b/src/meshfem3D/model_crust_1_0.f90
@@ -265,8 +265,7 @@
write(IMAIN,*)
! allocates temporary array
- allocate(bnd(CRUST1_NP,CRUST1_NLA,CRUST1_NLO), &
- stat=ier)
+ allocate(bnd(CRUST1_NP,CRUST1_NLA,CRUST1_NLO),stat=ier)
if( ier /= 0 ) call exit_MPI(0,'error allocating crustal arrays in read routine')
! initializes
@@ -372,8 +371,7 @@
close(77)
! checks min/max
- if(h_moho_min == HUGEVAL .or. h_moho_max == -HUGEVAL) &
- stop 'incorrect moho depths in read_crust_1_0_model'
+ if(h_moho_min == HUGEVAL .or. h_moho_max == -HUGEVAL) stop 'incorrect moho depths in read_crust_1_0_model'
! debug: file output for smoothed data
open(77,file='tmp-crust1.0-smooth.dat',status='unknown')
@@ -401,8 +399,7 @@
close(77)
! checks min/max
- if(h_moho_min == HUGEVAL .or. h_moho_max == -HUGEVAL) &
- stop 'incorrect moho depths in read_crust_1_0_model'
+ if(h_moho_min == HUGEVAL .or. h_moho_max == -HUGEVAL) stop 'incorrect moho depths in read_crust_1_0_model'
! frees memory
deallocate(ths,thc)
@@ -435,15 +432,6 @@
double precision :: lat,lon
double precision,dimension(CRUST1_NP) :: rho,thick,velp,vels
-
- !-------------------------------
- ! CAP smoothing - extension range in degree
- ! note: using a smaller range, e.g. 0.5 degrees, leads to undefined Jacobian error at different places.
- ! this is probably due to stretching elements below sharp gradients, especially with deep moho values.
- ! so far, the only thing that works is to smooth out values and take special care of the Andes...
- ! TODO: one could try to adapt this degree range to the simulation resolution in future
- double precision,parameter :: CAP_DEFAULT_DEGREE = 1.0d0
-
! work-around to avoid Jacobian problems when stretching mesh elements;
! one could also try to slightly change the shape of the doubling element bricks (which cause the problem)...
!
@@ -469,18 +457,17 @@
stop 'error in latitude/longitude range in crust1.0'
endif
- ! makes sure lat/lon are within range
+ ! makes sure lat/lon are within crust1.0 range
if(lat==90.0d0) lat=89.9999d0
if(lat==-90.0d0) lat=-89.9999d0
if(lon==180.0d0) lon=179.9999d0
if(lon==-180.0d0) lon=-179.9999d0
- ! sets up smoothing points
- ! by default uses CAP smoothing with 1 degree
- cap_degree = CAP_DEFAULT_DEGREE
+ ! sets up smoothing points based on cap smoothing
+ cap_degree = CAP_SMOOTHING_DEGREE_DEFAULT
! checks if inside/outside of critical region for mesh stretching
- if( SMOOTH_CRUST ) then
+ if( SMOOTH_CRUST_EVEN_MORE ) then
dist = dsqrt( (lon-LON_CRITICAL_ANDES)**2 + (lat-LAT_CRITICAL_ANDES )**2 )
if( dist < CRITICAL_RANGE ) then
! increases cap smoothing degree
@@ -493,8 +480,7 @@
endif
endif
- ! gets smoothing points and weights
- ! (see routine in model_crust.f90)
+ ! gets smoothing points and weights (see routine in model_crust_2_0.f90)
call CAP_vardegree(lon,lat,xlon,xlat,weight,cap_degree,NTHETA,NPHI)
! initializes
@@ -590,4 +576,3 @@
end subroutine get_crust_1_0_structure
-
diff --git a/src/meshfem3D/model_crust.f90 b/src/meshfem3D/model_crust_2_0.f90
similarity index 97%
rename from src/meshfem3D/model_crust.f90
rename to src/meshfem3D/model_crust_2_0.f90
index 7e9843e..6c85721 100644
--- a/src/meshfem3D/model_crust.f90
+++ b/src/meshfem3D/model_crust_2_0.f90
@@ -46,7 +46,7 @@
! reads and smooths crust2.0 model
!--------------------------------------------------------------------------------------------------
- module model_crust_par
+ module model_crust_2_0_par
! crustal_model_constants
! crustal model parameters for crust2.0
@@ -59,18 +59,18 @@
character(len=2) :: abbreviation(NCAP_CRUST/2,NCAP_CRUST)
character(len=2) :: code(NKEYS_CRUST)
- end module model_crust_par
+ end module model_crust_2_0_par
!
!--------------------------------------------------------------------------------------------------
!
- subroutine model_crust_broadcast(myrank)
+ subroutine model_crust_2_0_broadcast(myrank)
! standard routine to setup model
use constants
- use model_crust_par
+ use model_crust_2_0_par
implicit none
@@ -85,7 +85,7 @@
stat=ier)
if( ier /= 0 ) call exit_MPI(myrank,'error allocating crustal arrays')
- ! the variables read are declared and stored in structure model_crust_par
+ ! the variables read are declared and stored in structure model_crust_2_0_par
if(myrank == 0) call read_crust_model()
! broadcast the information read on the master to the nodes
@@ -97,7 +97,7 @@
call bcast_all_ch_array2(abbreviation,NCAP_CRUST/2,NCAP_CRUST,2)
call bcast_all_ch_array(code,NKEYS_CRUST,2)
- end subroutine model_crust_broadcast
+ end subroutine model_crust_2_0_broadcast
!
!-------------------------------------------------------------------------------------------------
@@ -106,7 +106,7 @@
subroutine model_crust(lat,lon,x,vp,vs,rho,moho,found_crust,elem_in_crust)
use constants
- use model_crust_par
+ use model_crust_2_0_par
implicit none
@@ -204,7 +204,7 @@
subroutine read_crust_model()
use constants
- use model_crust_par
+ use model_crust_2_0_par
implicit none
@@ -276,7 +276,7 @@
! The cap is rotated to the North Pole.
use constants
- use model_crust_par,only: NLAYERS_CRUST,NKEYS_CRUST,NCAP_CRUST
+ use model_crust_2_0_par,only: NLAYERS_CRUST,NKEYS_CRUST,NCAP_CRUST
implicit none
@@ -334,12 +334,11 @@
if(lon==180.0d0) lon=179.9999d0
if(lon==-180.0d0) lon=-179.9999d0
- ! sets up smoothing points
- ! by default uses CAP smoothing with 1 degree
- cap_degree = 1.0d0
+ ! sets up smoothing points based on cap smoothing
+ cap_degree = CAP_SMOOTHING_DEGREE_DEFAULT
! checks if inside/outside of critical region for mesh stretching
- if( SMOOTH_CRUST ) then
+ if( SMOOTH_CRUST_EVEN_MORE ) then
dist = dsqrt( (lon-LON_CRITICAL_ANDES)**2 + (lat-LAT_CRITICAL_ANDES )**2 )
if( dist < CRITICAL_RANGE ) then
! increases cap smoothing degree
@@ -449,7 +448,7 @@
subroutine get_crust_structure(ikey,vptyp,vstyp,rhtyp,thtp,thlr,velocp,velocs,dens)
- use model_crust_par,only: NLAYERS_CRUST,NKEYS_CRUST
+ use model_crust_2_0_par,only: NLAYERS_CRUST,NKEYS_CRUST
implicit none
diff --git a/src/meshfem3D/model_crustmaps.f90 b/src/meshfem3D/model_crustmaps.f90
index 4d9dd2b..a917177 100644
--- a/src/meshfem3D/model_crustmaps.f90
+++ b/src/meshfem3D/model_crustmaps.f90
@@ -369,7 +369,8 @@
num_points = 1
! checks if inside/outside of critical region for mesh stretching
- if( SMOOTH_CRUST ) then
+ if( SMOOTH_CRUST_EVEN_MORE ) then
+
dist = dsqrt( (lon-LON_CRITICAL_EUROPE)**2 + (lat-LAT_CRITICAL_EUROPE )**2 )
if( dist < CRITICAL_RANGE_EUROPE ) then
! sets up smoothing points
@@ -388,6 +389,7 @@
call CAP_vardegree(lon,lat,xlon,xlat,weight,cap_degree,NTHETA,NPHI)
num_points = NTHETA*NPHI
endif
+
dist = dsqrt( (lon-LON_CRITICAL_ANDES)**2 + (lat-LAT_CRITICAL_ANDES )**2 )
if( dist < CRITICAL_RANGE_ANDES ) then
! sets up smoothing points
@@ -406,6 +408,7 @@
call CAP_vardegree(lon,lat,xlon,xlat,weight,cap_degree,NTHETA,NPHI)
num_points = NTHETA*NPHI
endif
+
endif
! initializes
diff --git a/src/meshfem3D/rules.mk b/src/meshfem3D/rules.mk
index c3aafe5..8c3ee5d 100644
--- a/src/meshfem3D/rules.mk
+++ b/src/meshfem3D/rules.mk
@@ -85,8 +85,8 @@ meshfem3D_OBJECTS = \
$O/model_aniso_mantle.check.o \
$O/model_atten3D_QRFSI12.check.o \
$O/model_attenuation.check.o \
- $O/model_crust.check.o \
$O/model_crust_1_0.check.o \
+ $O/model_crust_2_0.check.o \
$O/model_crustmaps.check.o \
$O/model_eucrust.check.o \
$O/model_epcrust.check.o \
@@ -132,8 +132,8 @@ meshfem3D_MODULES = \
$(FC_MODDIR)/model_ak135_par.$(FC_MODEXT) \
$(FC_MODDIR)/model_aniso_mantle_par.$(FC_MODEXT) \
$(FC_MODDIR)/model_atten3d_qrfsi12_par.$(FC_MODEXT) \
- $(FC_MODDIR)/model_crust_par.$(FC_MODEXT) \
$(FC_MODDIR)/model_crust_1_0_par.$(FC_MODEXT) \
+ $(FC_MODDIR)/model_crust_2_0_par.$(FC_MODEXT) \
$(FC_MODDIR)/model_crustmaps_par.$(FC_MODEXT) \
$(FC_MODDIR)/model_epcrust_par.$(FC_MODEXT) \
$(FC_MODDIR)/model_eucrust_par.$(FC_MODEXT) \
More information about the CIG-COMMITS
mailing list