45 #if (defined(ALLOW_GRDCHK) || defined(ALLOW_OPENAD)) 46 public :: read_ad_data
55 opt_maske, opt_n_cts, opt_kc_cts, &
56 opt_h_cold, opt_h_temp, opt_h, &
57 opt_temp_r, opt_omega_t, opt_age_t, &
58 opt_temp_c, opt_age_c, opt_omega_c, &
59 opt_flag_temp_age_only)
64 #if (NETCDF==2) /* time-slice file in NetCDF format */ 71 character(len=100),
intent(in) :: filename
72 real(dp),
intent(inout) :: z_sl
74 logical,
optional,
intent(in) :: opt_flag_temp_age_only
75 integer(i1b),
optional,
dimension(0:JMAX,0:IMAX), &
76 intent(inout) :: opt_maske, opt_n_cts
77 integer(i4b),
optional,
dimension(0:JMAX,0:IMAX), &
78 intent(inout) :: opt_kc_cts
79 real(dp),
optional,
dimension(0:JMAX,0:IMAX), &
80 intent(inout) :: opt_H_cold, opt_H_temp, opt_H
81 real(dp),
optional,
dimension(0:KRMAX,0:JMAX,0:IMAX), &
82 intent(inout) :: opt_temp_r
83 real(dp),
optional,
dimension(0:KTMAX,0:JMAX,0:IMAX), &
84 intent(inout) :: opt_omega_t, opt_age_t
85 real(dp),
optional,
dimension(0:KCMAX,0:JMAX,0:IMAX), &
86 intent(inout) :: opt_temp_c, opt_age_c, opt_omega_c
88 integer(i4b) :: i, j, kc, kt, kr, n
89 character(len=256) :: filename_with_path
90 logical :: flag_temp_age_only
92 #if (NETCDF==1) /* time-slice file in native binary format */ 95 character(len=256) :: ch_attr_title, ch_attr_institution, ch_attr_source, &
96 ch_attr_history, ch_attr_references
98 #elif (NETCDF==2) /* time-slice file in NetCDF format */ 100 integer(i4b) :: ncid, ncv
105 errormsg =
' >>> read_erg_nc: Parameter NETCDF must be either 1 or 2!' 109 integer(i1b),
dimension(0:IMAX,0:JMAX) :: maske_conv, maske_old_conv, &
111 flag_shelfy_stream_x_conv, &
112 flag_shelfy_stream_y_conv, &
113 flag_shelfy_stream_conv, &
114 flag_grounding_line_1_conv, &
115 flag_grounding_line_2_conv, &
116 flag_calving_front_1_conv, &
117 flag_calving_front_2_conv, &
118 flag_grounded_front_a_1_conv, &
119 flag_grounded_front_a_2_conv, &
120 flag_grounded_front_b_1_conv, &
121 flag_grounded_front_b_2_conv
122 integer(i4b),
dimension(0:IMAX,0:JMAX) :: kc_cts_conv
123 real(sp) :: time_conv, dummy_conv, z_sl_conv, &
124 V_tot_conv, A_grounded_conv, A_floating_conv, &
126 xi_conv(0:imax), eta_conv(0:jmax), &
127 sigma_level_c_conv(0:kcmax), sigma_level_t_conv(0:ktmax), &
128 sigma_level_r_conv(0:krmax)
129 real(sp),
dimension(0:IMAX,0:JMAX) :: lambda_conv, phi_conv, &
130 lon_conv, lat_conv, &
131 temp_s_conv, accum_conv, &
132 snowfall_conv, rainfall_conv, pdd_conv, &
133 as_perp_conv, as_perp_apl_conv, smb_corr_conv, &
135 zs_conv, zm_conv, zb_conv, zl_conv, zl0_conv, &
136 H_cold_conv, H_temp_conv, H_conv, &
137 Q_bm_conv, Q_tld_conv, &
140 vx_m_sia_conv, vy_m_sia_conv, vx_m_ssa_conv, vy_m_ssa_conv, &
141 dzs_dtau_conv, dzm_dtau_conv, dzb_dtau_conv, dzl_dtau_conv, &
142 dH_c_dtau_conv, dH_t_dtau_conv, dH_dtau_conv, &
143 vx_b_g_conv, vy_b_g_conv, vz_b_conv, vh_b_conv, &
144 vx_s_g_conv, vy_s_g_conv, vz_s_conv, vh_s_conv, &
145 vx_m_g_conv, vy_m_g_conv, vh_m_conv, &
146 temp_b_conv, temph_b_conv, &
147 tau_b_driving_conv, tau_b_drag_conv, &
148 p_b_w_conv, q_w_conv, q_w_x_conv, q_w_y_conv, H_w_conv, &
149 q_gl_g_conv, q_cf_g_conv, &
150 ratio_sl_x_conv, ratio_sl_y_conv, &
151 vis_ave_g_conv, vis_int_g_conv
152 real(sp),
dimension(0:IMAX,0:JMAX,0:KCMAX) :: vx_c_conv, vy_c_conv, vz_c_conv, &
153 temp_c_conv, age_c_conv, &
154 enth_c_conv, omega_c_conv, &
156 real(sp),
dimension(0:IMAX,0:JMAX,0:KTMAX) :: vx_t_conv, vy_t_conv, vz_t_conv, &
157 omega_t_conv, age_t_conv, &
160 real(sp),
dimension(0:IMAX,0:JMAX,0:KRMAX) :: temp_r_conv
162 #if (DISC>0) /* Ice discharge parameterisation */ 163 integer(i1b),
dimension(0:IMAX,0:JMAX) :: mask_mar_conv
164 real(sp),
dimension(0:IMAX,0:JMAX) :: dis_perp_conv, &
165 cst_dist_conv, cos_grad_tc_conv
168 real(dp) :: visc_min, visc_max, visc_init
169 logical :: flag_vis_ave_g
173 if (
present(opt_flag_temp_age_only))
then 174 flag_temp_age_only = opt_flag_temp_age_only
176 flag_temp_age_only = .false.
179 filename_with_path = trim(anfdatpath)//
'/'//trim(filename)
181 #if (NETCDF==1) /* time-slice file in native binary format */ 183 #if !defined(ALLOW_OPENAD) /* Normal */ 184 open(unit=11, iostat=ios, file=trim(filename_with_path), status=
'old', &
187 open(unit=11, iostat=ios, file=trim(filename_with_path))
188 #endif /* Normal vs. OpenAD */ 191 errormsg =
' >>> read_erg_nc: Error when opening the time-slice file!' 195 read(unit=11) ch_attr_title
196 read(unit=11) ch_attr_institution
197 read(unit=11) ch_attr_source
198 read(unit=11) ch_attr_history
199 read(unit=11) ch_attr_references
201 read(unit=11) time_conv
202 read(unit=11) dummy_conv
203 read(unit=11) z_sl_conv
205 read(unit=11) dummy_conv
206 read(unit=11) dummy_conv
207 read(unit=11) dummy_conv
208 read(unit=11) dummy_conv
210 read(unit=11) xi_conv
211 read(unit=11) eta_conv
212 read(unit=11) sigma_level_c_conv
213 read(unit=11) sigma_level_t_conv
214 read(unit=11) sigma_level_r_conv
215 read(unit=11) lon_conv
216 read(unit=11) lat_conv
217 read(unit=11) lambda_conv
218 read(unit=11) phi_conv
219 read(unit=11) temp_s_conv
220 read(unit=11) accum_conv
221 read(unit=11) snowfall_conv
222 read(unit=11) rainfall_conv
223 read(unit=11) pdd_conv
224 read(unit=11) as_perp_conv
225 read(unit=11) as_perp_apl_conv
226 read(unit=11) smb_corr_conv
228 #if (DISC>0) /* Ice discharge parameterisation */ 230 read(unit=11) dis_perp_conv
231 read(unit=11) cst_dist_conv
232 read(unit=11) cos_grad_tc_conv
233 read(unit=11) mask_mar_conv
237 read(unit=11) q_geo_conv
238 read(unit=11) maske_conv
239 read(unit=11) maske_old_conv
240 read(unit=11) n_cts_conv
241 read(unit=11) kc_cts_conv
242 read(unit=11) zs_conv
243 read(unit=11) zm_conv
244 read(unit=11) zb_conv
245 read(unit=11) zl_conv
246 read(unit=11) zl0_conv
247 read(unit=11) h_cold_conv
248 read(unit=11) h_temp_conv
250 read(unit=11) h_r_conv
251 read(unit=11) q_bm_conv
252 read(unit=11) q_tld_conv
253 read(unit=11) am_perp_conv
254 read(unit=11) qx_conv
255 read(unit=11) qy_conv
256 read(unit=11) vx_m_sia_conv
257 read(unit=11) vy_m_sia_conv
258 read(unit=11) vx_m_ssa_conv
259 read(unit=11) vy_m_ssa_conv
260 read(unit=11) dzs_dtau_conv
261 read(unit=11) dzm_dtau_conv
262 read(unit=11) dzb_dtau_conv
263 read(unit=11) dzl_dtau_conv
264 read(unit=11) dh_c_dtau_conv
265 read(unit=11) dh_t_dtau_conv
266 read(unit=11) dh_dtau_conv
267 read(unit=11) vx_b_g_conv
268 read(unit=11) vy_b_g_conv
269 read(unit=11) vz_b_conv
270 read(unit=11) vh_b_conv
271 read(unit=11) vx_s_g_conv
272 read(unit=11) vy_s_g_conv
273 read(unit=11) vz_s_conv
274 read(unit=11) vh_s_conv
275 read(unit=11) vx_m_g_conv
276 read(unit=11) vy_m_g_conv
277 read(unit=11) vh_m_conv
278 read(unit=11) temp_b_conv
279 read(unit=11) temph_b_conv
280 read(unit=11) tau_b_driving_conv
281 read(unit=11) tau_b_drag_conv
282 read(unit=11) p_b_w_conv
283 read(unit=11) q_w_conv
284 read(unit=11) q_w_x_conv
285 read(unit=11) q_w_y_conv
286 read(unit=11) h_w_conv
287 read(unit=11) q_gl_g_conv
288 read(unit=11) q_cf_g_conv
289 read(unit=11) ratio_sl_x_conv
290 read(unit=11) ratio_sl_y_conv
291 read(unit=11) flag_shelfy_stream_x_conv
292 read(unit=11) flag_shelfy_stream_y_conv
293 read(unit=11) flag_shelfy_stream_conv
294 read(unit=11) flag_grounding_line_1_conv
295 read(unit=11) flag_grounding_line_2_conv
296 read(unit=11) flag_calving_front_1_conv
297 read(unit=11) flag_calving_front_2_conv
298 read(unit=11) flag_grounded_front_a_1_conv
299 read(unit=11) flag_grounded_front_a_2_conv
300 read(unit=11) flag_grounded_front_b_1_conv
301 read(unit=11) flag_grounded_front_b_2_conv
302 read(unit=11) vis_ave_g_conv; flag_vis_ave_g = .true.
304 read(unit=11) vis_int_g_conv
306 read(unit=11) vx_c_conv
307 read(unit=11) vy_c_conv
308 read(unit=11) vz_c_conv
309 read(unit=11) vx_t_conv
310 read(unit=11) vy_t_conv
311 read(unit=11) vz_t_conv
312 read(unit=11) temp_c_conv
313 read(unit=11) omega_t_conv
314 read(unit=11) temp_r_conv
315 read(unit=11) enth_c_conv
316 read(unit=11) enth_t_conv
317 read(unit=11) omega_c_conv
318 read(unit=11) enh_c_conv
319 read(unit=11) enh_t_conv
320 read(unit=11) age_c_conv
321 read(unit=11) age_t_conv
323 close(unit=11,status=
'keep')
325 #elif (NETCDF==2) /* time-slice file in NetCDF format */ 327 call check( nf90_open(trim(filename_with_path), nf90_nowrite, ncid) )
329 call check( nf90_inq_varid(ncid,
'time', ncv) )
330 call check( nf90_get_var(ncid, ncv, time_conv) )
332 if ( nf90_inq_varid(ncid,
'delta_ts', ncv) == nf90_noerr )
then 333 call check( nf90_get_var(ncid, ncv, dummy_conv) )
334 else if ( nf90_inq_varid(ncid,
'glac_index', ncv) == nf90_noerr )
then 335 call check( nf90_get_var(ncid, ncv, dummy_conv) )
340 call check( nf90_inq_varid(ncid,
'z_sl', ncv) )
341 call check( nf90_get_var(ncid, ncv, z_sl_conv) )
343 call check( nf90_inq_varid(ncid,
'x', ncv) )
344 call check( nf90_get_var(ncid, ncv, xi_conv) )
346 call check( nf90_inq_varid(ncid,
'y', ncv) )
347 call check( nf90_get_var(ncid, ncv, eta_conv) )
349 call check( nf90_inq_varid(ncid,
'sigma_level_c', ncv) )
350 call check( nf90_get_var(ncid, ncv, sigma_level_c_conv) )
352 call check( nf90_inq_varid(ncid,
'sigma_level_t', ncv) )
353 call check( nf90_get_var(ncid, ncv, sigma_level_t_conv) )
355 call check( nf90_inq_varid(ncid,
'sigma_level_r', ncv) )
356 call check( nf90_get_var(ncid, ncv, sigma_level_r_conv) )
358 call check( nf90_inq_varid(ncid,
'lon', ncv) )
359 call check( nf90_get_var(ncid, ncv, lon_conv) )
361 call check( nf90_inq_varid(ncid,
'lat', ncv) )
362 call check( nf90_get_var(ncid, ncv, lat_conv) )
364 call check( nf90_inq_varid(ncid,
'lambda', ncv) )
365 call check( nf90_get_var(ncid, ncv, lambda_conv) )
367 call check( nf90_inq_varid(ncid,
'phi', ncv) )
368 call check( nf90_get_var(ncid, ncv, phi_conv) )
370 call check( nf90_inq_varid(ncid,
'temp_s', ncv) )
371 call check( nf90_get_var(ncid, ncv, temp_s_conv) )
373 call check( nf90_inq_varid(ncid,
'prec', ncv) )
374 call check( nf90_get_var(ncid, ncv, accum_conv) )
376 call check( nf90_inq_varid(ncid,
'snowfall', ncv) )
377 call check( nf90_get_var(ncid, ncv, snowfall_conv) )
379 call check( nf90_inq_varid(ncid,
'rainfall', ncv) )
380 call check( nf90_get_var(ncid, ncv, rainfall_conv) )
382 call check( nf90_inq_varid(ncid,
'pdd', ncv) )
383 call check( nf90_get_var(ncid, ncv, pdd_conv) )
385 call check( nf90_inq_varid(ncid,
'as_perp', ncv) )
386 call check( nf90_get_var(ncid, ncv, as_perp_conv) )
388 call check( nf90_inq_varid(ncid,
'as_perp_apl', ncv) )
389 call check( nf90_get_var(ncid, ncv, as_perp_apl_conv) )
391 call check( nf90_inq_varid(ncid,
'smb_corr', ncv) )
392 call check( nf90_get_var(ncid, ncv, smb_corr_conv) )
394 #if (DISC>0) /* Ice discharge parameterisation */ 396 call check( nf90_inq_varid(ncid,
'dis_perp', ncv) )
397 call check( nf90_get_var(ncid, ncv, dis_perp_conv) )
399 call check( nf90_inq_varid(ncid,
'cst_dist', ncv) )
400 call check( nf90_get_var(ncid, ncv, cst_dist_conv) )
402 call check( nf90_inq_varid(ncid,
'cos_grad_tc', ncv) )
403 call check( nf90_get_var(ncid, ncv, cos_grad_tc_conv) )
405 call check( nf90_inq_varid(ncid,
'mask_mar', ncv) )
406 call check( nf90_get_var(ncid, ncv, mask_mar_conv) )
410 call check( nf90_inq_varid(ncid,
'q_geo', ncv) )
411 call check( nf90_get_var(ncid, ncv, q_geo_conv) )
413 call check( nf90_inq_varid(ncid,
'maske', ncv) )
414 call check( nf90_get_var(ncid, ncv, maske_conv) )
416 call check( nf90_inq_varid(ncid,
'maske_old', ncv) )
417 call check( nf90_get_var(ncid, ncv, maske_old_conv) )
419 call check( nf90_inq_varid(ncid,
'n_cts', ncv) )
420 call check( nf90_get_var(ncid, ncv, n_cts_conv) )
422 call check( nf90_inq_varid(ncid,
'kc_cts', ncv) )
423 call check( nf90_get_var(ncid, ncv, kc_cts_conv) )
425 call check( nf90_inq_varid(ncid,
'zs', ncv) )
426 call check( nf90_get_var(ncid, ncv, zs_conv) )
428 call check( nf90_inq_varid(ncid,
'zm', ncv) )
429 call check( nf90_get_var(ncid, ncv, zm_conv) )
431 call check( nf90_inq_varid(ncid,
'zb', ncv) )
432 call check( nf90_get_var(ncid, ncv, zb_conv) )
434 call check( nf90_inq_varid(ncid,
'zl', ncv) )
435 call check( nf90_get_var(ncid, ncv, zl_conv) )
437 call check( nf90_inq_varid(ncid,
'zl0', ncv) )
438 call check( nf90_get_var(ncid, ncv, zl0_conv) )
440 call check( nf90_inq_varid(ncid,
'H_cold', ncv) )
441 call check( nf90_get_var(ncid, ncv, h_cold_conv) )
443 call check( nf90_inq_varid(ncid,
'H_temp', ncv) )
444 call check( nf90_get_var(ncid, ncv, h_temp_conv) )
446 call check( nf90_inq_varid(ncid,
'H', ncv) )
447 call check( nf90_get_var(ncid, ncv, h_conv) )
449 call check( nf90_inq_varid(ncid,
'H_R', ncv) )
450 call check( nf90_get_var(ncid, ncv, h_r_conv) )
452 call check( nf90_inq_varid(ncid,
'Q_bm', ncv) )
453 call check( nf90_get_var(ncid, ncv, q_bm_conv) )
455 call check( nf90_inq_varid(ncid,
'Q_tld', ncv) )
456 call check( nf90_get_var(ncid, ncv, q_tld_conv) )
458 call check( nf90_inq_varid(ncid,
'am_perp', ncv) )
459 call check( nf90_get_var(ncid, ncv, am_perp_conv) )
461 call check( nf90_inq_varid(ncid,
'qx', ncv) )
462 call check( nf90_get_var(ncid, ncv, qx_conv) )
464 call check( nf90_inq_varid(ncid,
'qy', ncv) )
465 call check( nf90_get_var(ncid, ncv, qy_conv) )
467 if ( nf90_inq_varid(ncid,
'vx_m_sia', ncv) == nf90_noerr )
then 468 call check( nf90_get_var(ncid, ncv, vx_m_sia_conv) )
470 write(6,
'(/1x,a)')
'>>> read_erg_nc: Variable vx_m_sia' 471 write(6,
'(1x,a)')
' not available in read file *.nc.' 472 vx_m_sia_conv = 0.0_sp
475 if ( nf90_inq_varid(ncid,
'vy_m_sia', ncv) == nf90_noerr )
then 476 call check( nf90_get_var(ncid, ncv, vy_m_sia_conv) )
478 write(6,
'(/1x,a)')
'>>> read_erg_nc: Variable vy_m_sia' 479 write(6,
'(1x,a)')
' not available in read file *.nc.' 480 vy_m_sia_conv = 0.0_sp
483 if ( nf90_inq_varid(ncid,
'vx_m_ssa', ncv) == nf90_noerr )
then 484 call check( nf90_get_var(ncid, ncv, vx_m_ssa_conv) )
486 write(6,
'(/1x,a)')
'>>> read_erg_nc: Variable vx_m_ssa' 487 write(6,
'(1x,a)')
' not available in read file *.nc.' 488 vx_m_ssa_conv = 0.0_sp
491 if ( nf90_inq_varid(ncid,
'vy_m_ssa', ncv) == nf90_noerr )
then 492 call check( nf90_get_var(ncid, ncv, vy_m_ssa_conv) )
494 write(6,
'(/1x,a)')
'>>> read_erg_nc: Variable vy_m_ssa' 495 write(6,
'(1x,a)')
' not available in read file *.nc.' 496 vy_m_ssa_conv = 0.0_sp
499 call check( nf90_inq_varid(ncid,
'dzs_dt', ncv) )
500 call check( nf90_get_var(ncid, ncv, dzs_dtau_conv) )
502 call check( nf90_inq_varid(ncid,
'dzm_dt', ncv) )
503 call check( nf90_get_var(ncid, ncv, dzm_dtau_conv) )
505 call check( nf90_inq_varid(ncid,
'dzb_dt', ncv) )
506 call check( nf90_get_var(ncid, ncv, dzb_dtau_conv) )
508 call check( nf90_inq_varid(ncid,
'dzl_dt', ncv) )
509 call check( nf90_get_var(ncid, ncv, dzl_dtau_conv) )
511 call check( nf90_inq_varid(ncid,
'dH_c_dt', ncv) )
512 call check( nf90_get_var(ncid, ncv, dh_c_dtau_conv) )
514 call check( nf90_inq_varid(ncid,
'dH_t_dt', ncv) )
515 call check( nf90_get_var(ncid, ncv, dh_t_dtau_conv) )
517 call check( nf90_inq_varid(ncid,
'dH_dt', ncv) )
518 call check( nf90_get_var(ncid, ncv, dh_dtau_conv) )
520 call check( nf90_inq_varid(ncid,
'vx_b_g', ncv) )
521 call check( nf90_get_var(ncid, ncv, vx_b_g_conv) )
523 call check( nf90_inq_varid(ncid,
'vy_b_g', ncv) )
524 call check( nf90_get_var(ncid, ncv, vy_b_g_conv) )
526 call check( nf90_inq_varid(ncid,
'vz_b', ncv) )
527 call check( nf90_get_var(ncid, ncv, vz_b_conv) )
529 call check( nf90_inq_varid(ncid,
'vh_b', ncv) )
530 call check( nf90_get_var(ncid, ncv, vh_b_conv) )
532 call check( nf90_inq_varid(ncid,
'vx_s_g', ncv) )
533 call check( nf90_get_var(ncid, ncv, vx_s_g_conv) )
535 call check( nf90_inq_varid(ncid,
'vy_s_g', ncv) )
536 call check( nf90_get_var(ncid, ncv, vy_s_g_conv) )
538 call check( nf90_inq_varid(ncid,
'vz_s', ncv) )
539 call check( nf90_get_var(ncid, ncv, vz_s_conv) )
541 call check( nf90_inq_varid(ncid,
'vh_s', ncv) )
542 call check( nf90_get_var(ncid, ncv, vh_s_conv) )
544 call check( nf90_inq_varid(ncid,
'vx_m_g', ncv) )
545 call check( nf90_get_var(ncid, ncv, vx_m_g_conv) )
547 call check( nf90_inq_varid(ncid,
'vy_m_g', ncv) )
548 call check( nf90_get_var(ncid, ncv, vy_m_g_conv) )
550 call check( nf90_inq_varid(ncid,
'vh_m', ncv) )
551 call check( nf90_get_var(ncid, ncv, vh_m_conv) )
553 call check( nf90_inq_varid(ncid,
'temp_b', ncv) )
554 call check( nf90_get_var(ncid, ncv, temp_b_conv) )
556 call check( nf90_inq_varid(ncid,
'temph_b', ncv) )
557 call check( nf90_get_var(ncid, ncv, temph_b_conv) )
559 call check( nf90_inq_varid(ncid,
'tau_b_driving', ncv) )
560 call check( nf90_get_var(ncid, ncv, tau_b_driving_conv) )
562 call check( nf90_inq_varid(ncid,
'tau_b_drag', ncv) )
563 call check( nf90_get_var(ncid, ncv, tau_b_drag_conv) )
565 call check( nf90_inq_varid(ncid,
'p_b_w', ncv) )
566 call check( nf90_get_var(ncid, ncv, p_b_w_conv) )
568 call check( nf90_inq_varid(ncid,
'q_w', ncv) )
569 call check( nf90_get_var(ncid, ncv, q_w_conv) )
571 call check( nf90_inq_varid(ncid,
'q_w_x', ncv) )
572 call check( nf90_get_var(ncid, ncv, q_w_x_conv) )
574 call check( nf90_inq_varid(ncid,
'q_w_y', ncv) )
575 call check( nf90_get_var(ncid, ncv, q_w_y_conv) )
577 call check( nf90_inq_varid(ncid,
'H_w', ncv) )
578 call check( nf90_get_var(ncid, ncv, h_w_conv) )
580 call check( nf90_inq_varid(ncid,
'q_gl_g', ncv) )
581 call check( nf90_get_var(ncid, ncv, q_gl_g_conv) )
583 call check( nf90_inq_varid(ncid,
'q_cf_g', ncv) )
584 call check( nf90_get_var(ncid, ncv, q_cf_g_conv) )
586 call check( nf90_inq_varid(ncid,
'ratio_sl_x', ncv) )
587 call check( nf90_get_var(ncid, ncv, ratio_sl_x_conv) )
589 call check( nf90_inq_varid(ncid,
'ratio_sl_y', ncv) )
590 call check( nf90_get_var(ncid, ncv, ratio_sl_y_conv) )
592 call check( nf90_inq_varid(ncid,
'flag_shelfy_stream_x', ncv) )
593 call check( nf90_get_var(ncid, ncv, flag_shelfy_stream_x_conv) )
595 call check( nf90_inq_varid(ncid,
'flag_shelfy_stream_y', ncv) )
596 call check( nf90_get_var(ncid, ncv, flag_shelfy_stream_y_conv) )
598 call check( nf90_inq_varid(ncid,
'flag_shelfy_stream', ncv) )
599 call check( nf90_get_var(ncid, ncv, flag_shelfy_stream_conv) )
601 call check( nf90_inq_varid(ncid,
'flag_grounding_line_1', ncv) )
602 call check( nf90_get_var(ncid, ncv, flag_grounding_line_1_conv) )
604 call check( nf90_inq_varid(ncid,
'flag_grounding_line_2', ncv) )
605 call check( nf90_get_var(ncid, ncv, flag_grounding_line_2_conv) )
607 call check( nf90_inq_varid(ncid,
'flag_calving_front_1', ncv) )
608 call check( nf90_get_var(ncid, ncv, flag_calving_front_1_conv) )
610 call check( nf90_inq_varid(ncid,
'flag_calving_front_2', ncv) )
611 call check( nf90_get_var(ncid, ncv, flag_calving_front_2_conv) )
613 call check( nf90_inq_varid(ncid,
'flag_grounded_front_a_1', ncv) )
614 call check( nf90_get_var(ncid, ncv, flag_grounded_front_a_1_conv) )
616 call check( nf90_inq_varid(ncid,
'flag_grounded_front_a_2', ncv) )
617 call check( nf90_get_var(ncid, ncv, flag_grounded_front_a_2_conv) )
619 call check( nf90_inq_varid(ncid,
'flag_grounded_front_b_1', ncv) )
620 call check( nf90_get_var(ncid, ncv, flag_grounded_front_b_1_conv) )
622 call check( nf90_inq_varid(ncid,
'flag_grounded_front_b_2', ncv) )
623 call check( nf90_get_var(ncid, ncv, flag_grounded_front_b_2_conv) )
625 if ( nf90_inq_varid(ncid,
'vis_ave_g', ncv) == nf90_noerr )
then 626 call check( nf90_get_var(ncid, ncv, vis_ave_g_conv) )
627 flag_vis_ave_g = .true.
629 write(6,
'(/1x,a)')
'>>> read_erg_nc: Variable vis_ave_g' 630 write(6,
'(1x,a)')
' not available in read file *.nc.' 631 vis_ave_g_conv = 0.0_sp
632 flag_vis_ave_g = .false.
635 call check( nf90_inq_varid(ncid,
'vis_int_g', ncv) )
636 call check( nf90_get_var(ncid, ncv, vis_int_g_conv) )
638 call check( nf90_inq_varid(ncid,
'vx_c', ncv) )
639 call check( nf90_get_var(ncid, ncv, vx_c_conv) )
641 call check( nf90_inq_varid(ncid,
'vy_c', ncv) )
642 call check( nf90_get_var(ncid, ncv, vy_c_conv) )
644 call check( nf90_inq_varid(ncid,
'vz_c', ncv) )
645 call check( nf90_get_var(ncid, ncv, vz_c_conv) )
647 call check( nf90_inq_varid(ncid,
'vx_t', ncv) )
648 call check( nf90_get_var(ncid, ncv, vx_t_conv) )
650 call check( nf90_inq_varid(ncid,
'vy_t', ncv) )
651 call check( nf90_get_var(ncid, ncv, vy_t_conv) )
653 call check( nf90_inq_varid(ncid,
'vz_t', ncv) )
654 call check( nf90_get_var(ncid, ncv, vz_t_conv) )
656 call check( nf90_inq_varid(ncid,
'temp_c', ncv) )
657 call check( nf90_get_var(ncid, ncv, temp_c_conv) )
659 call check( nf90_inq_varid(ncid,
'omega_t', ncv) )
660 call check( nf90_get_var(ncid, ncv, omega_t_conv) )
662 call check( nf90_inq_varid(ncid,
'temp_r', ncv) )
663 call check( nf90_get_var(ncid, ncv, temp_r_conv) )
665 call check( nf90_inq_varid(ncid,
'enth_c', ncv) )
666 call check( nf90_get_var(ncid, ncv, enth_c_conv) )
668 call check( nf90_inq_varid(ncid,
'enth_t', ncv) )
669 call check( nf90_get_var(ncid, ncv, enth_t_conv) )
671 call check( nf90_inq_varid(ncid,
'omega_c', ncv) )
672 call check( nf90_get_var(ncid, ncv, omega_c_conv) )
674 call check( nf90_inq_varid(ncid,
'enh_c', ncv) )
675 call check( nf90_get_var(ncid, ncv, enh_c_conv) )
677 call check( nf90_inq_varid(ncid,
'enh_t', ncv) )
678 call check( nf90_get_var(ncid, ncv, enh_t_conv) )
680 call check( nf90_inq_varid(ncid,
'age_c', ncv) )
681 call check( nf90_get_var(ncid, ncv, age_c_conv) )
683 call check( nf90_inq_varid(ncid,
'age_t', ncv) )
684 call check( nf90_get_var(ncid, ncv, age_t_conv) )
686 call check( nf90_close(ncid) )
692 if (.not.flag_temp_age_only)
then 694 z_sl =
real(z_sl_conv,
dp)
696 h_r =
real(h_r_conv,
dp)
699 xi(i) =
real(xi_conv(i),dp)
703 eta(j) =
real(eta_conv(j),dp)
709 maske(j,i) = maske_conv(i,j)
711 n_cts(j,i) = n_cts_conv(i,j)
712 kc_cts(j,i) = kc_cts_conv(i,j)
713 zs(j,i) =
real(zs_conv(i,j),dp)
714 zm(j,i) =
real(zm_conv(i,j),dp)
715 zb(j,i) =
real(zb_conv(i,j),dp)
716 zl(j,i) =
real(zl_conv(i,j),dp)
718 h_c(j,i) =
real(H_cold_conv(i,j),dp)
719 h_t(j,i) =
real(H_temp_conv(i,j),dp)
720 #elif (CALCMOD==0 || CALCMOD==2 || CALCMOD==3 || CALCMOD==-1) 721 h_c(j,i) =
real(H_conv(i,j),dp)
724 q_bm(j,i) =
real(Q_bm_conv(i,j),dp)/year_sec
725 q_tld(j,i) =
real(Q_tld_conv(i,j),dp)/year_sec
726 am_perp(j,i) =
real(am_perp_conv(i,j),dp)/year_sec
727 qx(j,i) =
real(qx_conv(i,j),dp)/year_sec
728 qy(j,i) =
real(qy_conv(i,j),dp)/year_sec
729 vx_m_sia(j,i) =
real(vx_m_sia_conv(i,j),dp)/year_sec
730 vy_m_sia(j,i) =
real(vy_m_sia_conv(i,j),dp)/year_sec
731 vx_m_ssa(j,i) =
real(vx_m_ssa_conv(i,j),dp)/year_sec
732 vy_m_ssa(j,i) =
real(vy_m_ssa_conv(i,j),dp)/year_sec
733 dzs_dtau(j,i) =
real(dzs_dtau_conv(i,j),dp)/year_sec
734 dzm_dtau(j,i) =
real(dzm_dtau_conv(i,j),dp)/year_sec
735 dzb_dtau(j,i) =
real(dzb_dtau_conv(i,j),dp)/year_sec
736 dzl_dtau(j,i) =
real(dzl_dtau_conv(i,j),dp)/year_sec
737 dh_c_dtau(j,i) =
real(dH_c_dtau_conv(i,j),dp)/year_sec
738 dh_t_dtau(j,i) =
real(dH_t_dtau_conv(i,j),dp)/year_sec
739 vx_b_g(j,i) =
real(vx_b_g_conv(i,j),dp)/year_sec
740 vy_b_g(j,i) =
real(vy_b_g_conv(i,j),dp)/year_sec
741 vz_b(j,i) =
real(vz_b_conv(i,j),dp)/year_sec
742 vx_s_g(j,i) =
real(vx_s_g_conv(i,j),dp)/year_sec
743 vy_s_g(j,i) =
real(vy_s_g_conv(i,j),dp)/year_sec
744 vz_s(j,i) =
real(vz_s_conv(i,j),dp)/year_sec
745 temp_b(j,i) =
real(temp_b_conv(i,j),dp)
746 temph_b(j,i) =
real(temph_b_conv(i,j),dp)
747 p_b_w(j,i) =
real(p_b_w_conv(i,j),dp)
748 q_w(j,i) =
real(q_w_conv(i,j),dp)/year_sec
749 q_w_x(j,i) =
real(q_w_x_conv(i,j),dp)/year_sec
750 q_w_y(j,i) =
real(q_w_y_conv(i,j),dp)/year_sec
751 h_w(j,i) =
real(H_w_conv(i,j),dp)
752 ratio_sl_x(j,i) =
real(ratio_sl_x_conv(i,j),dp)
753 ratio_sl_y(j,i) =
real(ratio_sl_y_conv(i,j),dp)
755 if (flag_shelfy_stream_x_conv(i,j) == 1_i1b)
then 761 if (flag_shelfy_stream_y_conv(i,j) == 1_i1b)
then 767 if (flag_shelfy_stream_conv(i,j) == 1_i1b)
then 773 if (flag_grounding_line_1_conv(i,j) == 1_i1b)
then 779 if (flag_grounding_line_2_conv(i,j) == 1_i1b)
then 785 if (flag_calving_front_1_conv(i,j) == 1_i1b)
then 791 if (flag_calving_front_2_conv(i,j) == 1_i1b)
then 797 if (flag_grounded_front_a_1_conv(i,j) == 1_i1b)
then 803 if (flag_grounded_front_a_2_conv(i,j) == 1_i1b)
then 809 if (flag_grounded_front_b_1_conv(i,j) == 1_i1b)
then 815 if (flag_grounded_front_b_2_conv(i,j) == 1_i1b)
then 821 vis_ave_g(j,i) =
real(vis_ave_g_conv(i,j),dp)
822 vis_int_g(j,i) =
real(vis_int_g_conv(i,j),dp)
825 temp_r(kr,j,i) =
real(temp_r_conv(i,j,kr),dp)
829 vx_t(kt,j,i) =
real(vx_t_conv(i,j,kt),dp)/year_sec
830 vy_t(kt,j,i) =
real(vy_t_conv(i,j,kt),dp)/year_sec
831 vz_t(kt,j,i) =
real(vz_t_conv(i,j,kt),dp)/year_sec
832 omega_t(kt,j,i) =
real(omega_t_conv(i,j,kt),dp)
833 age_t(kt,j,i) =
real(age_t_conv(i,j,kt),dp)*year_sec
834 enth_t(kt,j,i) =
real(enth_t_conv(i,j,kt),dp)
835 enh_t(kt,j,i) =
real(enh_t_conv(i,j,kt),dp)
839 vx_c(kc,j,i) =
real(vx_c_conv(i,j,kc),dp)/year_sec
840 vy_c(kc,j,i) =
real(vy_c_conv(i,j,kc),dp)/year_sec
841 vz_c(kc,j,i) =
real(vz_c_conv(i,j,kc),dp)/year_sec
842 temp_c(kc,j,i) =
real(temp_c_conv(i,j,kc),dp)
843 age_c(kc,j,i) =
real(age_c_conv(i,j,kc),dp)*year_sec
844 enth_c(kc,j,i) =
real(enth_c_conv(i,j,kc),dp)
845 omega_c(kc,j,i) =
real(omega_c_conv(i,j,kc),dp)
846 enh_c(kc,j,i) =
real(enh_c_conv(i,j,kc),dp)
852 if (.not.flag_vis_ave_g)
then 854 #if (defined(VISC_MIN) && defined(VISC_MAX)) 858 visc_min = 1.0e+10_dp
859 visc_max = 1.0e+25_dp
862 #if (defined(VISC_INIT_SSA)) 863 visc_init = visc_init_ssa
865 visc_init = 1.0e+15_dp
871 if ((
maske(j,i)==0_i1b).or.(
maske(j,i)==3_i1b))
then 892 if (
present(opt_maske) &
893 .and.
present(opt_n_cts) &
894 .and.
present(opt_kc_cts) &
895 .and.
present(opt_h_cold) &
896 .and.
present(opt_h_temp) &
897 .and.
present(opt_h) &
898 .and.
present(opt_temp_r) &
899 .and.
present(opt_omega_t) &
900 .and.
present(opt_age_t) &
901 .and.
present(opt_temp_c) &
902 .and.
present(opt_age_c) &
903 .and.
present(opt_omega_c) &
909 opt_maske(j,i) = maske_conv(i,j)
910 opt_n_cts(j,i) = n_cts_conv(i,j)
911 opt_kc_cts(j,i) = kc_cts_conv(i,j)
913 opt_h_cold(j,i) =
real(H_cold_conv(i,j),dp)
914 opt_h_temp(j,i) =
real(H_temp_conv(i,j),dp)
915 opt_h(j,i) =
real(H_conv(i,j),dp)
918 opt_temp_r(kr,j,i) =
real(temp_r_conv(i,j,kr),dp)
922 opt_omega_t(kt,j,i) =
real(omega_t_conv(i,j,kt),dp)
923 opt_age_t(kt,j,i) =
real(age_t_conv(i,j,kt),dp)*year_sec
927 opt_temp_c(kc,j,i) =
real(temp_c_conv(i,j,kc),dp)
928 opt_age_c(kc,j,i) =
real(age_c_conv(i,j,kc),dp)*year_sec
929 opt_omega_c(kc,j,i) =
real(omega_c_conv(i,j,kc),dp)
937 errormsg =
' >>> read_erg_nc: optional argument(s) missing!' 962 character(len=100),
intent(in) :: target_topo_dat_name
969 integer(i1b),
dimension(0:IMAX,0:JMAX) :: maske_conv
970 real(sp),
dimension(0:IMAX,0:JMAX) :: zs_conv, zb_conv, zl_conv, H_conv
971 character(len=256) :: target_topo_dat_path
972 character(len=256) :: filename_with_path
975 integer(i4b) :: ncid, ncv
980 character(len=64),
parameter :: thisroutine =
'read_target_topo_nc' 982 #if defined(ALLOW_OPENAD) /* OpenAD */ 991 #if (defined(TARGET_TOPO_DAT_PATH)) 992 target_topo_dat_path = trim(target_topo_dat_path)
994 target_topo_dat_path =
'dummy' 995 errormsg =
' >>> read_target_topo_nc: TARGET_TOPO_DAT_PATH must be defined!' 1000 = trim(target_topo_dat_path)//
'/'//trim(target_topo_dat_name)
1002 call check( nf90_open(trim(filename_with_path), nf90_nowrite, ncid), &
1005 call check( nf90_inq_varid(ncid,
'maske', ncv), thisroutine )
1006 call check( nf90_get_var(ncid, ncv, maske_conv), thisroutine )
1008 call check( nf90_inq_varid(ncid,
'zs', ncv), thisroutine )
1009 call check( nf90_get_var(ncid, ncv, zs_conv), thisroutine )
1011 call check( nf90_inq_varid(ncid,
'zb', ncv), thisroutine )
1012 call check( nf90_get_var(ncid, ncv, zb_conv), thisroutine )
1014 call check( nf90_inq_varid(ncid,
'zl', ncv), thisroutine )
1015 call check( nf90_get_var(ncid, ncv, zl_conv), thisroutine )
1017 call check( nf90_inq_varid(ncid,
'H', ncv), thisroutine )
1018 call check( nf90_get_var(ncid, ncv, h_conv), thisroutine )
1020 call check( nf90_close(ncid), thisroutine )
1022 #else /* NETCDF != 2, not allowed */ 1030 errormsg =
' >>> read_target_topo_nc: Parameter NETCDF must be set to 2!' 1037 #if !defined(ALLOW_OPENAD) /* Normal */ 1042 h_target =
real(transpose(H_conv) ,dp)
1050 h_target(j,i) =
real(H_conv(i,j) ,dp)
1053 #endif /* Normal vs. OpenAD */ 1066 integer(i4b) :: n, n_data_kei_max
1068 real(dp) :: r_val, d_dummy
1069 character :: ch_dummy
1070 character(len=256) :: filename_with_path
1072 n_data_kei_max = 10000
1078 filename_with_path = trim(inpath)//
'/general/kei.dat' 1080 #if !defined(ALLOW_OPENAD) /* Normal */ 1081 open(unit=11, iostat=ios, file=trim(filename_with_path), status=
'old')
1083 open(unit=11, iostat=ios, file=trim(filename_with_path))
1084 #endif /* Normal vs. OpenAD */ 1087 errormsg =
' >>> read_kei: Error when opening the kei file!' 1091 read(unit=11, fmt=
'(a)') ch_dummy
1092 read(unit=11, fmt=
'(15x,f5.0)')
kei_r_max 1094 read(unit=11, fmt=
'(a)') ch_dummy
1099 errormsg =
' >>> read_kei: Array kei too small!' 1103 read(unit=11, fmt=
'(a)') ch_dummy
1104 read(unit=11, fmt=
'(a)') ch_dummy
1105 read(unit=11, fmt=
'(a)') ch_dummy
1108 read(unit=11, fmt=*) d_dummy,
kei(n)
1111 close(unit=11, status=
'keep')
1125 integer(i4b),
parameter :: n_unit=31
1128 character(len=256) :: filename_with_path
1130 filename_with_path = trim(inpath)//
'/'//trim(
ch_domain_short)//
'/'// &
1131 trim(phys_para_file)
1133 #if !defined(ALLOW_OPENAD) /* Normal */ 1134 open(n_unit, iostat=ios, file=trim(filename_with_path), status=
'old')
1136 open(n_unit, iostat=ios, file=trim(filename_with_path))
1137 #endif /* Normal vs. OpenAD */ 1140 errormsg =
' >>> read_phys_para: Error when opening the phys_para file!' 1144 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */ 1146 #else /* Martian polar caps */ 1155 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */ 1159 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */ 1176 #if (!defined(NMARS) && !defined(SMARS)) /* not Martian polar caps */ 1178 #if (!defined(GRL)) /* not Greenland */ 1181 #else /* Greenland */ 1204 close(n_unit, status=
'keep')
1215 integer(i4b),
intent(in) :: n_unit
1216 real(dp),
intent(out) :: para
1220 #if !defined(ALLOW_OPENAD) /* Normal */ 1223 character(len=1024) :: txt
1227 character :: ch_dummy
1228 character(len=*),
parameter :: fmt1=
'(a)', fmt3=
'(15x)' 1229 logical :: first_read
1233 #endif /* Normal vs. OpenAD */ 1237 do while (check ==
'%')
1239 #if !defined(ALLOW_OPENAD) /* Normal */ 1241 read(unit=n_unit, fmt=
'(a)') txt
1245 if (check /=
'%')
then 1246 i0 = index(txt,
'=') + 1
1247 read(txt(i0:), *) para
1250 #else /* OpenAD: cannot differentiate through index function */ 1252 if (first_read)
then 1253 read(unit=n_unit, fmt=trim(fmt1), advance=
'no') check
1256 read(unit=n_unit, fmt=trim(fmt1)) ch_dummy
1257 read(unit=n_unit, fmt=trim(fmt1), advance=
'no') check
1260 if (check /=
'%')
then 1261 read(unit=n_unit, fmt=trim(fmt3), advance=
'no')
1262 read(unit=n_unit, fmt=*) para
1265 #endif /* Normal vs. OpenAD */ 1271 #if (defined(ALLOW_GRDCHK) || defined(ALLOW_OPENAD)) 1275 subroutine read_ad_data()
1282 integer(i4b) :: ios, i, j, k, kc, kt, KDATA
1283 logical :: debug = .true.
1285 #if (defined(AGE_COST)) /* Only designed for AGE_COST right now */ 1291 kdata = kcmax + ktmax
1296 open(unit=90, iostat=ios, &
1297 file=
'subroutines/openad/gridded_age_obs.dat', status=
'old')
1298 open(unit=91, iostat=ios, &
1299 file=
'subroutines/openad/gridded_age_unc.dat', status=
'old')
1305 read(unit=90, fmt=
'(es11.4)', advance=
"no") age_data(k,j,i)
1306 read(unit=91, fmt=
'(es11.4)', advance=
"no") age_unc(k,j,i)
1308 read(unit=90, fmt=
'(es11.4)' ) age_data(k,j,i)
1309 read(unit=91, fmt=
'(es11.4)' ) age_unc(k,j,i)
1322 age_data(k,j,i) = age_data(k,j,i) * year_sec
1323 age_unc(k,j,i) = age_unc(k,j,i) * year_sec
1325 age_c (k,j,i) = age_data(k,j,i)
1332 #endif /* Only designed for AGE_COST right now */ 1334 end subroutine read_ad_data
1336 #endif /* inclusion of OpenAD only routine */ 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.
integer(i1b), 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) 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.
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) age_c
age_c(kc,j,i): Age in the upper (kc) ice domain
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.
subroutine error(error_message)
Main routine of error_m: Writing of error messages and stopping execution.
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.
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.
logical, dimension(0:jmax, 0:imax) flag_grounding_line_1
flag_grounding_line_1(j,i): Grounding line flag. .true.: grounding line point (grounded ice point wit...
real(dp), dimension(0:jmax, 0:imax) vis_ave_g
vis_ave_g(j,i): Depth-averaged viscosity of the SIA/SSA, at (i,j)
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)
subroutine, public read_erg_nc(z_sl, filename, opt_maske, opt_n_cts, opt_kc_cts, opt_H_cold, opt_H_temp, opt_H, opt_temp_r, opt_omega_t, opt_age_t, opt_temp_c, opt_age_c, opt_omega_c, opt_flag_temp_age_only)
Reading of time-slice files in native binary or NetCDF format.
real(dp), dimension(0:jmax, 0:imax) p_b_w
p_b_w(j,i): Basal water pressure
integer(i1b), 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) vx_m_ssa
vx_m_ssa(j,i): Mean (depth-averaged) SSA velocity in x-direction, at (i+1/2,j)
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)
logical, dimension(0:jmax, 0:imax) flag_shelfy_stream_x
flag_shelfy_stream_x(j,i): Shelfy stream flag in x-direction, at (i+1/2,j). .true.: shelfy stream point .false.: otherwise
logical, dimension(0:jmax, 0:imax) flag_grounded_front_b_1
flag_grounded_front_b_1(j,i): Marine-terminating grounded front flag. .true.: grounded front point (g...
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(i1b), dimension(0:jmax, 0:imax) maske_target
maske_target(j,i): Target topography (ice-land-ocean mask)
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).
real(dp), dimension(0:jmax, 0:imax) vis_int_g
vis_int_g(j,i): Depth-integrated viscosity of the SIA/SSA, at (i,j)
subroutine read_phys_para_value(para, n_unit)
Reading of a value of a physical parameter from the phys_para file.
real(dp), dimension(0:jmax, 0:imax) ratio_sl_y
ratio_sl_y(j,i): Ratio of basal to surface velocity (slip ratio) in y-direction, at (i...
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:jmax, 0:imax) vy_m_ssa
vy_m_ssa(j,i): Mean (depth-averaged) SSA velocity in y-direction, at (i,j+1/2)
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), dimension(0:kcmax, 0:jmax, 0:imax) age_c_neu
(.)_neu: New value of quantity (.) computed during an integration step
real(dp) rho_c_r
RHO_C_R: Density times specific heat of the lithosphere.
Reading of several input data.
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.
integer(i4b), dimension(0:jmax, 0:imax) kc_cts
kc_cts(j,i): Position kc of the CTS (for COLD, ENTC, ENTM)
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)
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.
Writing of error messages and stopping execution.
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
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), parameter eps_dp
eps_dp: Small number to double-precision accuracy
real(dp) kappa_c
KAPPA_C: Heat conductivity of crustal material (dust)
real(dp) rho_i
RHO_I: Density of ice.
real(dp) nue
NUE: Water diffusivity in ice.
real(dp) c_c
C_C: Specific heat of crustal material (dust)
real(dp) beta1_lt_0
BETA1_LT_0: Degree-day factor for snow at low summer temperatures.
logical, dimension(0:jmax, 0:imax) flag_shelfy_stream
flag_shelfy_stream(j,i): Shelfy stream flag on the main grid. .true.: grounded ice, and at least one neighbour on the staggered grid is a shelfy stream point .false.: otherwise
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)
real(dp), dimension(0:jmax, 0:imax) vy_m_sia
vy_m_sia(j,i): Mean (depth-averaged) SIA velocity in y-direction, at (i,j+1/2)
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.
integer(i1b), 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) q_bm
Q_bm(j,i): Basal melting rate.
logical, dimension(0:jmax, 0:imax) flag_grounded_front_b_2
flag_grounded_front_b_2(j,i): Marine-terminating grounded front flag. .true.: grounded front point (o...
real(dp), dimension(0:jmax, 0:imax) q_w
q_w(j,i): Scalar volume flux of the basal water
real(dp), dimension(-190:10) c
C(n): Tabulated values for the specific heat of ice.
real(dp), dimension(0:jmax, 0:imax) q_w_y
q_w_y(j,i): Scalar volume flux of the basal water in y-direction
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.
real(dp), dimension(0:jmax, 0:imax) temph_b
temph_b(j,i): Basal temperature relative to the pressure melting point
logical, dimension(0:jmax, 0:imax) flag_grounded_front_a_2
flag_grounded_front_a_2(j,i): Land-terminating grounded front flag. .true.: grounded front point (ice...
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.
logical, dimension(0:jmax, 0:imax) flag_calving_front_1
flag_calving_front_1(j,i): Calving front flag. .true.: calving front point (floating ice point with a...
real(dp), dimension(0:jmax, 0:imax) ratio_sl_x
ratio_sl_x(j,i): Ratio of basal to surface velocity (slip ratio) in x-direction, at (i+1/2...
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.
logical, dimension(0:jmax, 0:imax) flag_calving_front_2
flag_calving_front_2(j,i): Calving front flag. .true.: calving front point (ocean point with at least...
logical, dimension(0:jmax, 0:imax) flag_shelfy_stream_y
flag_shelfy_stream_y(j,i): Shelfy stream flag in y-direction, at (i,j+1/2). .true.: shelfy stream point .false.: otherwise
real(dp) l
L: Latent heat of ice.
integer(i4b) n_data_kei
n_data_kei: Number of tabulated values of the kei function
real(dp), dimension(0:jmax, 0:imax) q_w_x
q_w_x(j,i): Scalar volume flux of the basal water in x-direction
character(len=256) errormsg
errormsg: Error-message string
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
logical, dimension(0:jmax, 0:imax) flag_grounding_line_2
flag_grounding_line_2(j,i): Grounding line flag. .true.: grounding line point (floating ice point wit...
real(dp), dimension(0:jmax, 0:imax) vx_m_sia
vx_m_sia(j,i): Mean (depth-averaged) SIA velocity in x-direction, at (i+1/2,j)
real(dp), dimension(0:jmax, 0:imax) dzb_dtau
dzb_dtau(j,i): Derivative of zb by tau (time)
subroutine, public read_target_topo_nc(target_topo_dat_name)
Reading of the target-topography file (in NetCDF format).
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.
logical, dimension(0:jmax, 0:imax) flag_grounded_front_a_1
flag_grounded_front_a_1(j,i): Land-terminating grounded front flag. .true.: grounded front point (gro...
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))