SICOPOLIS V5-dev  Revision 1105
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, as_perp_conv, as_perp_apl_conv, &
95  q_geo_conv, &
96  zs_conv, zm_conv, zb_conv, zl_conv, zl0_conv, &
97  H_cold_conv, H_temp_conv, H_conv, &
98  Q_bm_conv, Q_tld_conv, &
99  am_perp_conv, &
100  qx_conv, qy_conv, &
101  dzs_dtau_conv, dzm_dtau_conv, dzb_dtau_conv, dzl_dtau_conv, &
102  dH_c_dtau_conv, dH_t_dtau_conv, dH_dtau_conv, &
103  vx_b_g_conv, vy_b_g_conv, vz_b_conv, vh_b_conv, &
104  vx_s_g_conv, vy_s_g_conv, vz_s_conv, vh_s_conv, &
105  vx_m_g_conv, vy_m_g_conv, vh_m_conv, &
106  temp_b_conv, temph_b_conv, &
107  tau_b_driving_conv, tau_b_drag_conv, &
108  p_b_w_conv, H_w_conv, q_gl_g_conv, q_cf_g_conv
109  real(sp), dimension(0:IMAX,0:JMAX,0:KCMAX) :: vx_c_conv, vy_c_conv, vz_c_conv, &
110  temp_c_conv, age_c_conv, &
111  enth_c_conv, omega_c_conv, &
112  enh_c_conv
113  real(sp), dimension(0:IMAX,0:JMAX,0:KTMAX) :: vx_t_conv, vy_t_conv, vz_t_conv, &
114  omega_t_conv, age_t_conv, &
115  enth_t_conv, &
116  enh_t_conv
117  real(sp), dimension(0:IMAX,0:JMAX,0:KRMAX) :: temp_r_conv
118 
119 #if (DISC>0) /* Ice discharge parameterisation */
120  integer(i2b), dimension(0:IMAX,0:JMAX) :: mask_mar_conv
121  real(sp), dimension(0:IMAX,0:JMAX) :: dis_perp_conv, &
122  cst_dist_conv, cos_grad_tc_conv
123 #endif
124 
125 !-------- Read data from time-slice file of previous simulation --------
126 
127  filename_with_path = trim(anfdatpath)//'/'//trim(filename)
128 
129 #if (NETCDF==1) /* time-slice file in native binary format */
130 
131  open(unit=11, iostat=ios, file=trim(filename_with_path), status='old', &
132  form='unformatted')
133 
134  if (ios /= 0) stop ' >>> read_erg_nc: Error when opening the time-slice file!'
135 
136  read(unit=11) ch_attr_title
137  read(unit=11) ch_attr_institution
138  read(unit=11) ch_attr_source
139  read(unit=11) ch_attr_history
140  read(unit=11) ch_attr_references
141 
142  read(unit=11) time_conv
143  read(unit=11) dummy_conv ! this is either delta_ts or glac_index; not needed
144  read(unit=11) z_sl_conv
145 
146  read(unit=11) dummy_conv ! this is V_tot; not needed
147  read(unit=11) dummy_conv ! this is V_af; not needed
148  read(unit=11) dummy_conv ! this is A_grounded; not needed
149  read(unit=11) dummy_conv ! this is A_floating; not needed
150 
151  read(unit=11) xi_conv
152  read(unit=11) eta_conv
153  read(unit=11) sigma_level_c_conv
154  read(unit=11) sigma_level_t_conv
155  read(unit=11) sigma_level_r_conv
156  read(unit=11) lon_conv
157  read(unit=11) lat_conv
158  read(unit=11) lambda_conv
159  read(unit=11) phi_conv
160  read(unit=11) temp_s_conv
161  read(unit=11) accum_conv
162  read(unit=11) as_perp_conv
163  read(unit=11) as_perp_apl_conv
164 
165 #if (DISC>0) /* Ice discharge parameterisation */
166 
167  read(unit=11) dis_perp_conv
168  read(unit=11) cst_dist_conv
169  read(unit=11) cos_grad_tc_conv
170  read(unit=11) mask_mar_conv
171 
172 #endif
173 
174  read(unit=11) q_geo_conv
175  read(unit=11) maske_conv
176  read(unit=11) maske_old_conv
177  read(unit=11) n_cts_conv
178  read(unit=11) kc_cts_conv
179  read(unit=11) zs_conv
180  read(unit=11) zm_conv
181  read(unit=11) zb_conv
182  read(unit=11) zl_conv
183  read(unit=11) zl0_conv
184  read(unit=11) h_cold_conv
185  read(unit=11) h_temp_conv
186  read(unit=11) h_conv
187  read(unit=11) h_r_conv
188  read(unit=11) q_bm_conv
189  read(unit=11) q_tld_conv
190  read(unit=11) am_perp_conv
191  read(unit=11) qx_conv
192  read(unit=11) qy_conv
193  read(unit=11) dzs_dtau_conv
194  read(unit=11) dzm_dtau_conv
195  read(unit=11) dzb_dtau_conv
196  read(unit=11) dzl_dtau_conv
197  read(unit=11) dh_c_dtau_conv
198  read(unit=11) dh_t_dtau_conv
199  read(unit=11) dh_dtau_conv
200  read(unit=11) vx_b_g_conv
201  read(unit=11) vy_b_g_conv
202  read(unit=11) vz_b_conv
203  read(unit=11) vh_b_conv
204  read(unit=11) vx_s_g_conv
205  read(unit=11) vy_s_g_conv
206  read(unit=11) vz_s_conv
207  read(unit=11) vh_s_conv
208  read(unit=11) vx_m_g_conv
209  read(unit=11) vy_m_g_conv
210  read(unit=11) vh_m_conv
211  read(unit=11) temp_b_conv
212  read(unit=11) temph_b_conv
213  read(unit=11) tau_b_driving_conv
214  read(unit=11) tau_b_drag_conv
215  read(unit=11) p_b_w_conv
216  read(unit=11) h_w_conv
217  read(unit=11) q_gl_g_conv
218  read(unit=11) q_cf_g_conv
219 
220  read(unit=11) vx_c_conv
221  read(unit=11) vy_c_conv
222  read(unit=11) vz_c_conv
223  read(unit=11) vx_t_conv
224  read(unit=11) vy_t_conv
225  read(unit=11) vz_t_conv
226  read(unit=11) temp_c_conv
227  read(unit=11) omega_t_conv
228  read(unit=11) temp_r_conv
229  read(unit=11) enth_c_conv
230  read(unit=11) enth_t_conv
231  read(unit=11) omega_c_conv
232  read(unit=11) enh_c_conv
233  read(unit=11) enh_t_conv
234  read(unit=11) age_c_conv
235  read(unit=11) age_t_conv
236 
237  close(unit=11,status='keep')
238 
239 #elif (NETCDF==2) /* time-slice file in NetCDF format */
240 
241  call check( nf90_open(trim(filename_with_path), nf90_nowrite, ncid) )
242 
243  call check( nf90_inq_varid(ncid, 'time', ncv) )
244  call check( nf90_get_var(ncid, ncv, time_conv) )
245 
246  if ( nf90_inq_varid(ncid, 'delta_ts', ncv) == nf90_noerr ) then
247  call check( nf90_get_var(ncid, ncv, dummy_conv) )
248  else if ( nf90_inq_varid(ncid, 'glac_index', ncv) == nf90_noerr ) then
249  call check( nf90_get_var(ncid, ncv, dummy_conv) )
250  else
251  dummy_conv = 0.0_sp
252  end if
253 
254  call check( nf90_inq_varid(ncid, 'z_sl', ncv) )
255  call check( nf90_get_var(ncid, ncv, z_sl_conv) )
256 
257  call check( nf90_inq_varid(ncid, 'x', ncv) )
258  call check( nf90_get_var(ncid, ncv, xi_conv) )
259 
260  call check( nf90_inq_varid(ncid, 'y', ncv) )
261  call check( nf90_get_var(ncid, ncv, eta_conv) )
262 
263  call check( nf90_inq_varid(ncid, 'sigma_level_c', ncv) )
264  call check( nf90_get_var(ncid, ncv, sigma_level_c_conv) )
265 
266  call check( nf90_inq_varid(ncid, 'sigma_level_t', ncv) )
267  call check( nf90_get_var(ncid, ncv, sigma_level_t_conv) )
268 
269  call check( nf90_inq_varid(ncid, 'sigma_level_r', ncv) )
270  call check( nf90_get_var(ncid, ncv, sigma_level_r_conv) )
271 
272  call check( nf90_inq_varid(ncid, 'lon', ncv) )
273  call check( nf90_get_var(ncid, ncv, lon_conv) )
274 
275  call check( nf90_inq_varid(ncid, 'lat', ncv) )
276  call check( nf90_get_var(ncid, ncv, lat_conv) )
277 
278  call check( nf90_inq_varid(ncid, 'lambda', ncv) )
279  call check( nf90_get_var(ncid, ncv, lambda_conv) )
280 
281  call check( nf90_inq_varid(ncid, 'phi', ncv) )
282  call check( nf90_get_var(ncid, ncv, phi_conv) )
283 
284  call check( nf90_inq_varid(ncid, 'temp_s', ncv) )
285  call check( nf90_get_var(ncid, ncv, temp_s_conv) )
286 
287  if ( nf90_inq_varid(ncid, 'prec', ncv) == nf90_noerr ) then
288  call check( nf90_get_var(ncid, ncv, accum_conv) )
289  else
290  write(6,'(/1x,a)') '>>> read_erg_nc: Variable prec'
291  write(6, '(1x,a)') ' not available in read file *.nc.'
292  accum_conv = 0.0_sp
293  end if
294 
295  call check( nf90_inq_varid(ncid, 'as_perp', ncv) )
296  call check( nf90_get_var(ncid, ncv, as_perp_conv) )
297 
298  if ( nf90_inq_varid(ncid, 'as_perp_apl', ncv) == nf90_noerr ) then
299  call check( nf90_get_var(ncid, ncv, as_perp_apl_conv) )
300  else
301  write(6,'(/1x,a)') '>>> read_erg_nc: Variable as_perp_apl'
302  write(6, '(1x,a)') ' not available in read file *.nc.'
303  as_perp_apl_conv = 0.0_sp
304  end if
305 
306 #if (DISC>0) /* Ice discharge parameterisation */
307 
308  if ( nf90_inq_varid(ncid, 'dis_perp', ncv) == nf90_noerr ) then
309  call check( nf90_get_var(ncid, ncv, dis_perp_conv) )
310  else
311  write(6,'(/1x,a)') '>>> read_erg_nc: Variable dis_perp'
312  write(6, '(1x,a)') ' not available in read file *.nc.'
313  dis_perp_conv = 0.0_sp
314  end if
315 
316  if ( nf90_inq_varid(ncid, 'cst_dist', ncv) == nf90_noerr ) then
317  call check( nf90_get_var(ncid, ncv, cst_dist_conv) )
318  else
319  write(6,'(/1x,a)') '>>> read_erg_nc: Variable cst_dist'
320  write(6, '(1x,a)') ' not available in read file *.nc.'
321  cst_dist_conv = real(no_value_pos_1,sp)
322  end if
323 
324  if ( nf90_inq_varid(ncid, 'cos_grad_tc', ncv) == nf90_noerr ) then
325  call check( nf90_get_var(ncid, ncv, cos_grad_tc_conv) )
326  else
327  write(6,'(/1x,a)') '>>> read_erg_nc: Variable cos_grad_tc'
328  write(6, '(1x,a)') ' not available in read file *.nc.'
329  cos_grad_tc_conv = 1.0_sp
330  end if
331 
332  if ( nf90_inq_varid(ncid, 'mask_mar', ncv) == nf90_noerr ) then
333  call check( nf90_get_var(ncid, ncv, mask_mar_conv) )
334  else
335  write(6,'(/1x,a)') '>>> read_erg_nc: Variable mask_mar'
336  write(6, '(1x,a)') ' not available in read file *.nc.'
337  mask_mar_conv = 0_i2b
338  end if
339 
340 #endif
341 
342  if ( nf90_inq_varid(ncid, 'q_geo', ncv) == nf90_noerr ) then
343  call check( nf90_get_var(ncid, ncv, q_geo_conv) )
344  else
345  write(6,'(/1x,a)') '>>> read_erg_nc: Variable q_geo'
346  write(6, '(1x,a)') ' not available in read file *.nc.'
347  q_geo_conv = 0.0_sp
348  end if
349 
350  call check( nf90_inq_varid(ncid, 'maske', ncv) )
351  call check( nf90_get_var(ncid, ncv, maske_conv) )
352 
353  if ( nf90_inq_varid(ncid, 'maske_old', ncv) == nf90_noerr ) then
354  call check( nf90_get_var(ncid, ncv, maske_old_conv) )
355  else
356  write(6,'(/1x,a)') '>>> read_erg_nc: Variable maske_old'
357  write(6, '(1x,a)') ' not available in read file *.nc.'
358  maske_old_conv = maske_conv
359  end if
360 
361  call check( nf90_inq_varid(ncid, 'n_cts', ncv) )
362  call check( nf90_get_var(ncid, ncv, n_cts_conv) )
363 
364  call check( nf90_inq_varid(ncid, 'kc_cts', ncv) )
365  call check( nf90_get_var(ncid, ncv, kc_cts_conv) )
366 
367  call check( nf90_inq_varid(ncid, 'zs', ncv) )
368  call check( nf90_get_var(ncid, ncv, zs_conv) )
369 
370  call check( nf90_inq_varid(ncid, 'zm', ncv) )
371  call check( nf90_get_var(ncid, ncv, zm_conv) )
372 
373  call check( nf90_inq_varid(ncid, 'zb', ncv) )
374  call check( nf90_get_var(ncid, ncv, zb_conv) )
375 
376  call check( nf90_inq_varid(ncid, 'zl', ncv) )
377  call check( nf90_get_var(ncid, ncv, zl_conv) )
378 
379  if ( nf90_inq_varid(ncid, 'zl0', ncv) == nf90_noerr ) then
380  call check( nf90_get_var(ncid, ncv, zl0_conv) )
381  else
382  write(6,'(/1x,a)') '>>> read_erg_nc: Variable zl0'
383  write(6, '(1x,a)') ' not available in read file *.nc.'
384  zl0_conv = 0.0_sp
385  end if
386 
387  call check( nf90_inq_varid(ncid, 'H_cold', ncv) )
388  call check( nf90_get_var(ncid, ncv, h_cold_conv) )
389 
390  call check( nf90_inq_varid(ncid, 'H_temp', ncv) )
391  call check( nf90_get_var(ncid, ncv, h_temp_conv) )
392 
393  call check( nf90_inq_varid(ncid, 'H', ncv) )
394  call check( nf90_get_var(ncid, ncv, h_conv) )
395 
396  call check( nf90_inq_varid(ncid, 'H_R', ncv) )
397  call check( nf90_get_var(ncid, ncv, h_r_conv) )
398 
399  call check( nf90_inq_varid(ncid, 'Q_bm', ncv) )
400  call check( nf90_get_var(ncid, ncv, q_bm_conv) )
401 
402  call check( nf90_inq_varid(ncid, 'Q_tld', ncv) )
403  call check( nf90_get_var(ncid, ncv, q_tld_conv) )
404 
405  call check( nf90_inq_varid(ncid, 'am_perp', ncv) )
406  call check( nf90_get_var(ncid, ncv, am_perp_conv) )
407 
408  call check( nf90_inq_varid(ncid, 'qx', ncv) )
409  call check( nf90_get_var(ncid, ncv, qx_conv) )
410 
411  call check( nf90_inq_varid(ncid, 'qy', ncv) )
412  call check( nf90_get_var(ncid, ncv, qy_conv) )
413 
414  call check( nf90_inq_varid(ncid, 'dzs_dt', ncv) )
415  call check( nf90_get_var(ncid, ncv, dzs_dtau_conv) )
416 
417  call check( nf90_inq_varid(ncid, 'dzm_dt', ncv) )
418  call check( nf90_get_var(ncid, ncv, dzm_dtau_conv) )
419 
420  call check( nf90_inq_varid(ncid, 'dzb_dt', ncv) )
421  call check( nf90_get_var(ncid, ncv, dzb_dtau_conv) )
422 
423  call check( nf90_inq_varid(ncid, 'dzl_dt', ncv) )
424  call check( nf90_get_var(ncid, ncv, dzl_dtau_conv) )
425 
426  call check( nf90_inq_varid(ncid, 'dH_c_dt', ncv) )
427  call check( nf90_get_var(ncid, ncv, dh_c_dtau_conv) )
428 
429  call check( nf90_inq_varid(ncid, 'dH_t_dt', ncv) )
430  call check( nf90_get_var(ncid, ncv, dh_t_dtau_conv) )
431 
432  call check( nf90_inq_varid(ncid, 'dH_dt', ncv) )
433  call check( nf90_get_var(ncid, ncv, dh_dtau_conv) )
434 
435  call check( nf90_inq_varid(ncid, 'vx_b_g', ncv) )
436  call check( nf90_get_var(ncid, ncv, vx_b_g_conv) )
437 
438  call check( nf90_inq_varid(ncid, 'vy_b_g', ncv) )
439  call check( nf90_get_var(ncid, ncv, vy_b_g_conv) )
440 
441  call check( nf90_inq_varid(ncid, 'vz_b', ncv) )
442  call check( nf90_get_var(ncid, ncv, vz_b_conv) )
443 
444  call check( nf90_inq_varid(ncid, 'vh_b', ncv) )
445  call check( nf90_get_var(ncid, ncv, vh_b_conv) )
446 
447  call check( nf90_inq_varid(ncid, 'vx_s_g', ncv) )
448  call check( nf90_get_var(ncid, ncv, vx_s_g_conv) )
449 
450  call check( nf90_inq_varid(ncid, 'vy_s_g', ncv) )
451  call check( nf90_get_var(ncid, ncv, vy_s_g_conv) )
452 
453  call check( nf90_inq_varid(ncid, 'vz_s', ncv) )
454  call check( nf90_get_var(ncid, ncv, vz_s_conv) )
455 
456  call check( nf90_inq_varid(ncid, 'vh_s', ncv) )
457  call check( nf90_get_var(ncid, ncv, vh_s_conv) )
458 
459  if ( nf90_inq_varid(ncid, 'vx_m_g', ncv) == nf90_noerr ) then
460  call check( nf90_get_var(ncid, ncv, vx_m_g_conv) )
461  else
462  write(6,'(/1x,a)') '>>> read_erg_nc: Variable vx_m_g'
463  write(6, '(1x,a)') ' not available in read file *.nc.'
464  vx_m_g_conv = 0.0_sp
465  end if
466 
467  if ( nf90_inq_varid(ncid, 'vy_m_g', ncv) == nf90_noerr ) then
468  call check( nf90_get_var(ncid, ncv, vy_m_g_conv) )
469  else
470  write(6,'(/1x,a)') '>>> read_erg_nc: Variable vy_m_g'
471  write(6, '(1x,a)') ' not available in read file *.nc.'
472  vy_m_g_conv = 0.0_sp
473  end if
474 
475  if ( nf90_inq_varid(ncid, 'vh_m', ncv) == nf90_noerr ) then
476  call check( nf90_get_var(ncid, ncv, vh_m_conv) )
477  else
478  write(6,'(/1x,a)') '>>> read_erg_nc: Variable vh_m'
479  write(6, '(1x,a)') ' not available in read file *.nc.'
480  vh_m_conv = 0.0_sp
481  end if
482 
483  call check( nf90_inq_varid(ncid, 'temp_b', ncv) )
484  call check( nf90_get_var(ncid, ncv, temp_b_conv) )
485 
486  call check( nf90_inq_varid(ncid, 'temph_b', ncv) )
487  call check( nf90_get_var(ncid, ncv, temph_b_conv) )
488 
489  if ( nf90_inq_varid(ncid, 'tau_b_driving', ncv) == nf90_noerr ) then
490  call check( nf90_get_var(ncid, ncv, tau_b_driving_conv) )
491  else
492  write(6,'(/1x,a)') '>>> read_erg_nc: Variable tau_b_driving'
493  write(6, '(1x,a)') ' not available in read file *.nc.'
494  tau_b_driving_conv = 0.0_sp
495  end if
496 
497  if ( nf90_inq_varid(ncid, 'tau_b_drag', ncv) == nf90_noerr ) then
498  call check( nf90_get_var(ncid, ncv, tau_b_drag_conv) )
499  else
500  write(6,'(/1x,a)') '>>> read_erg_nc: Variable tau_b_drag'
501  write(6, '(1x,a)') ' not available in read file *.nc.'
502  tau_b_drag_conv = 0.0_sp
503  end if
504 
505  call check( nf90_inq_varid(ncid, 'p_b_w', ncv) )
506  call check( nf90_get_var(ncid, ncv, p_b_w_conv) )
507 
508  call check( nf90_inq_varid(ncid, 'H_w', ncv) )
509  call check( nf90_get_var(ncid, ncv, h_w_conv) )
510 
511  call check( nf90_inq_varid(ncid, 'q_gl_g', ncv) )
512  call check( nf90_get_var(ncid, ncv, q_gl_g_conv) )
513 
514  if ( nf90_inq_varid(ncid, 'q_cf_g', ncv) == nf90_noerr ) then
515  call check( nf90_get_var(ncid, ncv, q_cf_g_conv) )
516  else
517  write(6,'(/1x,a)') '>>> read_erg_nc: Variable q_cf_g'
518  write(6, '(1x,a)') ' not available in read file *.nc.'
519  q_cf_g_conv = 0.0_sp
520  end if
521 
522  call check( nf90_inq_varid(ncid, 'vx_c', ncv) )
523  call check( nf90_get_var(ncid, ncv, vx_c_conv) )
524 
525  call check( nf90_inq_varid(ncid, 'vy_c', ncv) )
526  call check( nf90_get_var(ncid, ncv, vy_c_conv) )
527 
528  call check( nf90_inq_varid(ncid, 'vz_c', ncv) )
529  call check( nf90_get_var(ncid, ncv, vz_c_conv) )
530 
531  call check( nf90_inq_varid(ncid, 'vx_t', ncv) )
532  call check( nf90_get_var(ncid, ncv, vx_t_conv) )
533 
534  call check( nf90_inq_varid(ncid, 'vy_t', ncv) )
535  call check( nf90_get_var(ncid, ncv, vy_t_conv) )
536 
537  call check( nf90_inq_varid(ncid, 'vz_t', ncv) )
538  call check( nf90_get_var(ncid, ncv, vz_t_conv) )
539 
540  call check( nf90_inq_varid(ncid, 'temp_c', ncv) )
541  call check( nf90_get_var(ncid, ncv, temp_c_conv) )
542 
543  call check( nf90_inq_varid(ncid, 'omega_t', ncv) )
544  call check( nf90_get_var(ncid, ncv, omega_t_conv) )
545 
546  call check( nf90_inq_varid(ncid, 'temp_r', ncv) )
547  call check( nf90_get_var(ncid, ncv, temp_r_conv) )
548 
549  call check( nf90_inq_varid(ncid, 'enth_c', ncv) )
550  call check( nf90_get_var(ncid, ncv, enth_c_conv) )
551 
552  call check( nf90_inq_varid(ncid, 'enth_t', ncv) )
553  call check( nf90_get_var(ncid, ncv, enth_t_conv) )
554 
555  call check( nf90_inq_varid(ncid, 'omega_c', ncv) )
556  call check( nf90_get_var(ncid, ncv, omega_c_conv) )
557 
558  if ( nf90_inq_varid(ncid, 'enh_c', ncv) == nf90_noerr ) then
559  call check( nf90_get_var(ncid, ncv, enh_c_conv) )
560  else
561  write(6,'(/1x,a)') '>>> read_erg_nc: Variable enh_c'
562  write(6, '(1x,a)') ' not available in read file *.nc.'
563  enh_c_conv = 1.0_sp
564  end if
565 
566  if ( nf90_inq_varid(ncid, 'enh_t', ncv) == nf90_noerr ) then
567  call check( nf90_get_var(ncid, ncv, enh_t_conv) )
568  else
569  write(6,'(/1x,a)') '>>> read_erg_nc: Variable enh_t'
570  write(6, '(1x,a)') ' not available in read file *.nc.'
571  enh_t_conv = 1.0_sp
572  end if
573 
574  call check( nf90_inq_varid(ncid, 'age_c', ncv) )
575  call check( nf90_get_var(ncid, ncv, age_c_conv) )
576 
577  call check( nf90_inq_varid(ncid, 'age_t', ncv) )
578  call check( nf90_get_var(ncid, ncv, age_t_conv) )
579 
580  call check( nf90_close(ncid) )
581 
582 #endif
583 
584 !-------- Convert data to real*8 and years to seconds --------
585 
586  z_sl = real(z_sl_conv,dp)
587 
588  h_r = real(h_r_conv,dp)
589 
590  do i=0, imax
591  xi(i) = real(xi_conv(i),dp)
592  end do
593 
594  do j=0, jmax
595  eta(j) = real(eta_conv(j),dp)
596  end do
597 
598  do i=0, imax
599  do j=0, jmax
600 
601  maske(j,i) = maske_conv(i,j)
602  maske_old(j,i) = maske_old_conv(i,j)
603  n_cts(j,i) = n_cts_conv(i,j)
604  kc_cts(j,i) = kc_cts_conv(i,j)
605  zs(j,i) = real(zs_conv(i,j),dp)
606  zm(j,i) = real(zm_conv(i,j),dp)
607  zb(j,i) = real(zb_conv(i,j),dp)
608  zl(j,i) = real(zl_conv(i,j),dp)
609 #if (CALCMOD==1)
610  h_c(j,i) = real(H_cold_conv(i,j),dp)
611  h_t(j,i) = real(H_temp_conv(i,j),dp)
612 #elif (CALCMOD==0 || CALCMOD==2 || CALCMOD==3 || CALCMOD==-1)
613  h_c(j,i) = real(H_conv(i,j),dp)
614  h_t(j,i) = 0.0_dp
615 #endif
616  q_bm(j,i) = real(Q_bm_conv(i,j),dp)/year_sec
617  q_tld(j,i) = real(Q_tld_conv(i,j),dp)/year_sec
618  am_perp(j,i) = real(am_perp_conv(i,j),dp)/year_sec
619  qx(j,i) = real(qx_conv(i,j),dp)/year_sec
620  qy(j,i) = real(qy_conv(i,j),dp)/year_sec
621  dzs_dtau(j,i) = real(dzs_dtau_conv(i,j),dp)/year_sec
622  dzm_dtau(j,i) = real(dzm_dtau_conv(i,j),dp)/year_sec
623  dzb_dtau(j,i) = real(dzb_dtau_conv(i,j),dp)/year_sec
624  dzl_dtau(j,i) = real(dzl_dtau_conv(i,j),dp)/year_sec
625  dh_c_dtau(j,i) = real(dH_c_dtau_conv(i,j),dp)/year_sec
626  dh_t_dtau(j,i) = real(dH_t_dtau_conv(i,j),dp)/year_sec
627  vx_b_g(j,i) = real(vx_b_g_conv(i,j),dp)/year_sec
628  vy_b_g(j,i) = real(vy_b_g_conv(i,j),dp)/year_sec
629  vz_b(j,i) = real(vz_b_conv(i,j),dp)/year_sec
630  vx_s_g(j,i) = real(vx_s_g_conv(i,j),dp)/year_sec
631  vy_s_g(j,i) = real(vy_s_g_conv(i,j),dp)/year_sec
632  vz_s(j,i) = real(vz_s_conv(i,j),dp)/year_sec
633  temp_b(j,i) = real(temp_b_conv(i,j),dp)
634  temph_b(j,i) = real(temph_b_conv(i,j),dp)
635  p_b_w(j,i) = real(p_b_w_conv(i,j),dp)
636  h_w(j,i) = real(H_w_conv(i,j),dp)
637 
638  do kr=0, krmax
639  temp_r(kr,j,i) = real(temp_r_conv(i,j,kr),dp)
640  end do
641 
642  do kt=0, ktmax
643  vx_t(kt,j,i) = real(vx_t_conv(i,j,kt),dp)/year_sec
644  vy_t(kt,j,i) = real(vy_t_conv(i,j,kt),dp)/year_sec
645  vz_t(kt,j,i) = real(vz_t_conv(i,j,kt),dp)/year_sec
646  omega_t(kt,j,i) = real(omega_t_conv(i,j,kt),dp)
647  age_t(kt,j,i) = real(age_t_conv(i,j,kt),dp)*year_sec
648  enth_t(kt,j,i) = real(enth_t_conv(i,j,kt),dp)
649  enh_t(kt,j,i) = real(enh_t_conv(i,j,kt),dp)
650  end do
651 
652  do kc=0, kcmax
653  vx_c(kc,j,i) = real(vx_c_conv(i,j,kc),dp)/year_sec
654  vy_c(kc,j,i) = real(vy_c_conv(i,j,kc),dp)/year_sec
655  vz_c(kc,j,i) = real(vz_c_conv(i,j,kc),dp)/year_sec
656  temp_c(kc,j,i) = real(temp_c_conv(i,j,kc),dp)
657  age_c(kc,j,i) = real(age_c_conv(i,j,kc),dp)*year_sec
658  enth_c(kc,j,i) = real(enth_c_conv(i,j,kc),dp)
659  omega_c(kc,j,i) = real(omega_c_conv(i,j,kc),dp)
660  enh_c(kc,j,i) = real(enh_c_conv(i,j,kc),dp)
661  end do
662 
663  end do
664  end do
665 
666  end subroutine read_erg_nc
667 
668 !-------------------------------------------------------------------------------
669 !> Reading of the target-topography file (in NetCDF format).
670 !<------------------------------------------------------------------------------
671  subroutine read_target_topo_nc(target_topo_dat_name)
673  use sico_variables_m, only : maske_target, &
675  h_target
676 
677 #if (NETCDF==2)
678  use netcdf
679  use nc_check_m
680 #endif
681 
682  implicit none
683 
684  character(len=100), intent(in) :: target_topo_dat_name
685 
686 ! Return variables
687 ! (defined as global variables in module sico_variables_m):
688 !
689 ! maske_target, zs_target, zb_target, zl_target, H_target
690 
691  integer(i2b), dimension(0:IMAX,0:JMAX) :: maske_conv
692  real(sp), dimension(0:IMAX,0:JMAX) :: zs_conv, zb_conv, zl_conv, H_conv
693  character(len=256) :: target_topo_dat_path
694  character(len=256) :: filename_with_path
695 
696 #if (NETCDF==2)
697  integer(i4b) :: ncid, ncv
698  ! ncid: ID of the output file
699  ! ncv: Variable ID
700 #endif
701 
702  character(len=64), parameter :: thisroutine = 'read_target_topo_nc'
703 
704 !-------- Read target-topography data
705 ! (from a time-slice file of a previous simulation) --------
706 
707 #if (NETCDF==2)
708 
709 #if (defined(TARGET_TOPO_DAT_PATH))
710  target_topo_dat_path = trim(target_topo_dat_path)
711 #else
712  target_topo_dat_path = 'dummy'
713  stop ' >>> read_target_topo_nc: TARGET_TOPO_DAT_PATH must be defined!'
714 #endif
715 
716  filename_with_path &
717  = trim(target_topo_dat_path)//'/'//trim(target_topo_dat_name)
718 
719  call check( nf90_open(trim(filename_with_path), nf90_nowrite, ncid), &
720  thisroutine )
721 
722  call check( nf90_inq_varid(ncid, 'maske', ncv), thisroutine )
723  call check( nf90_get_var(ncid, ncv, maske_conv), thisroutine )
724 
725  call check( nf90_inq_varid(ncid, 'zs', ncv), thisroutine )
726  call check( nf90_get_var(ncid, ncv, zs_conv), thisroutine )
727 
728  call check( nf90_inq_varid(ncid, 'zb', ncv), thisroutine )
729  call check( nf90_get_var(ncid, ncv, zb_conv), thisroutine )
730 
731  call check( nf90_inq_varid(ncid, 'zl', ncv), thisroutine )
732  call check( nf90_get_var(ncid, ncv, zl_conv), thisroutine )
733 
734  call check( nf90_inq_varid(ncid, 'H', ncv), thisroutine )
735  call check( nf90_get_var(ncid, ncv, h_conv), thisroutine )
736 
737  call check( nf90_close(ncid), thisroutine )
738 
739 #else /* NETCDF != 2, not allowed */
740 
741  maske_conv = 0_i2b ! dummy values
742  zs_conv = 0.0_sp ! dummy values
743  zb_conv = 0.0_sp ! dummy values
744  zl_conv = 0.0_sp ! dummy values
745  h_conv = 0.0_sp ! dummy values
746 
747  stop ' >>> read_target_topo_nc: Parameter NETCDF must be set to 2!'
748 
749 #endif
750 
751 !-------- Convert data to double precision --------
752 
753  maske_target = transpose(maske_conv)
754  zs_target = real(transpose(zs_conv),dp)
755  zb_target = real(transpose(zb_conv),dp)
756  zl_target = real(transpose(zl_conv),dp)
757  h_target = real(transpose(H_conv) ,dp)
758 
759  end subroutine read_target_topo_nc
760 
761 !-------------------------------------------------------------------------------
762 !> Reading of the tabulated kei function.
763 !<------------------------------------------------------------------------------
764  subroutine read_kei()
767 
768  implicit none
769 
770  integer(i4b) :: n, n_data_kei_max
771  integer(i4b) :: ios
772  real(dp) :: r_val, d_dummy
773  character :: ch_dummy
774  character(len=256) :: filename_with_path
775 
776  n_data_kei_max = 10000
777 
778  kei = 0.0_dp
779 
780 !-------- Reading of tabulated values --------
781 
782  filename_with_path = trim(inpath)//'/general/kei.dat'
783 
784  open(unit=11, iostat=ios, file=trim(filename_with_path), status='old')
785 
786  if (ios /= 0) stop ' >>> read_kei: Error when opening the kei file!'
787 
788  read(unit=11, fmt='(a)') ch_dummy
789  read(unit=11, fmt='(15x,f5.0)') kei_r_max
790  read(unit=11, fmt='(15x,f5.0)') kei_r_incr
791  read(unit=11, fmt='(a)') ch_dummy
792 
794 
795  if (n_data_kei > n_data_kei_max) stop ' >>> read_kei: Array kei too small!'
796 
797  read(unit=11, fmt='(a)') ch_dummy
798  read(unit=11, fmt='(a)') ch_dummy
799  read(unit=11, fmt='(a)') ch_dummy
800 
801  do n=-n_data_kei, n_data_kei
802  read(unit=11, fmt=*) d_dummy, kei(n)
803  end do
804 
805  close(unit=11, status='keep')
806 
807  end subroutine read_kei
808 
809 !-------------------------------------------------------------------------------
810 !> Reading of physical parameters.
811 !<------------------------------------------------------------------------------
812  subroutine read_phys_para()
814  use sico_variables_m
815  use sico_vars_m
816 
817  implicit none
818 
819  integer(i4b), parameter :: n_unit=31
820  integer(i4b) :: ios
821  integer(i4b) :: n
822  character(len=256) :: filename_with_path
823 
824  filename_with_path = trim(inpath)//'/'//trim(ch_domain_short)//'/'// &
825  trim(phys_para_file)
826 
827  open(n_unit, iostat=ios, file=trim(filename_with_path), status='old')
828 
829  if (ios /= 0) &
830  stop ' >>> read_phys_para: Error when opening the phys_para file!'
831 
832 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */
833  call read_phys_para_value(rho, n_unit)
834 #else /* Martian polar caps */
835  call read_phys_para_value(rho_i, n_unit)
836 #endif
837  call read_phys_para_value(rho_w, n_unit)
838  call read_phys_para_value(rho_sw, n_unit)
839  call read_phys_para_value(l, n_unit)
840  call read_phys_para_value(g, n_unit)
841  call read_phys_para_value(nue, n_unit)
842  call read_phys_para_value(beta, n_unit)
843 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */
844  call read_phys_para_value(delta_tm_sw, n_unit)
845 #endif
846  call read_phys_para_value(omega_max, n_unit)
847 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
848  call read_phys_para_value(rho_c, n_unit)
849  call read_phys_para_value(kappa_c, n_unit)
850  call read_phys_para_value(c_c, n_unit)
851 #endif
852  call read_phys_para_value(h_r, n_unit)
853  call read_phys_para_value(rho_c_r, n_unit)
854  call read_phys_para_value(kappa_r, n_unit)
855  call read_phys_para_value(rho_a, n_unit)
856  call read_phys_para_value(r_t, n_unit)
857 
858  call read_phys_para_value(r, n_unit)
859  call read_phys_para_value(a, n_unit)
860  call read_phys_para_value(b, n_unit)
861  call read_phys_para_value(phi0, n_unit)
862  call read_phys_para_value(lambda0, n_unit)
863 
864 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */
865  call read_phys_para_value(s_stat_0, n_unit)
866 #if (!defined(GRL)) /* not Greenland */
867  call read_phys_para_value(beta1_0, n_unit)
868  call read_phys_para_value(beta2_0, n_unit)
869 #else /* Greenland */
870  call read_phys_para_value(beta1_lt_0, n_unit)
871  call read_phys_para_value(beta1_ht_0, n_unit)
872  call read_phys_para_value(beta2_lt_0, n_unit)
873  call read_phys_para_value(beta2_ht_0, n_unit)
874  call read_phys_para_value(phi_sep_0, n_unit)
875 #endif
876  call read_phys_para_value(pmax_0, n_unit)
877  call read_phys_para_value(mu_0, n_unit)
878 #endif
879 
880  do n=10, -190, -1
881  call read_phys_para_value(rf(n), n_unit)
882  end do
883 
884  do n=10, -190, -1
885  call read_phys_para_value(kappa(n), n_unit)
886  end do
887 
888  do n=10, -190, -1
889  call read_phys_para_value(c(n), n_unit)
890  end do
891 
892  close(n_unit, status='keep')
893 
894  end subroutine read_phys_para
895 
896 !-------------------------------------------------------------------------------
897 !> Reading of a value of a physical parameter from the phys_para file.
898 !<------------------------------------------------------------------------------
899  subroutine read_phys_para_value(para, n_unit)
901  implicit none
902 
903  integer(i4b), intent(in) :: n_unit
904  real(dp), intent(out) :: para
905 
906  character :: check, ch_dummy
907  character(len=*), parameter :: fmt1='(a)', fmt3='(15x)'
908  logical :: first_read
909 
910  first_read = .true.
911 
912  do
913  if (first_read) then
914  read(unit=n_unit, fmt=trim(fmt1), advance='no') check
915  first_read=.false.
916  else
917  read(unit=n_unit, fmt=trim(fmt1)) ch_dummy
918  read(unit=n_unit, fmt=trim(fmt1), advance='no') check
919  end if
920 
921  if (check /= '%') then ! no comment line
922  read(unit=n_unit, fmt=trim(fmt3), advance='no')
923  read(unit=n_unit, fmt=*) para
924  exit
925  end if
926  end do
927 
928  end subroutine read_phys_para_value
929 
930 !-------------------------------------------------------------------------------
931 
932 end module read_m
933 !
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:765
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:813
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:900
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:672
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))