Skip to content

Commit 4208388

Browse files
danielpeterdanielpeter
committed
updates comments; adds moho check for Berkeley model testing
1 parent 68a8f50 commit 4208388

File tree

3 files changed

+46
-23
lines changed

3 files changed

+46
-23
lines changed

‎src/meshfem3D/model_crust_berkeley.f90‎

Lines changed: 13 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -148,7 +148,7 @@ subroutine model_berkeley_crust(x,theta,phi,vp,vs,rho,moho,found_crust,elem_in_c
148148

149149
depth = (1.d0 - x) * EARTH_R_KM
150150

151-
call get_crust_val_csem(theta,phi,depth,rho,vp,vsv,vsh,moho_depth)
151+
call get_crust_val_csem(theta,phi,depth,rho,vp,vsv,vsh,moho_depth,moho_only)
152152

153153
! using crustal values
154154
if (USE_OLD_VERSION_FORMAT) then
@@ -226,7 +226,7 @@ subroutine model_berkeley_crust_aniso(x,theta,phi,vpv,vph,vsv,vsh,eta_aniso,rho,
226226

227227
depth = (1.d0 - x) * EARTH_R_KM
228228

229-
call get_crust_val_csem(theta,phi,depth,rho,vp,vsv,vsh,moho_depth)
229+
call get_crust_val_csem(theta,phi,depth,rho,vp,vsv,vsh,moho_depth,moho_only)
230230

231231
! using crustal values
232232
if (USE_OLD_VERSION_FORMAT) then
@@ -266,14 +266,15 @@ end subroutine model_berkeley_crust_aniso
266266
!--------------------------------------------------------------------------------------------------
267267
!
268268

269-
subroutine get_crust_val_csem(theta,phi,z,rho,vp,vsv,vsh,moho_depth)
269+
subroutine get_crust_val_csem(theta,phi,z,rho,vp,vsv,vsh,moho_depth,moho_only)
270270

271271
use model_crust_berkeley_par
272272

273273
implicit none
274274

275275
double precision,intent(in) :: theta,phi,z
276276
double precision,intent(out) :: rho,vp,vsv,vsh,moho_depth
277+
logical,intent(in) :: moho_only
277278

278279
! local parameters
279280
! 4-th order GLL positions
@@ -290,13 +291,22 @@ subroutine get_crust_val_csem(theta,phi,z,rho,vp,vsv,vsh,moho_depth)
290291
double precision,external :: moho_filtre
291292
double precision,external :: lagrange
292293

294+
! initialize
295+
rho = 0.d0
296+
vp = 0.d0
297+
vsv = 0.d0
298+
vsh = 0.d0
299+
293300
! get moho depth
294301
moho_depth = moho1D_depth - moho_filtre(theta,phi)
295302

296303
!debug
297304
!print *,"debug: [get_crust_val_csem] Moho depth:",moho_depth, &
298305
! "moho1D_depth:",moho1D_depth,"moho_filtre:",moho_filtre(theta,phi)
299306

307+
! check if anything further to do or return only moho
308+
if (moho_only) return
309+
300310
!
301311
! horizontal interpolation for all registered depths
302312
!

‎src/shared/get_model_parameters.F90‎

Lines changed: 29 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -1190,7 +1190,8 @@ subroutine get_model_parameters_radii()
11901190
RHO_TOP_OC,RHO_BOTTOM_OC,RHO_OCEANS
11911191

11921192
use shared_parameters, only: &
1193-
HONOR_1D_SPHERICAL_MOHO,CASE_3D,CRUSTAL,REFERENCE_1D_MODEL
1193+
HONOR_1D_SPHERICAL_MOHO,CASE_3D,CRUSTAL,REFERENCE_1D_MODEL, &
1194+
NCHUNKS,NEX_XI,NEX_ETA
11941195

11951196
! reference models
11961197
use model_prem_par
@@ -1439,33 +1440,33 @@ subroutine get_model_parameters_radii()
14391440
! (same as values in model_ccrem.f90)
14401441
CCREM_RSURFACE = R_PLANET
14411442
ROCEAN = CCREM_RSURFACE ! no ocean
1442-
RMIDDLE_CRUST = CCREM_RSURFACE - 20000.d0 ! depth = 20 km
1443-
RMOHO = CCREM_RSURFACE - 35000.d0 ! depth = 35 km
1444-
R80 = CCREM_RSURFACE - 80000.d00 ! depth = 80 km
1445-
R220 = CCREM_RSURFACE - 220000.d0 ! depth = 220 km
1446-
R400 = CCREM_RSURFACE - 410000.d0 ! depth = 410 km - CCREM depth 410km discontinuity
1447-
R600 = CCREM_RSURFACE - 600000.d0 ! depth = 600 km
1448-
R670 = CCREM_RSURFACE - 660000.d0 ! depth = 660 km - CCREM depth 660km discontinuity
1449-
R771 = CCREM_RSURFACE - 771000.d0 ! depth = 771 km (PREM)
1443+
RMIDDLE_CRUST = CCREM_RSURFACE - 20000.d0 ! depth = 20 km
1444+
RMOHO = CCREM_RSURFACE - 35000.d0 ! depth = 35 km
1445+
R80 = CCREM_RSURFACE - 80000.d00 ! depth = 80 km
1446+
R220 = CCREM_RSURFACE - 220000.d0 ! depth = 220 km
1447+
R400 = CCREM_RSURFACE - 410000.d0 ! depth = 410 km - CCREM depth 410km discontinuity
1448+
R600 = CCREM_RSURFACE - 600000.d0 ! depth = 600 km
1449+
R670 = CCREM_RSURFACE - 660000.d0 ! depth = 660 km - CCREM depth 660km discontinuity
1450+
R771 = CCREM_RSURFACE - 771000.d0 ! depth = 771 km (PREM)
14501451
RTOPDDOUBLEPRIME = CCREM_RSURFACE - 2741000.d0 ! depth = 2741 km (PREM)
1451-
RCMB = CCREM_RSURFACE - 2891000.d0 ! depth = 2891 km (PREM)
1452-
RICB = CCREM_RSURFACE - 5153500.d0 ! depth = 5153.5 km
1452+
RCMB = CCREM_RSURFACE - 2891000.d0 ! depth = 2891 km (PREM)
1453+
RICB = CCREM_RSURFACE - 5153500.d0 ! depth = 5153.5 km
14531454

14541455
RHO_TOP_OC = 9.9131d0 * 1000.d0 / RHOAV
14551456
RHO_BOTTOM_OC = 12.1478d0 * 1000.d0 / RHOAV
14561457

14571458
case (REFERENCE_MODEL_SEMUCB)
14581459
! Berkeley SEMUCB Model - discontinuities
14591460
ROCEAN = 6368000.d0
1460-
RMIDDLE_CRUST = 6356000.d0
1461-
RMOHO = 6341000.d0 ! moho depth = 30 km
1462-
R80 = 6291000.d0
1463-
R120 = -1.d0 ! no d120 discontinuity, set to fictitious value
1464-
R220 = 6151000.d0
1465-
R400 = 5961000.d0
1466-
R600 = 5771000.d0
1467-
R670 = 5721000.d0
1468-
R771 = 5600000.d0
1461+
RMIDDLE_CRUST = 6356000.d0 ! depth = 15 km
1462+
RMOHO = 6341000.d0 ! moho depth = 30 km
1463+
R80 = 6291000.d0 ! depth = 80 km
1464+
R120 = -1.d0 ! no d120 discontinuity, set to fictitious value
1465+
R220 = 6151000.d0 ! depth = 220 km
1466+
R400 = 5961000.d0 ! depth = 410 km
1467+
R600 = 5771000.d0 ! depth = 600 km
1468+
R670 = 5721000.d0 ! depth = 650 km
1469+
R771 = 5600000.d0 ! depth = 771 km
14691470
RTOPDDOUBLEPRIME = 3630000.d0
14701471
RCMB = 3480000.d0
14711472
RICB = 1221500.d0
@@ -1605,6 +1606,14 @@ subroutine get_model_parameters_radii()
16051606
! moves fictitious moho closer to this 1d moho depth
16061607
RMOHO_FICTITIOUS_IN_MESHER = R_PLANET - 29000.000 ! down at 29 km depth
16071608
R80_FICTITIOUS_IN_MESHER = R_PLANET - 130000.000 ! down at 130 km depth
1609+
1610+
! coarse global mesh modification (to allow for small testing examples)
1611+
! to avoid Jacobian errors around lat/lon ~ 80 deg / 100 deg due to stretching
1612+
if (NCHUNKS == 6 .and. min(NEX_XI,NEX_ETA) <= 48) then
1613+
RMOHO_FICTITIOUS_IN_MESHER = R_PLANET - 35000.000 ! down at 35 km depth
1614+
R80_FICTITIOUS_IN_MESHER = R_PLANET - 140000.000 ! down at 100 km depth
1615+
endif
1616+
16081617
endif
16091618
endif
16101619

‎src/shared/read_compute_parameters.f90‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -524,6 +524,10 @@ subroutine rcp_check_parameters()
524524
! for flat topography, NEX = 32 setting still okay
525525
nex_minimum = 32
526526
endif
527+
528+
! Berkeley coarse global mesh modification (to allow for small testing examples)
529+
if (REFERENCE_1D_MODEL == REFERENCE_MODEL_SEMUCB) nex_minimum = 32
530+
527531
! checks nex
528532
if (NEX_XI < nex_minimum) &
529533
stop 'NEX_XI must be greater to cut the sphere into slices with positive Jacobian'

0 commit comments

Comments
 (0)