GopherProxy

	trun.py - pism-exp-gsw - ice stream and sediment transport experiments
	git clone git://src.adamsgaard.dk/pism-exp-gsw
	Log
	Files
	Refs
	README
	LICENSE
	---
	trun.py (13935B)
	---
	1 #!/usr/bin/env python3
	2 import MISMIP
	3
	4 # This scripts generates bash scripts that run MISMIP experiments and ge…
	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, semi…
	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 self.Lz *= 2.0
	103
	104 def physics_options(self, input_file, step):
	105 "Options corresponding to modeling choices."
	106 config_filename = self.config(step)
	107
	108 #options = ["-energy none", # isothermal setup; allows selectin…
	109 #"-ssa_flow_law isothermal_glen", # isothermal setup
	110 #"-yield_stress constant",
	111 #"-tauc %e" % MISMIP.C(self.experiment),
	112 #"-pseudo_plastic",
	113 #"-gradient eta",
	114 #"-pseudo_plastic_q %e" % MISMIP.m(self.experiment),
	115 #"-pseudo_plastic_uthreshold %e" % MISMIP.secpera(),
	116 #"-front_retreat_file %s" % input_file, # prescribe t…
	117 #"-config_override %s" % config_filename,
	118 #"-ssa_method fd",
	119 #"-cfbc", # calving front boundary con…
	120 #"-part_grid", # sub-grid front motion para…
	121 #"-ssafd_ksp_rtol 1e-7",
	122 #"-ys 0",
	123 #"-ye %d" % MISMIP.run_length(self.experiment, step),
	124 #"-options_left",
	125 #]
	126
	127
	128 options = ["-stress_balance ssa+sia",
	129 #"-ssa_flow_law isothermal_glen", # isothermal setup
	130 #"-yield_stress constant",
	131 #"-tauc %e" % MISMIP.C(self.experiment),
	132 #"-pseudo_plastic",
	133 "-gradient eta",
	134 #"-pseudo_plastic_q %e" % MISMIP.m(self.experiment),
	135 #"-pseudo_plastic_uthreshold %e" % MISMIP.secpera(),
	136 #"-yield_stress mohr_coulomb",
	137 #"-till_flux",
	138 "-front_retreat_file %s" % input_file, # prescribe th…
	139 "-config_override %s" % config_filename,
	140 #"-ssa_method fd",
	141 "-cfbc", # PIK: calvin…
	142 "-part_grid", # PIK: sub-grid front mo…
	143 "-kill_icebergs", # PIK: https://pism-docs.org…
	144 "-subgl", # PIK: https://pism-…
	145 #"-ssafd_ksp_rtol 1e-7",
	146 "-ys 0",
	147 "-ye 1e3",
	148 #"-ye %d" % MISMIP.run_length(self.experiment, step),
	149 "-options_left",
	150 "-hydrology routing",
	151 "-plastic_phi 24", # Tulaczyk et al 2000
	152 "-till_cohesion 0", # Tulaczyk et al 2000
	153 "-geothermal_flux 70e-3", # Feldmann and Lever…
	154 "-sia_e 4.5", # Martin et al 2010
	155 "-ssa_e 0.512", # Martin et al 2010
	156 "-bed_def lc", # Lingle and Clark 1985, Bueler…
	157 "-stress_balance.sia.max_diffusivity 1e4",
	158 "-sea_level constant,delta_sl -ocean_delta_sl_file se…
	159 ]
	160
	161 #options = [,
	162 #"-yield_stress mohr_coulomb",
	163 #"-Lz 1e4",
	164 ##"-till_flux",
	165 #"-cfbc", # calving front boundary con…
	166 #"-part_grid", # sub-grid front motion para…
	167 #"-front_retreat_file %s" % input_file, # prescribe t…
	168 #"-ys 0",
	169 #"-ye %d" % MISMIP.run_length(self.experiment, step),
	170 #"-options_left",
	171 #]
	172
	173 return options
	174
	175 def config(self, step):
	176 '''Generates a config file containing flags and parameters
	177 for a particular experiment and step.
	178
	179 This takes care of flags and parameters that cannot be set usi…
	180 command-line options. (We try to use command-line options whenev…
	181 '''
	182
	183 filename = "MISMIP_conf_%s_A%s.nc" % (self.experiment, step)
	184
	185 nc = NC(filename, 'w', format="NETCDF3_CLASSIC")
	186
	187 var = nc.createVariable("pism_overrides", 'i')
	188
	189 attrs = {"geometry.update.use_basal_melt_rate": "no",
	190 "stress_balance.ssa.compute_surface_gradient_inward": "…
	191 "flow_law.isothermal_Glen.ice_softness": MISMIP.A(self.…
	192 "constants.ice.density": MISMIP.rho_i(),
	193 "constants.sea_water.density": MISMIP.rho_w(),
	194 "bootstrapping.defaults.geothermal_flux": 0.0,
	195 "stress_balance.ssa.Glen_exponent": MISMIP.n(),
	196 "constants.standard_gravity": MISMIP.g(),
	197 "ocean.sub_shelf_heat_flux_into_ice": 0.0,
	198 }
	199
	200 if self.model != 1:
	201 attrs["stress_balance.sia.bed_smoother.range"] = 0.0
	202
	203 for name, value in attrs.items():
	204 var.setncattr(name, value)
	205
	206 nc.close()
	207
	208 return filename
	209
	210 def bootstrap_options(self, step):
	211 boot_filename = "MISMIP_boot_%s_M%s_A%s.nc" % (self.experiment, …
	212
	213 import prepare
	214 prepare.pism_bootstrap_file(boot_filename, self.experiment, step…
	215 semianalytical_profile=self.semianal…
	216
	217 options = ["-i %s -bootstrap" % boot_filename,
	218 "-Mx %d" % self.Mx,
	219 "-My %d" % self.My,
	220 "-Mz %d" % self.Mz,
	221 "-Lz %d" % self.Lz]
	222
	223 return options, boot_filename
	224
	225 def output_options(self, step):
	226 output_file = self.output_filename(self.experiment, step)
	227 extra_file = "ex_" + output_file
	228 ts_file = "ts_" + output_file
	229
	230 options = ["-extra_file %s" % extra_file,
	231 "-extra_times 0:50:3e4",
	232 "-extra_vars thk,topg,velbar_mag,flux_mag,mask,dHdt,u…
	233 "-ts_file %s" % ts_file,
	234 "-ts_times 0:50:3e4",
	235 "-backup_size big",
	236 "-o %s" % output_file,
	237 "-o_order zyx",
	238 ]
	239
	240 return output_file, options
	241
	242 def output_filename(self, experiment, step):
	243 return "%s%s_%s_M%s_A%s.nc" % (self.initials, self.model, experi…
	244
	245 def options(self, step, input_file=None):
	246 '''Generates a string of PISM options corresponding to a MISMIP …
	247
	248 if input_file is None:
	249 input_options, input_file = self.bootstrap_options(step)
	250 else:
	251 input_options = ["-i %s" % input_file]
	252
	253 physics = self.physics_options(input_file, step)
	254
	255 output_file, output_options = self.output_options(step)
	256
	257 return output_file, (input_options + physics + output_options)
	258
	259 def run_step(self, step, input_file=None):
	260 out, opts = self.options(step, input_file)
	261 print('echo "# Step %s-%s"' % (self.experiment, step))
	262 print("%s %s" % (self.executable, ' '.join(opts)))
	263 print('echo "Done."\n')
	264
	265 return out
	266
	267 def run(self, step=None):
	268 print('echo "# Experiment %s"' % self.experiment)
	269
	270 if self.experiment in ('1a', '1b'):
	271 # bootstrap
	272 input_file = None
	273 # steps 1 to 9
	274 steps = list(range(1, 10))
	275
	276 if self.experiment in ('2a', '2b'):
	277 # start from step 9 of the corresponding experiment 1
	278 input_file = self.output_filename(self.experiment.replace("2…
	279 # steps 8 to 1
	280 steps = list(range(8, 0, -1))
	281
	282 if self.experiment == '3a':
	283 # bootstrap
	284 input_file = None
	285 # steps 1 to 13
	286 steps = list(range(1, 14))
	287
	288 if self.experiment == '3b':
	289 # bootstrap
	290 input_file = None
	291 # steps 1 to 15
	292 steps = list(range(1, 16))
	293
	294 if step is not None:
	295 input_file = None
	296 steps = [step]
	297
	298 for step in steps:
	299 input_file = self.run_step(step, input_file)
	300
	301
	302 def run_mismip(initials, executable, semianalytic):
	303 Mx = 601
	304 models = (1, 2)
	305 modes = (1, 2, 3)
	306 experiments = ('1a', '1b', '2a', '2b', '3a', '3b')
	307
	308 print(preamble)
	309
	310 for model in models:
	311 for mode in modes:
	312 for experiment in experiments:
	313 e = Experiment(experiment,
	314 initials=initials,
	315 executable=executable,
	316 model=model, mode=mode, Mx=Mx,
	317 semianalytic=semianalytic)
	318 e.run()
	319
	320
	321 if __name__ == "__main__":
	322 from optparse import OptionParser
	323
	324 parser = OptionParser()
	325
	326 parser.usage = "%prog [options]"
	327 parser.description = "Creates a script running MISMIP experiments."
	328 parser.add_option("--initials", dest="initials", type="string",
	329 help="Initials (3 letters)", default="ABC")
	330 parser.add_option("-e", "--experiment", dest="experiment", type="str…
	331 default='1a',
	332 help="MISMIP experiments (one of '1a', '1b', '2a',…
	333 parser.add_option("-s", "--step", dest="step", type="int",
	334 help="MISMIP step number")
	335 parser.add_option("-u", "--uniform_thickness",
	336 action="store_false", dest="semianalytic", default…
	337 help="Use uniform 10 m ice thickness")
	338 parser.add_option("-a", "--all",
	339 action="store_true", dest="all", default=False,
	340 help="Run all experiments")
	341 parser.add_option("-m", "--mode", dest="mode", type="int", default=1,
	342 help="MISMIP grid mode")
	343 parser.add_option("--Mx", dest="Mx", type="int", default=601,
	344 help="Custom grid size; use with --mode=3")
	345 parser.add_option("--model", dest="model", type="int", default=1,
	346 help="Models: 1 - SSA only; 2 - SIA+SSA")
	347 parser.add_option("--executable", dest="executable", type="string",
	348 help="Executable to run, e.g. 'mpiexec -n 4 pismr'…
	349
	350 (opts, args) = parser.parse_args()
	351
	352 if opts.all:
	353 run_mismip(opts.initials, opts.executable, opts.semianalytic)
	354 exit(0)
	355
	356 def escape(arg):
	357 if arg.find(" ") >= 0:
	358 parts = arg.split("=")
	359 return "%s=\"%s\"" % (parts[0], ' '.join(parts[1:]))
	360 else:
	361 return arg
	362
	363 arg_list = [escape(a) for a in sys.argv]
	364
	365 print("#!/bin/bash")
	366 print("# This script was created by examples/mismip/run.py. The comm…
	367 print("# %s" % (' '.join(arg_list)))
	368
	369 if opts.executable is None:
	370 print(preamble)
	371
	372 e = Experiment(opts.experiment,
	373 initials=opts.initials,
	374 executable=opts.executable,
	375 model=opts.model,
	376 mode=opts.mode,
	377 Mx=opts.Mx,
	378 semianalytic=opts.semianalytic)
	379
	380 if opts.step is not None:
	381 e.run(opts.step)
	382 else:
	383 e.run()