SICOPOLIS V5-dev  Revision 1173
read_m.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Module : r e a d _ m
4 !
5 !> @file
6 !!
7 !! Reading of several input data.
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2009-2017 Ralf Greve
12 !!
13 !! @section License
14 !!
15 !! This file is part of SICOPOLIS.
16 !!
17 !! SICOPOLIS is free software: you can redistribute it and/or modify
18 !! it under the terms of the GNU General Public License as published by
19 !! the Free Software Foundation, either version 3 of the License, or
20 !! (at your option) any later version.
21 !!
22 !! SICOPOLIS is distributed in the hope that it will be useful,
23 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
24 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 !! GNU General Public License for more details.
26 !!
27 !! You should have received a copy of the GNU General Public License
28 !! along with SICOPOLIS. If not, see <http://www.gnu.org/licenses/>.
29 !<
30 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
31 
32 !-------------------------------------------------------------------------------
33 !> Reading of several input data.
34 !<------------------------------------------------------------------------------
35 module read_m
36 
37  use sico_types_m
38 
39  implicit none
40 
41  private
43 
44 contains
45 
46 !-------------------------------------------------------------------------------
47 !> Reading of time-slice files in native binary or NetCDF format.
48 !<------------------------------------------------------------------------------
49  subroutine read_erg_nc(z_sl, filename)
50 
52  use sico_vars_m
53 
54 #if (NETCDF==2) /* time-slice file in NetCDF format */
55  use netcdf
56  use nc_check_m
57 #endif
58 
59  implicit none
60 
61  character(len=100), intent(in) :: filename
62 
63  real(dp), intent(out) :: z_sl
64 
65  integer(i4b) :: i, j, kc, kt, kr, n
66  character(len=256) :: filename_with_path
67 
68 #if (NETCDF==1) /* time-slice file in native binary format */
69 
70  integer(i4b) :: ios
71  character(len=256) :: ch_attr_title, ch_attr_institution, ch_attr_source, &
72  ch_attr_history, ch_attr_references
73 
74 #elif (NETCDF==2) /* time-slice file in NetCDF format */
75 
76  integer(i4b) :: ncid, ncv
77  ! ncid: ID of the output file
78  ! ncv: Variable ID
79 
80 #else
81  stop ' >>> read_erg_nc: Parameter NETCDF must be either 1 or 2!'
82 #endif
83 
84  integer(i2b), dimension(0:IMAX,0:JMAX) :: maske_conv, maske_old_conv, &
85  n_cts_conv, kc_cts_conv
86  real(sp) :: time_conv, dummy_conv, z_sl_conv, &
87  V_tot_conv, A_grounded_conv, A_floating_conv, &
88  H_R_conv, &
89  xi_conv(0:imax), eta_conv(0:jmax), &
90  sigma_level_c_conv(0:kcmax), sigma_level_t_conv(0:ktmax), &
91  sigma_level_r_conv(0:krmax)
92  real(sp), dimension(0:IMAX,0:JMAX) :: lambda_conv, phi_conv, &
93  lon_conv, lat_conv, &
94  temp_s_conv, accum_conv, &
95  as_perp_conv, as_perp_apl_conv, smb_corr_conv, &
96  q_geo_conv, &
97  zs_conv, zm_conv, zb_conv, zl_conv, zl0_conv, &
98  H_cold_conv, H_temp_conv, H_conv, &
99  Q_bm_conv, Q_tld_conv, &
100  am_perp_conv, &
101  qx_conv, qy_conv, &
102  dzs_dtau_conv, dzm_dtau_conv, dzb_dtau_conv, dzl_dtau_conv, &
103  dH_c_dtau_conv, dH_t_dtau_conv, dH_dtau_conv, &
104  vx_b_g_conv, vy_b_g_conv, vz_b_conv, vh_b_conv, &
105  vx_s_g_conv, vy_s_g_conv, vz_s_conv, vh_s_conv, &
106  vx_m_g_conv, vy_m_g_conv, vh_m_conv, &
107  temp_b_conv, temph_b_conv, &
108  tau_b_driving_conv, tau_b_drag_conv, &
109  p_b_w_conv, H_w_conv, q_gl_g_conv, q_cf_g_conv
110  real(sp), dimension(0:IMAX,0:JMAX,0:KCMAX) :: vx_c_conv, vy_c_conv, vz_c_conv, &
111  temp_c_conv, age_c_conv, &
112  enth_c_conv, omega_c_conv, &
113  enh_c_conv
114  real(sp), dimension(0:IMAX,0:JMAX,0:KTMAX) :: vx_t_conv, vy_t_conv, vz_t_conv, &
115  omega_t_conv, age_t_conv, &
116  enth_t_conv, &
117  enh_t_conv
118  real(sp), dimension(0:IMAX,0:JMAX,0:KRMAX) :: temp_r_conv
119 
120 #if (DISC>0) /* Ice discharge parameterisation */
121  integer(i2b), dimension(0:IMAX,0:JMAX) :: mask_mar_conv
122  real(sp), dimension(0:IMAX,0:JMAX) :: dis_perp_conv, &
123  cst_dist_conv, cos_grad_tc_conv
124 #endif
125 
126 !-------- Read data from time-slice file of previous simulation --------
127 
128  filename_with_path = trim(anfdatpath)//'/'//trim(filename)
129 
130 #if (NETCDF==1) /* time-slice file in native binary format */
131 
132  open(unit=11, iostat=ios, file=trim(filename_with_path), status='old', &
133  form='unformatted')
134 
135  if (ios /= 0) stop ' >>> read_erg_nc: Error when opening the time-slice file!'
136 
137  read(unit=11) ch_attr_title
138  read(unit=11) ch_attr_institution
139  read(unit=11) ch_attr_source
140  read(unit=11) ch_attr_history
141  read(unit=11) ch_attr_references
142 
143  read(unit=11) time_conv
144  read(unit=11) dummy_conv ! this is either delta_ts or glac_index; not needed
145  read(unit=11) z_sl_conv
146 
147  read(unit=11) dummy_conv ! this is V_tot; not needed
148  read(unit=11) dummy_conv ! this is V_af; not needed
149  read(unit=11) dummy_conv ! this is A_grounded; not needed
150  read(unit=11) dummy_conv ! this is A_floating; not needed
151 
152  read(unit=11) xi_conv
153  read(unit=11) eta_conv
154  read(unit=11) sigma_level_c_conv
155  read(unit=11) sigma_level_t_conv
156  read(unit=11) sigma_level_r_conv
157  read(unit=11) lon_conv
158  read(unit=11) lat_conv
159  read(unit=11) lambda_conv
160  read(unit=11) phi_conv
161  read(unit=11) temp_s_conv
162  read(unit=11) accum_conv
163  read(unit=11) as_perp_conv
164  read(unit=11) as_perp_apl_conv
165  read(unit=11) smb_corr_conv
166 
167 #if (DISC>0) /* Ice discharge parameterisation */
168 
169  read(unit=11) dis_perp_conv
170  read(unit=11) cst_dist_conv
171  read(unit=11) cos_grad_tc_conv
172  read(unit=11) mask_mar_conv
173 
174 #endif
175 
176  read(unit=11) q_geo_conv
177  read(unit=11) maske_conv
178  read(unit=11) maske_old_conv
179  read(unit=11) n_cts_conv
180  read(unit=11) kc_cts_conv
181  read(unit=11) zs_conv
182  read(unit=11) zm_conv
183  read(unit=11) zb_conv
184  read(unit=11) zl_conv
185  read(unit=11) zl0_conv
186  read(unit=11) h_cold_conv
187  read(unit=11) h_temp_conv
188  read(unit=11) h_conv
189  read(unit=11) h_r_conv
190  read(unit=11) q_bm_conv
191  read(unit=11) q_tld_conv
192  read(unit=11) am_perp_conv
193  read(unit=11) qx_conv
194  read(unit=11) qy_conv
195  read(unit=11) dzs_dtau_conv
196  read(unit=11) dzm_dtau_conv
197  read(unit=11) dzb_dtau_conv
198  read(unit=11) dzl_dtau_conv
199  read(unit=11) dh_c_dtau_conv
200  read(unit=11) dh_t_dtau_conv
201  read(unit=11) dh_dtau_conv
202  read(unit=11) vx_b_g_conv
203  read(unit=11) vy_b_g_conv
204  read(unit=11) vz_b_conv
205  read(unit=11) vh_b_conv
206  read(unit=11) vx_s_g_conv
207  read(unit=11) vy_s_g_conv
208  read(unit=11) vz_s_conv
209  read(unit=11) vh_s_conv
210  read(unit=11) vx_m_g_conv
211  read(unit=11) vy_m_g_conv
212  read(unit=11) vh_m_conv
213  read(unit=11) temp_b_conv
214  read(unit=11) temph_b_conv
215  read(unit=11) tau_b_driving_conv
216  read(unit=11) tau_b_drag_conv
217  read(unit=11) p_b_w_conv
218  read(unit=11) h_w_conv
219  read(unit=11) q_gl_g_conv
220  read(unit=11) q_cf_g_conv
221 
222  read(unit=11) vx_c_conv
223  read(unit=11) vy_c_conv
224  read(unit=11) vz_c_conv
225  read(unit=11) vx_t_conv
226  read(unit=11) vy_t_conv
227  read(unit=11) vz_t_conv
228  read(unit=11) temp_c_conv
229  read(unit=11) omega_t_conv
230  read(unit=11) temp_r_conv
231  read(unit=11) enth_c_conv
232  read(unit=11) enth_t_conv
233  read(unit=11) omega_c_conv
234  read(unit=11) enh_c_conv
235  read(unit=11) enh_t_conv
236  read(unit=11) age_c_conv
237  read(unit=11) age_t_conv
238 
239  close(unit=11,status='keep')
240 
241 #elif (NETCDF==2) /* time-slice file in NetCDF format */
242 
243  call check( nf90_open(trim(filename_with_path), nf90_nowrite, ncid) )
244 
245  call check( nf90_inq_varid(ncid, 'time', ncv) )
246  call check( nf90_get_var(ncid, ncv, time_conv) )
247 
248  if ( nf90_inq_varid(ncid, 'delta_ts', ncv) == nf90_noerr ) then
249  call check( nf90_get_var(ncid, ncv, dummy_conv) )
250  else if ( nf90_inq_varid(ncid, 'glac_index', ncv) == nf90_noerr ) then
251  call check( nf90_get_var(ncid, ncv, dummy_conv) )
252  else
253  dummy_conv = 0.0_sp
254  end if
255 
256  call check( nf90_inq_varid(ncid, 'z_sl', ncv) )
257  call check( nf90_get_var(ncid, ncv, z_sl_conv) )
258 
259  call check( nf90_inq_varid(ncid, 'x', ncv) )
260  call check( nf90_get_var(ncid, ncv, xi_conv) )
261 
262  call check( nf90_inq_varid(ncid, 'y', ncv) )
263  call check( nf90_get_var(ncid, ncv, eta_conv) )
264 
265  call check( nf90_inq_varid(ncid, 'sigma_level_c', ncv) )
266  call check( nf90_get_var(ncid, ncv, sigma_level_c_conv) )
267 
268  call check( nf90_inq_varid(ncid, 'sigma_level_t', ncv) )
269  call check( nf90_get_var(ncid, ncv, sigma_level_t_conv) )
270 
271  call check( nf90_inq_varid(ncid, 'sigma_level_r', ncv) )
272  call check( nf90_get_var(ncid, ncv, sigma_level_r_conv) )
273 
274  call check( nf90_inq_varid(ncid, 'lon', ncv) )
275  call check( nf90_get_var(ncid, ncv, lon_conv) )
276 
277  call check( nf90_inq_varid(ncid, 'lat', ncv) )
278  call check( nf90_get_var(ncid, ncv, lat_conv) )
279 
280  call check( nf90_inq_varid(ncid, 'lambda', ncv) )
281  call check( nf90_get_var(ncid, ncv, lambda_conv) )
282 
283  call check( nf90_inq_varid(ncid, 'phi', ncv) )
284  call check( nf90_get_var(ncid, ncv, phi_conv) )
285 
286  call check( nf90_inq_varid(ncid, 'temp_s', ncv) )
287  call check( nf90_get_var(ncid, ncv, temp_s_conv) )
288 
289  if ( nf90_inq_varid(ncid, 'prec', ncv) == nf90_noerr ) then
290  call check( nf90_get_var(ncid, ncv, accum_conv) )
291  else
292  write(6,'(/1x,a)') '>>> read_erg_nc: Variable prec'
293  write(6, '(1x,a)') ' not available in read file *.nc.'
294  accum_conv = 0.0_sp
295  end if
296 
297  call check( nf90_inq_varid(ncid, 'as_perp', ncv) )
298  call check( nf90_get_var(ncid, ncv, as_perp_conv) )
299 
300  if ( nf90_inq_varid(ncid, 'as_perp_apl', ncv) == nf90_noerr ) then
301  call check( nf90_get_var(ncid, ncv, as_perp_apl_conv) )
302  else
303  write(6,'(/1x,a)') '>>> read_erg_nc: Variable as_perp_apl'
304  write(6, '(1x,a)') ' not available in read file *.nc.'
305  as_perp_apl_conv = 0.0_sp
306  end if
307 
308  if ( nf90_inq_varid(ncid, 'smb_corr', ncv) == nf90_noerr ) then
309  call check( nf90_get_var(ncid, ncv, smb_corr_conv) )
310  else
311  write(6,'(/1x,a)') '>>> read_erg_nc: Variable smb_corr'
312  write(6, '(1x,a)') ' not available in read file *.nc.'
313  smb_corr_conv = 0.0_sp
314  end if
315 
316 #if (DISC>0) /* Ice discharge parameterisation */
317 
318  if ( nf90_inq_varid(ncid, 'dis_perp', ncv) == nf90_noerr ) then
319  call check( nf90_get_var(ncid, ncv, dis_perp_conv) )
320  else
321  write(6,'(/1x,a)') '>>> read_erg_nc: Variable dis_perp'
322  write(6, '(1x,a)') ' not available in read file *.nc.'
323  dis_perp_conv = 0.0_sp
324  end if
325 
326  if ( nf90_inq_varid(ncid, 'cst_dist', ncv) == nf90_noerr ) then
327  call check( nf90_get_var(ncid, ncv, cst_dist_conv) )
328  else
329  write(6,'(/1x,a)') '>>> read_erg_nc: Variable cst_dist'
330  write(6, '(1x,a)') ' not available in read file *.nc.'
331  cst_dist_conv = real(no_value_pos_1,sp)
332  end if
333 
334  if ( nf90_inq_varid(ncid, 'cos_grad_tc', ncv) == nf90_noerr ) then
335  call check( nf90_get_var(ncid, ncv, cos_grad_tc_conv) )
336  else
337  write(6,'(/1x,a)') '>>> read_erg_nc: Variable cos_grad_tc'
338  write(6, '(1x,a)') ' not available in read file *.nc.'
339  cos_grad_tc_conv = 1.0_sp
340  end if
341 
342  if ( nf90_inq_varid(ncid, 'mask_mar', ncv) == nf90_noerr ) then
343  call check( nf90_get_var(ncid, ncv, mask_mar_conv) )
344  else
345  write(6,'(/1x,a)') '>>> read_erg_nc: Variable mask_mar'
346  write(6, '(1x,a)') ' not available in read file *.nc.'
347  mask_mar_conv = 0_i2b
348  end if
349 
350 #endif
351 
352  if ( nf90_inq_varid(ncid, 'q_geo', ncv) == nf90_noerr ) then
353  call check( nf90_get_var(ncid, ncv, q_geo_conv) )
354  else
355  write(6,'(/1x,a)') '>>> read_erg_nc: Variable q_geo'
356  write(6, '(1x,a)') ' not available in read file *.nc.'
357  q_geo_conv = 0.0_sp
358  end if
359 
360  call check( nf90_inq_varid(ncid, 'maske', ncv) )
361  call check( nf90_get_var(ncid, ncv, maske_conv) )
362 
363  if ( nf90_inq_varid(ncid, 'maske_old', ncv) == nf90_noerr ) then
364  call check( nf90_get_var(ncid, ncv, maske_old_conv) )
365  else
366  write(6,'(/1x,a)') '>>> read_erg_nc: Variable maske_old'
367  write(6, '(1x,a)') ' not available in read file *.nc.'
368  maske_old_conv = maske_conv
369  end if
370 
371  call check( nf90_inq_varid(ncid, 'n_cts', ncv) )
372  call check( nf90_get_var(ncid, ncv, n_cts_conv) )
373 
374  call check( nf90_inq_varid(ncid, 'kc_cts', ncv) )
375  call check( nf90_get_var(ncid, ncv, kc_cts_conv) )
376 
377  call check( nf90_inq_varid(ncid, 'zs', ncv) )
378  call check( nf90_get_var(ncid, ncv, zs_conv) )
379 
380  call check( nf90_inq_varid(ncid, 'zm', ncv) )
381  call check( nf90_get_var(ncid, ncv, zm_conv) )
382 
383  call check( nf90_inq_varid(ncid, 'zb', ncv) )
384  call check( nf90_get_var(ncid, ncv, zb_conv) )
385 
386  call check( nf90_inq_varid(ncid, 'zl', ncv) )
387  call check( nf90_get_var(ncid, ncv, zl_conv) )
388 
389  if ( nf90_inq_varid(ncid, 'zl0', ncv) == nf90_noerr ) then
390  call check( nf90_get_var(ncid, ncv, zl0_conv) )
391  else
392  write(6,'(/1x,a)') '>>> read_erg_nc: Variable zl0'
393  write(6, '(1x,a)') ' not available in read file *.nc.'
394  zl0_conv = 0.0_sp
395  end if
396 
397  call check( nf90_inq_varid(ncid, 'H_cold', ncv) )
398  call check( nf90_get_var(ncid, ncv, h_cold_conv) )
399 
400  call check( nf90_inq_varid(ncid, 'H_temp', ncv) )
401  call check( nf90_get_var(ncid, ncv, h_temp_conv) )
402 
403  call check( nf90_inq_varid(ncid, 'H', ncv) )
404  call check( nf90_get_var(ncid, ncv, h_conv) )
405 
406  call check( nf90_inq_varid(ncid, 'H_R', ncv) )
407  call check( nf90_get_var(ncid, ncv, h_r_conv) )
408 
409  call check( nf90_inq_varid(ncid, 'Q_bm', ncv) )
410  call check( nf90_get_var(ncid, ncv, q_bm_conv) )
411 
412  call check( nf90_inq_varid(ncid, 'Q_tld', ncv) )
413  call check( nf90_get_var(ncid, ncv, q_tld_conv) )
414 
415  call check( nf90_inq_varid(ncid, 'am_perp', ncv) )
416  call check( nf90_get_var(ncid, ncv, am_perp_conv) )
417 
418  call check( nf90_inq_varid(ncid, 'qx', ncv) )
419  call check( nf90_get_var(ncid, ncv, qx_conv) )
420 
421  call check( nf90_inq_varid(ncid, 'qy', ncv) )
422  call check( nf90_get_var(ncid, ncv, qy_conv) )
423 
424  call check( nf90_inq_varid(ncid, 'dzs_dt', ncv) )
425  call check( nf90_get_var(ncid, ncv, dzs_dtau_conv) )
426 
427  call check( nf90_inq_varid(ncid, 'dzm_dt', ncv) )
428  call check( nf90_get_var(ncid, ncv, dzm_dtau_conv) )
429 
430  call check( nf90_inq_varid(ncid, 'dzb_dt', ncv) )
431  call check( nf90_get_var(ncid, ncv, dzb_dtau_conv) )
432 
433  call check( nf90_inq_varid(ncid, 'dzl_dt', ncv) )
434  call check( nf90_get_var(ncid, ncv, dzl_dtau_conv) )
435 
436  call check( nf90_inq_varid(ncid, 'dH_c_dt', ncv) )
437  call check( nf90_get_var(ncid, ncv, dh_c_dtau_conv) )
438 
439  call check( nf90_inq_varid(ncid, 'dH_t_dt', ncv) )
440  call check( nf90_get_var(ncid, ncv, dh_t_dtau_conv) )
441 
442  call check( nf90_inq_varid(ncid, 'dH_dt', ncv) )
443  call check( nf90_get_var(ncid, ncv, dh_dtau_conv) )
444 
445  call check( nf90_inq_varid(ncid, 'vx_b_g', ncv) )
446  call check( nf90_get_var(ncid, ncv, vx_b_g_conv) )
447 
448  call check( nf90_inq_varid(ncid, 'vy_b_g', ncv) )
449  call check( nf90_get_var(ncid, ncv, vy_b_g_conv) )
450 
451  call check( nf90_inq_varid(ncid, 'vz_b', ncv) )
452  call check( nf90_get_var(ncid, ncv, vz_b_conv) )
453 
454  call check( nf90_inq_varid(ncid, 'vh_b', ncv) )
455  call check( nf90_get_var(ncid, ncv, vh_b_conv) )
456 
457  call check( nf90_inq_varid(ncid, 'vx_s_g', ncv) )
458  call check( nf90_get_var(ncid, ncv, vx_s_g_conv) )
459 
460  call check( nf90_inq_varid(ncid, 'vy_s_g', ncv) )
461  call check( nf90_get_var(ncid, ncv, vy_s_g_conv) )
462 
463  call check( nf90_inq_varid(ncid, 'vz_s', ncv) )
464  call check( nf90_get_var(ncid, ncv, vz_s_conv) )
465 
466  call check( nf90_inq_varid(ncid, 'vh_s', ncv) )
467  call check( nf90_get_var(ncid, ncv, vh_s_conv) )
468 
469  if ( nf90_inq_varid(ncid, 'vx_m_g', ncv) == nf90_noerr ) then
470  call check( nf90_get_var(ncid, ncv, vx_m_g_conv) )
471  else
472  write(6,'(/1x,a)') '>>> read_erg_nc: Variable vx_m_g'
473  write(6, '(1x,a)') ' not available in read file *.nc.'
474  vx_m_g_conv = 0.0_sp
475  end if
476 
477  if ( nf90_inq_varid(ncid, 'vy_m_g', ncv) == nf90_noerr ) then
478  call check( nf90_get_var(ncid, ncv, vy_m_g_conv) )
479  else
480  write(6,'(/1x,a)') '>>> read_erg_nc: Variable vy_m_g'
481  write(6, '(1x,a)') ' not available in read file *.nc.'
482  vy_m_g_conv = 0.0_sp
483  end if
484 
485  if ( nf90_inq_varid(ncid, 'vh_m', ncv) == nf90_noerr ) then
486  call check( nf90_get_var(ncid, ncv, vh_m_conv) )
487  else
488  write(6,'(/1x,a)') '>>> read_erg_nc: Variable vh_m'
489  write(6, '(1x,a)') ' not available in read file *.nc.'
490  vh_m_conv = 0.0_sp
491  end if
492 
493  call check( nf90_inq_varid(ncid, 'temp_b', ncv) )
494  call check( nf90_get_var(ncid, ncv, temp_b_conv) )
495 
496  call check( nf90_inq_varid(ncid, 'temph_b', ncv) )
497  call check( nf90_get_var(ncid, ncv, temph_b_conv) )
498 
499  if ( nf90_inq_varid(ncid, 'tau_b_driving', ncv) == nf90_noerr ) then
500  call check( nf90_get_var(ncid, ncv, tau_b_driving_conv) )
501  else
502  write(6,'(/1x,a)') '>>> read_erg_nc: Variable tau_b_driving'
503  write(6, '(1x,a)') ' not available in read file *.nc.'
504  tau_b_driving_conv = 0.0_sp
505  end if
506 
507  if ( nf90_inq_varid(ncid, 'tau_b_drag', ncv) == nf90_noerr ) then
508  call check( nf90_get_var(ncid, ncv, tau_b_drag_conv) )
509  else
510  write(6,'(/1x,a)') '>>> read_erg_nc: Variable tau_b_drag'
511  write(6, '(1x,a)') ' not available in read file *.nc.'
512  tau_b_drag_conv = 0.0_sp
513  end if
514 
515  call check( nf90_inq_varid(ncid, 'p_b_w', ncv) )
516  call check( nf90_get_var(ncid, ncv, p_b_w_conv) )
517 
518  call check( nf90_inq_varid(ncid, 'H_w', ncv) )
519  call check( nf90_get_var(ncid, ncv, h_w_conv) )
520 
521  call check( nf90_inq_varid(ncid, 'q_gl_g', ncv) )
522  call check( nf90_get_var(ncid, ncv, q_gl_g_conv) )
523 
524  if ( nf90_inq_varid(ncid, 'q_cf_g', ncv) == nf90_noerr ) then
525  call check( nf90_get_var(ncid, ncv, q_cf_g_conv) )
526  else
527  write(6,'(/1x,a)') '>>> read_erg_nc: Variable q_cf_g'
528  write(6, '(1x,a)') ' not available in read file *.nc.'
529  q_cf_g_conv = 0.0_sp
530  end if
531 
532  call check( nf90_inq_varid(ncid, 'vx_c', ncv) )
533  call check( nf90_get_var(ncid, ncv, vx_c_conv) )
534 
535  call check( nf90_inq_varid(ncid, 'vy_c', ncv) )
536  call check( nf90_get_var(ncid, ncv, vy_c_conv) )
537 
538  call check( nf90_inq_varid(ncid, 'vz_c', ncv) )
539  call check( nf90_get_var(ncid, ncv, vz_c_conv) )
540 
541  call check( nf90_inq_varid(ncid, 'vx_t', ncv) )
542  call check( nf90_get_var(ncid, ncv, vx_t_conv) )
543 
544  call check( nf90_inq_varid(ncid, 'vy_t', ncv) )
545  call check( nf90_get_var(ncid, ncv, vy_t_conv) )
546 
547  call check( nf90_inq_varid(ncid, 'vz_t', ncv) )
548  call check( nf90_get_var(ncid, ncv, vz_t_conv) )
549 
550  call check( nf90_inq_varid(ncid, 'temp_c', ncv) )
551  call check( nf90_get_var(ncid, ncv, temp_c_conv) )
552 
553  call check( nf90_inq_varid(ncid, 'omega_t', ncv) )
554  call check( nf90_get_var(ncid, ncv, omega_t_conv) )
555 
556  call check( nf90_inq_varid(ncid, 'temp_r', ncv) )
557  call check( nf90_get_var(ncid, ncv, temp_r_conv) )
558 
559  call check( nf90_inq_varid(ncid, 'enth_c', ncv) )
560  call check( nf90_get_var(ncid, ncv, enth_c_conv) )
561 
562  call check( nf90_inq_varid(ncid, 'enth_t', ncv) )
563  call check( nf90_get_var(ncid, ncv, enth_t_conv) )
564 
565  call check( nf90_inq_varid(ncid, 'omega_c', ncv) )
566  call check( nf90_get_var(ncid, ncv, omega_c_conv) )
567 
568  if ( nf90_inq_varid(ncid, 'enh_c', ncv) == nf90_noerr ) then
569  call check( nf90_get_var(ncid, ncv, enh_c_conv) )
570  else
571  write(6,'(/1x,a)') '>>> read_erg_nc: Variable enh_c'
572  write(6, '(1x,a)') ' not available in read file *.nc.'
573  enh_c_conv = 1.0_sp
574  end if
575 
576  if ( nf90_inq_varid(ncid, 'enh_t', ncv) == nf90_noerr ) then
577  call check( nf90_get_var(ncid, ncv, enh_t_conv) )
578  else
579  write(6,'(/1x,a)') '>>> read_erg_nc: Variable enh_t'
580  write(6, '(1x,a)') ' not available in read file *.nc.'
581  enh_t_conv = 1.0_sp
582  end if
583 
584  call check( nf90_inq_varid(ncid, 'age_c', ncv) )
585  call check( nf90_get_var(ncid, ncv, age_c_conv) )
586 
587  call check( nf90_inq_varid(ncid, 'age_t', ncv) )
588  call check( nf90_get_var(ncid, ncv, age_t_conv) )
589 
590  call check( nf90_close(ncid) )
591 
592 #endif
593 
594 !-------- Convert data to real*8 and years to seconds --------
595 
596  z_sl = real(z_sl_conv,dp)
597 
598  h_r = real(h_r_conv,dp)
599 
600  do i=0, imax
601  xi(i) = real(xi_conv(i),dp)
602  end do
603 
604  do j=0, jmax
605  eta(j) = real(eta_conv(j),dp)
606  end do
607 
608  do i=0, imax
609  do j=0, jmax
610 
611  maske(j,i) = maske_conv(i,j)
612  maske_old(j,i) = maske_old_conv(i,j)
613  n_cts(j,i) = n_cts_conv(i,j)
614  kc_cts(j,i) = kc_cts_conv(i,j)
615  zs(j,i) = real(zs_conv(i,j),dp)
616  zm(j,i) = real(zm_conv(i,j),dp)
617  zb(j,i) = real(zb_conv(i,j),dp)
618  zl(j,i) = real(zl_conv(i,j),dp)
619 #if (CALCMOD==1)
620  h_c(j,i) = real(H_cold_conv(i,j),dp)
621  h_t(j,i) = real(H_temp_conv(i,j),dp)
622 #elif (CALCMOD==0 || CALCMOD==2 || CALCMOD==3 || CALCMOD==-1)
623  h_c(j,i) = real(H_conv(i,j),dp)
624  h_t(j,i) = 0.0_dp
625 #endif
626  q_bm(j,i) = real(Q_bm_conv(i,j),dp)/year_sec
627  q_tld(j,i) = real(Q_tld_conv(i,j),dp)/year_sec
628  am_perp(j,i) = real(am_perp_conv(i,j),dp)/year_sec
629  qx(j,i) = real(qx_conv(i,j),dp)/year_sec
630  qy(j,i) = real(qy_conv(i,j),dp)/year_sec
631  dzs_dtau(j,i) = real(dzs_dtau_conv(i,j),dp)/year_sec
632  dzm_dtau(j,i) = real(dzm_dtau_conv(i,j),dp)/year_sec
633  dzb_dtau(j,i) = real(dzb_dtau_conv(i,j),dp)/year_sec
634  dzl_dtau(j,i) = real(dzl_dtau_conv(i,j),dp)/year_sec
635  dh_c_dtau(j,i) = real(dH_c_dtau_conv(i,j),dp)/year_sec
636  dh_t_dtau(j,i) = real(dH_t_dtau_conv(i,j),dp)/year_sec
637  vx_b_g(j,i) = real(vx_b_g_conv(i,j),dp)/year_sec
638  vy_b_g(j,i) = real(vy_b_g_conv(i,j),dp)/year_sec
639  vz_b(j,i) = real(vz_b_conv(i,j),dp)/year_sec
640  vx_s_g(j,i) = real(vx_s_g_conv(i,j),dp)/year_sec
641  vy_s_g(j,i) = real(vy_s_g_conv(i,j),dp)/year_sec
642  vz_s(j,i) = real(vz_s_conv(i,j),dp)/year_sec
643  temp_b(j,i) = real(temp_b_conv(i,j),dp)
644  temph_b(j,i) = real(temph_b_conv(i,j),dp)
645  p_b_w(j,i) = real(p_b_w_conv(i,j),dp)
646  h_w(j,i) = real(H_w_conv(i,j),dp)
647 
648  do kr=0, krmax
649  temp_r(kr,j,i) = real(temp_r_conv(i,j,kr),dp)
650  end do
651 
652  do kt=0, ktmax
653  vx_t(kt,j,i) = real(vx_t_conv(i,j,kt),dp)/year_sec
654  vy_t(kt,j,i) = real(vy_t_conv(i,j,kt),dp)/year_sec
655  vz_t(kt,j,i) = real(vz_t_conv(i,j,kt),dp)/year_sec
656  omega_t(kt,j,i) = real(omega_t_conv(i,j,kt),dp)
657  age_t(kt,j,i) = real(age_t_conv(i,j,kt),dp)*year_sec
658  enth_t(kt,j,i) = real(enth_t_conv(i,j,kt),dp)
659  enh_t(kt,j,i) = real(enh_t_conv(i,j,kt),dp)
660  end do
661 
662  do kc=0, kcmax
663  vx_c(kc,j,i) = real(vx_c_conv(i,j,kc),dp)/year_sec
664  vy_c(kc,j,i) = real(vy_c_conv(i,j,kc),dp)/year_sec
665  vz_c(kc,j,i) = real(vz_c_conv(i,j,kc),dp)/year_sec
666  temp_c(kc,j,i) = real(temp_c_conv(i,j,kc),dp)
667  age_c(kc,j,i) = real(age_c_conv(i,j,kc),dp)*year_sec
668  enth_c(kc,j,i) = real(enth_c_conv(i,j,kc),dp)
669  omega_c(kc,j,i) = real(omega_c_conv(i,j,kc),dp)
670  enh_c(kc,j,i) = real(enh_c_conv(i,j,kc),dp)
671  end do
672 
673  end do
674  end do
675 
676  end subroutine read_erg_nc
677 
678 !-------------------------------------------------------------------------------
679 !> Reading of the target-topography file (in NetCDF format).
680 !<------------------------------------------------------------------------------
681  subroutine read_target_topo_nc(target_topo_dat_name)
683  use sico_variables_m, only : maske_target, &
685  h_target
686 
687 #if (NETCDF==2)
688  use netcdf
689  use nc_check_m
690 #endif
691 
692  implicit none
693 
694  character(len=100), intent(in) :: target_topo_dat_name
695 
696 ! Return variables
697 ! (defined as global variables in module sico_variables_m):
698 !
699 ! maske_target, zs_target, zb_target, zl_target, H_target
700 
701  integer(i2b), dimension(0:IMAX,0:JMAX) :: maske_conv
702  real(sp), dimension(0:IMAX,0:JMAX) :: zs_conv, zb_conv, zl_conv, H_conv
703  character(len=256) :: target_topo_dat_path
704  character(len=256) :: filename_with_path
705 
706 #if (NETCDF==2)
707  integer(i4b) :: ncid, ncv
708  ! ncid: ID of the output file
709  ! ncv: Variable ID
710 #endif
711 
712  character(len=64), parameter :: thisroutine = 'read_target_topo_nc'
713 
714 !-------- Read target-topography data
715 ! (from a time-slice file of a previous simulation) --------
716 
717 #if (NETCDF==2)
718 
719 #if (defined(TARGET_TOPO_DAT_PATH))
720  target_topo_dat_path = trim(target_topo_dat_path)
721 #else
722  target_topo_dat_path = 'dummy'
723  stop ' >>> read_target_topo_nc: TARGET_TOPO_DAT_PATH must be defined!'
724 #endif
725 
726  filename_with_path &
727  = trim(target_topo_dat_path)//'/'//trim(target_topo_dat_name)
728 
729  call check( nf90_open(trim(filename_with_path), nf90_nowrite, ncid), &
730  thisroutine )
731 
732  call check( nf90_inq_varid(ncid, 'maske', ncv), thisroutine )
733  call check( nf90_get_var(ncid, ncv, maske_conv), thisroutine )
734 
735  call check( nf90_inq_varid(ncid, 'zs', ncv), thisroutine )
736  call check( nf90_get_var(ncid, ncv, zs_conv), thisroutine )
737 
738  call check( nf90_inq_varid(ncid, 'zb', ncv), thisroutine )
739  call check( nf90_get_var(ncid, ncv, zb_conv), thisroutine )
740 
741  call check( nf90_inq_varid(ncid, 'zl', ncv), thisroutine )
742  call check( nf90_get_var(ncid, ncv, zl_conv), thisroutine )
743 
744  call check( nf90_inq_varid(ncid, 'H', ncv), thisroutine )
745  call check( nf90_get_var(ncid, ncv, h_conv), thisroutine )
746 
747  call check( nf90_close(ncid), thisroutine )
748 
749 #else /* NETCDF != 2, not allowed */
750 
751  maske_conv = 0_i2b ! dummy values
752  zs_conv = 0.0_sp ! dummy values
753  zb_conv = 0.0_sp ! dummy values
754  zl_conv = 0.0_sp ! dummy values
755  h_conv = 0.0_sp ! dummy values
756 
757  stop ' >>> read_target_topo_nc: Parameter NETCDF must be set to 2!'
758 
759 #endif
760 
761 !-------- Convert data to double precision --------
762 
763  maske_target = transpose(maske_conv)
764  zs_target = real(transpose(zs_conv),dp)
765  zb_target = real(transpose(zb_conv),dp)
766  zl_target = real(transpose(zl_conv),dp)
767  h_target = real(transpose(H_conv) ,dp)
768 
769  end subroutine read_target_topo_nc
770 
771 !-------------------------------------------------------------------------------
772 !> Reading of the tabulated kei function.
773 !<------------------------------------------------------------------------------
774  subroutine read_kei()
777 
778  implicit none
779 
780  integer(i4b) :: n, n_data_kei_max
781  integer(i4b) :: ios
782  real(dp) :: r_val, d_dummy
783  character :: ch_dummy
784  character(len=256) :: filename_with_path
785 
786  n_data_kei_max = 10000
787 
788  kei = 0.0_dp
789 
790 !-------- Reading of tabulated values --------
791 
792  filename_with_path = trim(inpath)//'/general/kei.dat'
793 
794  open(unit=11, iostat=ios, file=trim(filename_with_path), status='old')
795 
796  if (ios /= 0) stop ' >>> read_kei: Error when opening the kei file!'
797 
798  read(unit=11, fmt='(a)') ch_dummy
799  read(unit=11, fmt='(15x,f5.0)') kei_r_max
800  read(unit=11, fmt='(15x,f5.0)') kei_r_incr
801  read(unit=11, fmt='(a)') ch_dummy
802 
804 
805  if (n_data_kei > n_data_kei_max) stop ' >>> read_kei: Array kei too small!'
806 
807  read(unit=11, fmt='(a)') ch_dummy
808  read(unit=11, fmt='(a)') ch_dummy
809  read(unit=11, fmt='(a)') ch_dummy
810 
811  do n=-n_data_kei, n_data_kei
812  read(unit=11, fmt=*) d_dummy, kei(n)
813  end do
814 
815  close(unit=11, status='keep')
816 
817  end subroutine read_kei
818 
819 !-------------------------------------------------------------------------------
820 !> Reading of physical parameters.
821 !<------------------------------------------------------------------------------
822  subroutine read_phys_para()
824  use sico_variables_m
825  use sico_vars_m
826 
827  implicit none
828 
829  integer(i4b), parameter :: n_unit=31
830  integer(i4b) :: ios
831  integer(i4b) :: n
832  character(len=256) :: filename_with_path
833 
834  filename_with_path = trim(inpath)//'/'//trim(ch_domain_short)//'/'// &
835  trim(phys_para_file)
836 
837  open(n_unit, iostat=ios, file=trim(filename_with_path), status='old')
838 
839  if (ios /= 0) &
840  stop ' >>> read_phys_para: Error when opening the phys_para file!'
841 
842 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */
843  call read_phys_para_value(rho, n_unit)
844 #else /* Martian polar caps */
845  call read_phys_para_value(rho_i, n_unit)
846 #endif
847  call read_phys_para_value(rho_w, n_unit)
848  call read_phys_para_value(rho_sw, n_unit)
849  call read_phys_para_value(l, n_unit)
850  call read_phys_para_value(g, n_unit)
851  call read_phys_para_value(nue, n_unit)
852  call read_phys_para_value(beta, n_unit)
853 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */
854  call read_phys_para_value(delta_tm_sw, n_unit)
855 #endif
856  call read_phys_para_value(omega_max, n_unit)
857 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
858  call read_phys_para_value(rho_c, n_unit)
859  call read_phys_para_value(kappa_c, n_unit)
860  call read_phys_para_value(c_c, n_unit)
861 #endif
862  call read_phys_para_value(h_r, n_unit)
863  call read_phys_para_value(rho_c_r, n_unit)
864  call read_phys_para_value(kappa_r, n_unit)
865  call read_phys_para_value(rho_a, n_unit)
866  call read_phys_para_value(r_t, n_unit)
867 
868  call read_phys_para_value(r, n_unit)
869  call read_phys_para_value(a, n_unit)
870  call read_phys_para_value(b, n_unit)
871  call read_phys_para_value(phi0, n_unit)
872  call read_phys_para_value(lambda0, n_unit)
873 
874 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */
875  call read_phys_para_value(s_stat_0, n_unit)
876 #if (!defined(GRL)) /* not Greenland */
877  call read_phys_para_value(beta1_0, n_unit)
878  call read_phys_para_value(beta2_0, n_unit)
879 #else /* Greenland */
880  call read_phys_para_value(beta1_lt_0, n_unit)
881  call read_phys_para_value(beta1_ht_0, n_unit)
882  call read_phys_para_value(beta2_lt_0, n_unit)
883  call read_phys_para_value(beta2_ht_0, n_unit)
884  call read_phys_para_value(phi_sep_0, n_unit)
885 #endif
886  call read_phys_para_value(pmax_0, n_unit)
887  call read_phys_para_value(mu_0, n_unit)
888 #endif
889 
890  do n=10, -190, -1
891  call read_phys_para_value(rf(n), n_unit)
892  end do
893 
894  do n=10, -190, -1
895  call read_phys_para_value(kappa(n), n_unit)
896  end do
897 
898  do n=10, -190, -1
899  call read_phys_para_value(c(n), n_unit)
900  end do
901 
902  close(n_unit, status='keep')
903 
904  end subroutine read_phys_para
905 
906 !-------------------------------------------------------------------------------
907 !> Reading of a value of a physical parameter from the phys_para file.
908 !<------------------------------------------------------------------------------
909  subroutine read_phys_para_value(para, n_unit)
911  implicit none
912 
913  integer(i4b), intent(in) :: n_unit
914  real(dp), intent(out) :: para
915 
916  character :: check, ch_dummy
917  character(len=*), parameter :: fmt1='(a)', fmt3='(15x)'
918  logical :: first_read
919 
920  first_read = .true.
921 
922  do
923  if (first_read) then
924  read(unit=n_unit, fmt=trim(fmt1), advance='no') check
925  first_read=.false.
926  else
927  read(unit=n_unit, fmt=trim(fmt1)) ch_dummy
928  read(unit=n_unit, fmt=trim(fmt1), advance='no') check
929  end if
930 
931  if (check /= '%') then ! no comment line
932  read(unit=n_unit, fmt=trim(fmt3), advance='no')
933  read(unit=n_unit, fmt=*) para
934  exit
935  end if
936  end do
937 
938  end subroutine read_phys_para_value
939 
940 !-------------------------------------------------------------------------------
941 
942 end module read_m
943 !
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) vx_c
vx_c(kc,j,i): Velocity in x-direction in the upper (kc) ice domain (at (i+1/2,j,kc)) ...
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) temp_c
temp_c(kc,j,i): Temperature in the upper (kc) ice domain
real(dp) kappa_r
KAPPA_R: Heat conductivity of the lithosphere.
real(dp) delta_tm_sw
DELTA_TM_SW: Melting point depression of sea water due to its average salinity.
subroutine, public read_erg_nc(z_sl, filename)
Reading of time-slice files in native binary or NetCDF format.
Definition: read_m.F90:50
real(dp) phi_sep_0
PHI_SEP_0: Separation latitude for the computation of the degree-day factors beta1 and beta2: Equator...
subroutine, public read_kei()
Reading of the tabulated kei function.
Definition: read_m.F90:775
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) age_c
age_c(kc,j,i): Age in the upper (kc) ice domain
integer(i2b), dimension(0:jmax, 0:imax) n_cts
n_cts(j,i): Mask for thermal conditions. -1: cold ice base, 0: temperate ice base with cold ice above...
real(dp) beta2_ht_0
BETA2_HT_0: Degree-day factor for ice at high summer temperatures.
real(dp) beta1_ht_0
BETA1_HT_0: Degree-day factor for snow at high summer temperatures.
real(dp), dimension(0:jmax, 0:imax) vx_b_g
(.)_g(j,i): Staggered-grid quantity (.) interpolated to grid point (i,j)
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) age_t
age_t(kt,j,i): Age in the lower (kt) ice domain
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) enh_c
enh_c(kc,j,i): Flow enhancement factor in the upper (kc) ice domain
real(dp) rho_w
RHO_W: Density of pure water.
real(dp), dimension(0:krmax, 0:jmax, 0:imax) temp_r
temp_r(kr,j,i): Temperature in the bedrock
real(dp) h_r
H_R: Thickness of the modelled lithosphere layer.
real(dp) a
A: Semi-major axis of the planet.
real(dp), dimension(0:jmax, 0:imax) vz_b
vz_b(j,i): Velocity in z-direction at the ice base, at (i,j)
real(dp) mu_0
MU_0: Firn-warming correction.
subroutine, public read_phys_para()
Reading of physical parameters.
Definition: read_m.F90:823
real(dp), dimension(0:jmax, 0:imax) vx_s_g
vx_s_g(j,i): Velocity in x-direction at the ice surface, at (i,j)
real(dp) omega_max
OMEGA_MAX: Threshold value for the water content of temperate ice.
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) enth_t
enth_t(kt,j,i): Enthalpy in the lower (kt) ice domain
real(dp), dimension(0:jmax, 0:imax) vz_s
vz_s(j,i): Velocity in z-direction at the ice surface, at (i,j)
real(dp), dimension(0:jmax, 0:imax) dh_t_dtau
dH_t_dtau(j,i): Derivative of H_t by tau (time)
real(dp), dimension(0:jmax, 0:imax) dh_c_dtau
dH_c_dtau(j,i): Derivative of H_c by tau (time)
integer(i2b), dimension(0:jmax, 0:imax) maske
maske(j,i): Ice-land-ocean mask. 0: grounded ice, 1: ice-free land, 2: ocean, 3: floating ice ...
real(dp), dimension(0:jmax, 0:imax) p_b_w
p_b_w(j,i): Basal water pressure
real(dp) s_stat_0
S_STAT_0: Standard deviation of the air termperature for the degree-day model.
real(dp), dimension(0:jmax, 0:imax) zl_target
zl_target(j,i): Target topography (lithosphere surface)
real(dp) beta2_lt_0
BETA2_LT_0: Degree-day factor for ice at low summer temperatures.
real(dp), dimension(0:jmax, 0:imax) zs_target
zs_target(j,i): Target topography (ice surface)
real(dp) b
B: Semi-minor axis of the planet.
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) vz_t
vz_t(kt,j,i): Velocity in z-direction in the lower (kt) ice domain (at (i,j,kt+1/2)) ...
real(dp), dimension(0:jmax, 0:imax) zl
zl(j,i): Coordinate z of the lithosphere surface
real(dp) g
G: Acceleration due to gravity.
real(dp) r
R: Radius of the planet.
integer(i2b), dimension(0:jmax, 0:imax) kc_cts
kc_cts(j,i): Position kc of the CTS (for COLD, ENTC, ENTM)
real(dp), dimension(0:jmax, 0:imax) q_tld
Q_tld(j,i): Water drainage rate from the temperate layer.
Declarations of global variables for SICOPOLIS (for the ANT domain).
Definition: sico_vars_m.F90:35
subroutine read_phys_para_value(para, n_unit)
Reading of a value of a physical parameter from the phys_para file.
Definition: read_m.F90:910
real(dp), dimension(0:jmax, 0:imax) temp_b
temp_b(j,i): Basal temperature
real(dp) beta1_0
BETA1_0: Degree-day factor for snow.
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) enth_c
enth_c(kc,j,i): Enthalpy in the upper (kc) ice domain
real(dp) rho_c_r
RHO_C_R: Density times specific heat of the lithosphere.
Reading of several input data.
Definition: read_m.F90:35
Declarations of kind types for SICOPOLIS.
real(dp), dimension(0:jmax, 0:imax) zs
zs(j,i): Coordinate z of the surface topography
real(dp) pmax_0
PMAX_0: Saturation factor for the formation of superimposed ice.
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) enh_t
enh_t(kt,j,i): Flow enhancement factor in the lower (kt) ice domain
real(dp) rho_c
RHO_C: Density of crustal material (dust)
Definition: sico_vars_m.F90:79
real(dp), dimension(0:jmax, 0:imax) h_c
H_c(j,i): Thickness of ice in the upper (kc) domain (thickness of the cold-ice layer for POLY...
real(dp), dimension(0:jmax, 0:imax) zb_target
zb_target(j,i): Target topography (ice base)
real(dp) kei_r_max
kei_r_max: Maximum value of the argument r of the tabulated kei function
real(dp), dimension(0:jmax, 0:imax) vy_b_g
(.)_g(j,i): Staggered-grid quantity (.) interpolated to grid point (i,j)
real(dp) beta2_0
BETA2_0: Degree-day factor for ice.
real(dp), dimension(-190:10) kappa
KAPPA(n): Tabulated values for the heat conductivity of ice.
real(dp), dimension(0:jmax, 0:imax) dzm_dtau
dzm_dtau(j,i): Derivative of zm by tau (time)
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) omega_c
omega_c(kc,j,i): Water content in the upper (kc) ice domain
integer(i2b), dimension(0:jmax, 0:imax) maske_old
maske_old(j,i): Old value of maske (at the previous time step)
real(dp), dimension(0:jmax, 0:imax) qx
qx(j,i): Volume flux in x-direction (depth-integrated vx, at (i+1/2,j))
real(dp) kappa_c
KAPPA_C: Heat conductivity of crustal material (dust)
Definition: sico_vars_m.F90:81
real(dp) rho_i
RHO_I: Density of ice.
Definition: sico_vars_m.F90:77
real(dp) nue
NUE: Water diffusivity in ice.
real(dp) c_c
C_C: Specific heat of crustal material (dust)
Definition: sico_vars_m.F90:83
real(dp) beta1_lt_0
BETA1_LT_0: Degree-day factor for snow at low summer temperatures.
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) vz_c
vz_c(kc,j,i): Velocity in z-direction in the upper (kc) ice domain (at (i,j,kc+1/2)) ...
real(dp) phi0
PHI0: Standard parallel of the stereographic projection.
real(dp), dimension(0:jmax, 0:imax) vy_s_g
vy_s_g(j,i): Velocity in x-direction at the ice surface, at (i,j)
NetCDF error capturing.
Definition: nc_check_m.F90:35
real(dp), dimension(0:jmax, 0:imax) zm
zm(j,i): Coordinate z of the bottom of the upper (kc) ice domain = top of the lower (kt) ice domain (...
real(dp), dimension(-190:10) rf
RF(n): Tabulated values for the rate factor of cold ice.
real(dp), dimension(0:jmax, 0:imax) q_bm
Q_bm(j,i): Basal melting rate.
real(dp), dimension(-190:10) c
C(n): Tabulated values for the specific heat of ice.
real(dp), dimension(-10000:10000) kei
kei(n): Tabulated values of the kei function (Kelvin function of zero order)
subroutine check(status, ch_calling_routine)
NetCDF error capturing.
Definition: nc_check_m.F90:46
real(dp), dimension(0:jmax, 0:imax) temph_b
temph_b(j,i): Basal temperature relative to the pressure melting point
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) vy_t
vy_t(kt,j,i): Velocity in y-direction in the lower (kt) ice domain (at (i,j+1/2,kt)) ...
real(dp), dimension(0:jmax, 0:imax) am_perp
am_perp(j,i): Ice volume flux across the z=zm interface
real(dp) beta
BETA: Clausius-Clapeyron gradient of ice.
real(dp), parameter no_value_pos_1
no_value_pos_1: Positive no-value parameter
real(dp), dimension(0:jmax) eta
eta(j): Coordinate eta (= y) of grid point j
real(dp), dimension(0:imax) xi
xi(i): Coordinate xi (= x) of grid point i
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) vy_c
vy_c(kc,j,i): Velocity in y-direction in the upper (kc) ice domain (at (i,j+1/2,kc)) ...
character(len=16) ch_domain_short
ch_domain_short: Short name of the computational domain
real(dp) r_t
R_T: Coefficient of the water-content dependence in the rate factor for temperate ice...
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) omega_t
omega_t(kt,j,i): Water content in the lower (kt) ice domain
real(dp) kei_r_incr
kei_r_incr: Increment of the argument r of the tabulated kei function
real(dp) rho_sw
RHO_SW: Density of sea water.
real(dp) l
L: Latent heat of ice.
integer, parameter sp
Single-precision reals.
integer(i4b) n_data_kei
n_data_kei: Number of tabulated values of the kei function
real(dp) rho
RHO: Density of ice.
real(dp) rho_a
RHO_A: Density of the asthenosphere.
real(dp), dimension(0:jmax, 0:imax) h_target
H_target(j,i): Target topography (ice thickness)
integer, parameter dp
Double-precision reals.
real(dp), dimension(0:jmax, 0:imax) zb
zb(j,i): Coordinate z of the ice base
real(dp), dimension(0:jmax, 0:imax) dzb_dtau
dzb_dtau(j,i): Derivative of zb by tau (time)
integer(i2b), dimension(0:jmax, 0:imax) maske_target
maske_target(j,i): Target topography (ice-land-ocean mask)
subroutine, public read_target_topo_nc(target_topo_dat_name)
Reading of the target-topography file (in NetCDF format).
Definition: read_m.F90:682
real(dp) lambda0
LAMBDA0: Reference longitude (central meridian) of the stereographic projection.
real(dp), dimension(0:jmax, 0:imax) h_t
H_t(j,i): Thickness of ice in the lower (kt) domain (thickness of the temperate layer for POLY...
real(dp), dimension(0:jmax, 0:imax) h_w
H_w(j,i): Thickness of the water column under the ice base.
Declarations of global variables for SICOPOLIS.
real(dp), dimension(0:jmax, 0:imax) dzs_dtau
dzs_dtau(j,i): Derivative of zs by tau (time)
real(dp), dimension(0:jmax, 0:imax) dzl_dtau
dzl_dtau(j,i): Derivative of zl by tau (time)
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) vx_t
vx_t(kt,j,i): Velocity in x-direction in the lower (kt) ice domain (at (i+1/2,j,kt)) ...
real(dp), dimension(0:jmax, 0:imax) qy
qy(j,i): Volume flux in y-direction (depth-integrated vy, at (i,j+1/2))