trun.py - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch)
(HTM) git clone git://src.adamsgaard.dk/pism
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
trun.py (11649B)
---
1 #!/usr/bin/env python3
2 import MISMIP
3
4 # This scripts generates bash scripts that run MISMIP experiments and generates
5 # all the necessary input files.
6 #
7 # Run run.py > my_new_mismip.sh and use that.
8
9 try:
10 from netCDF4 import Dataset as NC
11 except:
12 print("netCDF4 is not installed!")
13 sys.exit(1)
14
15 import sys
16
17 # The "standard" preamble used in many PISM scripts:
18 preamble = '''
19 if [ -n "${SCRIPTNAME:+1}" ] ; then
20 echo "[SCRIPTNAME=$SCRIPTNAME (already set)]"
21 echo ""
22 else
23 SCRIPTNAME="#(mismip.sh)"
24 fi
25
26 echo
27 echo "# =================================================================================="
28 echo "# MISMIP experiments"
29 echo "# =================================================================================="
30 echo
31
32 set -e # exit on error
33
34 NN=2 # default number of processors
35 if [ $# -gt 0 ] ; then # if user says "mismip.sh 8" then NN = 8
36 NN="$1"
37 fi
38
39 echo "$SCRIPTNAME NN = $NN"
40
41 # set MPIDO if using different MPI execution command, for example:
42 # $ export PISM_MPIDO="aprun -n "
43 if [ -n "${PISM_MPIDO:+1}" ] ; then # check if env var is already set
44 echo "$SCRIPTNAME PISM_MPIDO = $PISM_MPIDO (already set)"
45 else
46 PISM_MPIDO="mpiexec -n "
47 echo "$SCRIPTNAME PISM_MPIDO = $PISM_MPIDO"
48 fi
49
50 # check if env var PISM_DO was set (i.e. PISM_DO=echo for a 'dry' run)
51 if [ -n "${PISM_DO:+1}" ] ; then # check if env var DO is already set
52 echo "$SCRIPTNAME PISM_DO = $PISM_DO (already set)"
53 else
54 PISM_DO=""
55 fi
56
57 # prefix to pism (not to executables)
58 if [ -n "${PISM_BIN:+1}" ] ; then # check if env var is already set
59 echo "$SCRIPTNAME PISM_BIN = $PISM_BIN (already set)"
60 else
61 PISM_BIN="" # just a guess
62 echo "$SCRIPTNAME PISM_BIN = $PISM_BIN"
63 fi
64 '''
65
66
67 class Experiment:
68
69 "A MISMIP experiment."
70 experiment = ""
71 mode = 1
72 model = 1
73 semianalytic = True
74 Mx = 151
75 My = 3
76 Mz = 15
77 initials = "ABC"
78 executable = "$PISM_DO $PISM_MPIDO $NN ${PISM_BIN}pismr"
79
80 def __init__(self, experiment, model=1, mode=1, Mx=None, Mz=15, semianalytic=True,
81 initials="ABC", executable=None):
82 self.model = model
83 self.mode = mode
84 self.experiment = experiment
85 self.initials = initials
86 self.semianalytic = semianalytic
87
88 if executable:
89 self.executable = executable
90
91 if mode == 3:
92 self.Mx = Mx
93 else:
94 self.Mx = 2 * MISMIP.N(self.mode) + 1
95
96 self.My = 3
97
98 if self.experiment == "2b":
99 self.Lz = 7000
100 else:
101 self.Lz = 6000
102
103 def physics_options(self, input_file, step):
104 "Options corresponding to modeling choices."
105 config_filename = self.config(step)
106
107 options = ["-energy none", # isothermal setup; allows selecting cold-mode flow laws
108 "-ssa_flow_law isothermal_glen", # isothermal setup
109 "-yield_stress constant",
110 "-tauc %e" % MISMIP.C(self.experiment),
111 "-pseudo_plastic",
112 "-gradient eta",
113 "-pseudo_plastic_q %e" % MISMIP.m(self.experiment),
114 "-pseudo_plastic_uthreshold %e" % MISMIP.secpera(),
115 "-front_retreat_file %s" % input_file, # prescribe the maximum ice extent
116 "-config_override %s" % config_filename,
117 "-ssa_method fd",
118 "-cfbc", # calving front boundary conditions
119 "-part_grid", # sub-grid front motion parameterization
120 "-ssafd_ksp_rtol 1e-7",
121 "-ys 0",
122 "-ye %d" % MISMIP.run_length(self.experiment, step),
123 "-options_left",
124 ]
125
126 if self.model == 1:
127 options.extend(["-stress_balance ssa"])
128 else:
129 options.extend(["-stress_balance ssa+sia",
130 "-sia_flow_law isothermal_glen", # isothermal setup
131 ])
132
133 return options
134
135 def config(self, step):
136 '''Generates a config file containing flags and parameters
137 for a particular experiment and step.
138
139 This takes care of flags and parameters that *cannot* be set using
140 command-line options. (We try to use command-line options whenever we can.)
141 '''
142
143 filename = "MISMIP_conf_%s_A%s.nc" % (self.experiment, step)
144
145 nc = NC(filename, 'w', format="NETCDF3_CLASSIC")
146
147 var = nc.createVariable("pism_overrides", 'i')
148
149 attrs = {"ocean.always_grounded": "no",
150 "geometry.update.use_basal_melt_rate": "no",
151 "stress_balance.ssa.compute_surface_gradient_inward": "no",
152 "flow_law.isothermal_Glen.ice_softness": MISMIP.A(self.experiment, step),
153 "constants.ice.density": MISMIP.rho_i(),
154 "constants.sea_water.density": MISMIP.rho_w(),
155 "bootstrapping.defaults.geothermal_flux": 0.0,
156 "stress_balance.ssa.Glen_exponent": MISMIP.n(),
157 "constants.standard_gravity": MISMIP.g(),
158 "ocean.sub_shelf_heat_flux_into_ice": 0.0,
159 }
160
161 if self.model != 1:
162 attrs["stress_balance.sia.bed_smoother.range"] = 0.0
163
164 for name, value in attrs.items():
165 var.setncattr(name, value)
166
167 nc.close()
168
169 return filename
170
171 def bootstrap_options(self, step):
172 boot_filename = "MISMIP_boot_%s_M%s_A%s.nc" % (self.experiment, self.mode, step)
173
174 import prepare
175 prepare.pism_bootstrap_file(boot_filename, self.experiment, step, self.mode, N=self.Mx,
176 semianalytical_profile=self.semianalytic)
177
178 options = ["-i %s -bootstrap" % boot_filename,
179 "-Mx %d" % self.Mx,
180 "-My %d" % self.My,
181 "-Mz %d" % self.Mz,
182 "-Lz %d" % self.Lz]
183
184 return options, boot_filename
185
186 def output_options(self, step):
187 output_file = self.output_filename(self.experiment, step)
188 extra_file = "ex_" + output_file
189 ts_file = "ts_" + output_file
190
191 options = ["-extra_file %s" % extra_file,
192 "-extra_times 0:50:3e4",
193 "-extra_vars thk,topg,velbar_mag,flux_mag,mask,dHdt,usurf,hardav,velbase_mag,nuH,tauc,taud,taub,flux_divergence,cell_grounded_fraction",
194 "-ts_file %s" % ts_file,
195 "-ts_times 0:50:3e4",
196 "-o %s" % output_file,
197 "-o_order zyx",
198 ]
199
200 return output_file, options
201
202 def output_filename(self, experiment, step):
203 return "%s%s_%s_M%s_A%s.nc" % (self.initials, self.model, experiment, self.mode, step)
204
205 def options(self, step, input_file=None):
206 '''Generates a string of PISM options corresponding to a MISMIP experiment.'''
207
208 if input_file is None:
209 input_options, input_file = self.bootstrap_options(step)
210 else:
211 input_options = ["-i %s" % input_file]
212
213 physics = self.physics_options(input_file, step)
214
215 output_file, output_options = self.output_options(step)
216
217 return output_file, (input_options + physics + output_options)
218
219 def run_step(self, step, input_file=None):
220 out, opts = self.options(step, input_file)
221 print('echo "# Step %s-%s"' % (self.experiment, step))
222 print("%s %s" % (self.executable, ' '.join(opts)))
223 print('echo "Done."\n')
224
225 return out
226
227 def run(self, step=None):
228 print('echo "# Experiment %s"' % self.experiment)
229
230 if self.experiment in ('1a', '1b'):
231 # bootstrap
232 input_file = None
233 # steps 1 to 9
234 steps = list(range(1, 10))
235
236 if self.experiment in ('2a', '2b'):
237 # start from step 9 of the corresponding experiment 1
238 input_file = self.output_filename(self.experiment.replace("2", "1"), 9)
239 # steps 8 to 1
240 steps = list(range(8, 0, -1))
241
242 if self.experiment == '3a':
243 # bootstrap
244 input_file = None
245 # steps 1 to 13
246 steps = list(range(1, 14))
247
248 if self.experiment == '3b':
249 # bootstrap
250 input_file = None
251 # steps 1 to 15
252 steps = list(range(1, 16))
253
254 if step is not None:
255 input_file = None
256 steps = [step]
257
258 for step in steps:
259 input_file = self.run_step(step, input_file)
260
261
262 def run_mismip(initials, executable, semianalytic):
263 Mx = 601
264 models = (1, 2)
265 modes = (1, 2, 3)
266 experiments = ('1a', '1b', '2a', '2b', '3a', '3b')
267
268 print(preamble)
269
270 for model in models:
271 for mode in modes:
272 for experiment in experiments:
273 e = Experiment(experiment,
274 initials=initials,
275 executable=executable,
276 model=model, mode=mode, Mx=Mx,
277 semianalytic=semianalytic)
278 e.run()
279
280
281 if __name__ == "__main__":
282 from optparse import OptionParser
283
284 parser = OptionParser()
285
286 parser.usage = "%prog [options]"
287 parser.description = "Creates a script running MISMIP experiments."
288 parser.add_option("--initials", dest="initials", type="string",
289 help="Initials (3 letters)", default="ABC")
290 parser.add_option("-e", "--experiment", dest="experiment", type="string",
291 default='1a',
292 help="MISMIP experiments (one of '1a', '1b', '2a', '2b', '3a', '3b')")
293 parser.add_option("-s", "--step", dest="step", type="int",
294 help="MISMIP step number")
295 parser.add_option("-u", "--uniform_thickness",
296 action="store_false", dest="semianalytic", default=True,
297 help="Use uniform 10 m ice thickness")
298 parser.add_option("-a", "--all",
299 action="store_true", dest="all", default=False,
300 help="Run all experiments")
301 parser.add_option("-m", "--mode", dest="mode", type="int", default=1,
302 help="MISMIP grid mode")
303 parser.add_option("--Mx", dest="Mx", type="int", default=601,
304 help="Custom grid size; use with --mode=3")
305 parser.add_option("--model", dest="model", type="int", default=1,
306 help="Models: 1 - SSA only; 2 - SIA+SSA")
307 parser.add_option("--executable", dest="executable", type="string",
308 help="Executable to run, e.g. 'mpiexec -n 4 pismr'")
309
310 (opts, args) = parser.parse_args()
311
312 if opts.all:
313 run_mismip(opts.initials, opts.executable, opts.semianalytic)
314 exit(0)
315
316 def escape(arg):
317 if arg.find(" ") >= 0:
318 parts = arg.split("=")
319 return "%s=\"%s\"" % (parts[0], ' '.join(parts[1:]))
320 else:
321 return arg
322
323 arg_list = [escape(a) for a in sys.argv]
324
325 print("#!/bin/bash")
326 print("# This script was created by examples/mismip/run.py. The command was:")
327 print("# %s" % (' '.join(arg_list)))
328
329 if opts.executable is None:
330 print(preamble)
331
332 e = Experiment(opts.experiment,
333 initials=opts.initials,
334 executable=opts.executable,
335 model=opts.model,
336 mode=opts.mode,
337 Mx=opts.Mx,
338 semianalytic=opts.semianalytic)
339
340 if opts.step is not None:
341 e.run(opts.step)
342 else:
343 e.run()