tsg_create_flowline.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 --- tsg_create_flowline.py (2842B) --- 1 #!/usr/bin/env python3 2 # 3 # Copyright (C) 2011, 2014, 2018 Andy Aschwanden 4 # 5 # This file is part of PISM. 6 # 7 # PISM is free software; you can redistribute it and/or modify it under the 8 # terms of the GNU General Public License as published by the Free Software 9 # Foundation; either version 3 of the License, or (at your option) any later 10 # version. 11 # 12 # PISM is distributed in the hope that it will be useful, but WITHOUT ANY 13 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 14 # FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 15 # details. 16 # 17 # You should have received a copy of the GNU General Public License 18 # along with PISM; if not, write to the Free Software 19 # Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA 20 21 import numpy as np 22 23 try: 24 from netCDF4 import Dataset as CDF 25 except: 26 print("netCDF4 is not installed!") 27 sys.exit(1) 28 29 x, topg, thk = np.loadtxt('sg_35m_flowline.txt', unpack=True) 30 31 output = 'storglaciaren_flowline.nc' 32 33 # Write the data: 34 print("Writing the data to '%s'... " % output) 35 nc = CDF(output, "w") 36 nc.createDimension("x", size=len(x)) 37 38 x_var = nc.createVariable("x", 'f', dimensions=("x",)) 39 x_var.units = "m" 40 x_var[:] = x 41 42 topg_var = nc.createVariable("topg", 'f', dimensions=("x",)) 43 topg_var.units = "m" 44 topg_var.standard_name = "bedrock_altitude" 45 topg_var[:] = topg 46 47 thk_var = nc.createVariable("thk", 'f', dimensions=("x",)) 48 thk_var.units = "m" 49 thk_var.standard_name = "land_ice_thickness" 50 thk_var[:] = thk 51 52 usurf_var = nc.createVariable("usurf", 'f', dimensions=("x",)) 53 usurf_var.units = "m" 54 usurf_var.standard_name = "surface_altitude" 55 usurf_var[:] = topg + thk 56 57 qgeo = 0.042 58 bheatflx_var = nc.createVariable("bheatflx", 'f', dimensions=("x",)) 59 bheatflx_var.units = "W m-2" 60 bheatflx_var[:] = qgeo * np.ones_like(x) 61 62 # generate (somewhat) reasonable acab 63 acab_max = 2.5 # m/a 64 acab_min = -3.0 # m/a 65 np.linspace(acab_max, acab_min, 100) 66 acab = np.ones_like(x) 67 acab[:5] = 0 68 acab[5:105] = np.linspace(acab_max, acab_min, 100) 69 acab[105:] = acab_min 70 71 ice_density = 910.0 # kg m-3 72 acab_var = nc.createVariable("climatic_mass_balance", 'f', dimensions=("x",)) 73 acab_var.units = "kg m-2 year-1" 74 acab_var.standard_name = "land_ice_surface_specific_mass_balance_flux" 75 # convert from m/year ice equivalent into [kg m-2 year-1]: 76 acab_var[:] = acab * ice_density 77 78 # Set boundary conditions 79 # ------------------------------------------------------------------------------ 80 # 81 # (A) Surface temperature for temperature equation bc 82 Tma = -6.0 # degC, mean annual air temperature at Tarfala 83 zcts = 1400 # m a.s.l.; altitude where CTS is at the surface 84 85 86 artm = np.zeros_like(x) 87 artm[(topg + thk) < zcts] = Tma 88 artm_var = nc.createVariable("ice_surface_temp", 'f', dimensions=("x",)) 89 artm_var.units = "deg_C" 90 artm_var[:] = artm 91 92 nc.close()