source: issm/trunk-jpl/examples/shakti/moulin.py@ 25482

Last change on this file since 25482 was 25482, checked in by bdef, 5 years ago

NEW: py translation of Shakti example

File size: 2.1 KB
Line 
1import numpy as np
2from SetIceSheetBC import SetIceSheetBC
3from frictionshakti import frictionshakti
4from verbose import verbose
5
6#Start defining model parameters here
7# Set up bed topography and ice geometry for a tilted 500m thick slab
8md.geometry.base = .02 * md.mesh.x
9md.geometry.bed = md.geometry.base
10md.geometry.surface = .02 * md.mesh.x + 500
11md.geometry.thickness = md.geometry.surface - md.geometry.bed
12
13# Define ice sliding velocity (m / yr)
14md.initialization.vx = 1.0e-6 * md.constants.yts * np.ones(md.mesh.numberofvertices)
15md.initialization.vy = np.zeros(md.mesh.numberofvertices)
16md.initialization.pressure = np.zeros(md.mesh.numberofvertices)
17
18# Materials
19# Ice flow law parameter (note that the standard parameter A = B^(- 3))
20md.materials.rheology_B = 5e-25**(-1 / 3) * np.ones(md.mesh.numberofvertices)
21md.initialization.temperature = 273. * np.ones(md.mesh.numberofvertices)
22md.materials.rheology_n = 3. * np.ones(md.mesh.numberofelements)
23
24#Calving
25md.calving.calvingrate = np.zeros(md.mesh.numberofvertices)
26#md.calving.spclevelset = NaN(md.mesh.numberofvertices, 1)
27
28# Friction law and coefficient
29md.friction = frictionshakti(md)
30md.friction.coefficient = 20. * np.ones(md.mesh.numberofvertices)
31
32#Numerical parameters
33#md.stressbalance.viscosity_overshoot = 0.0
34md.masstransport.stabilization = 1.
35md.thermal.stabilization = 1.
36md.verbose = verbose(0)
37md.settings.waitonlock = 30
38md.stressbalance.restol = 0.05
39md.steadystate.reltol = 0.05
40md.stressbalance.reltol = 0.05
41md.stressbalance.abstol = np.nan
42md.timestepping.time_step = 1.
43md.timestepping.final_time = 3.
44
45#GIA:
46md.gia.lithosphere_thickness = 100. * np.ones(md.mesh.numberofvertices) # in km
47md.gia.mantle_viscosity = 1.0e21 * np.ones(md.mesh.numberofvertices) # in Pa.s
48md.materials.lithosphere_shear_modulus = 6.7e10 # in Pa
49md.materials.lithosphere_density = 3.32 # in g / cm^ - 3
50md.materials.mantle_shear_modulus = 1.45e11 # in Pa
51md.materials.mantle_density = 3.34 # in g / cm^ - 3
52
53#Boundary conditions:
54md = SetIceSheetBC(md)
55
56#Change name so that no test have the same name
57md.private.runtimename = True
Note: See TracBrowser for help on using the repository browser.