SICOPOLIS V5-dev  Revision 1368
boundary_m.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Module : b o u n d a r y _ m
4 !
5 !> @file
6 !!
7 !! Computation of the surface temperature (must be less than 0 deg C!)
8 !! and of the accumulation-ablation function.
9 !!
10 !! @section Copyright
11 !!
12 !! Copyright 2009-2018 Ralf Greve, Eduardo Flandez, Matthias Scheiter
13 !!
14 !! @section License
15 !!
16 !! This file is part of SICOPOLIS.
17 !!
18 !! SICOPOLIS is free software: you can redistribute it and/or modify
19 !! it under the terms of the GNU General Public License as published by
20 !! the Free Software Foundation, either version 3 of the License, or
21 !! (at your option) any later version.
22 !!
23 !! SICOPOLIS is distributed in the hope that it will be useful,
24 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
25 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
26 !! GNU General Public License for more details.
27 !!
28 !! You should have received a copy of the GNU General Public License
29 !! along with SICOPOLIS. If not, see <http://www.gnu.org/licenses/>.
30 !<
31 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
32 
33 !-------------------------------------------------------------------------------
34 !> Computation of the surface temperature (must be less than 0 deg C!)
35 !! and of the accumulation-ablation function.
36 !<------------------------------------------------------------------------------
37 module boundary_m
38 
39  use sico_types_m
41  use sico_vars_m
42  use error_m
43 
44  implicit none
45 
46  public
47 
48 contains
49 
50 !-------------------------------------------------------------------------------
51 !> Main routine of boundary_m:
52 !! Computation of the surface temperature (must be less than 0 deg C!)
53 !! and of the accumulation-ablation function.
54 !<------------------------------------------------------------------------------
55 subroutine boundary(time, dtime, dxi, deta, &
56  delta_ts, glac_index, z_sl, dzsl_dtau, z_mar)
57 
58 #if ((MARGIN==2) \
59  && (marine_ice_formation==2) \
60  && (marine_ice_calving==9))
62 #endif
63 
65 
66 implicit none
67 
68 real(dp), intent(in) :: time, dtime, dxi, deta
69 
70 real(dp), intent(out) :: delta_ts, glac_index, dzsl_dtau, z_mar
71 real(dp), intent(inout) :: z_sl
72 
73 ! Further return variables
74 ! (defined as global variables in module sico_variables_m):
75 !
76 ! accum(j,i), evap(j,i), runoff(j,i), as_perp(j,i),
77 ! calv_grounded(j,i), temp_s(j,i)
78 
79 integer(i4b) :: i, j
80 integer(i4b) :: i_gr, i_kl
81 real(dp) :: z_sl_old
82 real(dp) :: z_sl_min, t1, t2, t3, t4, t5, t6
83 real(dp) :: time_gr, time_kl, asp_a, asp_b
84 real(dp) :: ela_now
85 real(dp), dimension(0:JMAX,0:IMAX) :: ela_dist, asp_dist, grad_dist
86 real(dp), dimension(0:JMAX,0:IMAX) :: dist, dist2
87 logical, dimension(0:JMAX,0:IMAX) :: check_point
88 logical, save :: firstcall = .true.
89 
90 !-------- Initialization of intent(out) variables --------
91 
92 z_sl_old = z_sl
93 
94 delta_ts = 0.0_dp
95 glac_index = 0.0_dp
96 z_sl = 0.0_dp
97 dzsl_dtau = 0.0_dp
98 z_mar = 0.0_dp
99 
100 !-------- Surface-temperature deviation from present values --------
101 
102 #if (TSURFACE==1)
103 delta_ts = delta_ts0
104 ! ! Steady state with prescribed constant
105 ! ! air-temperature deviation
106 #elif (TSURFACE==3)
107 delta_ts = sine_amplit &
108  *cos(2.0_dp*pi*time/(sine_period*year_sec)) &
109  -sine_amplit
110 ! ! Sinusoidal air-temperature forcing
111 #elif (TSURFACE==4)
112 
113 ! ------ delta_ts from the GRIP record
114 
115 if (time/year_sec.lt.real(grip_time_min,dp)) then
116  delta_ts = griptemp(0)
117 else if (time/year_sec.lt.real(grip_time_max,dp)) then
118 
119  i_kl = floor(((time/year_sec) &
120  -real(grip_time_min,dp))/real(grip_time_stp,dp))
121  i_kl = max(i_kl, 0)
122 
123  i_gr = ceiling(((time/year_sec) &
124  -real(grip_time_min,dp))/real(grip_time_stp,dp))
125  i_gr = min(i_gr, ndata_grip)
126 
127  if (i_kl.eq.i_gr) then
128 
129  delta_ts = griptemp(i_kl)
130 
131  else
132 
133  time_kl = (grip_time_min + i_kl*grip_time_stp) *year_sec
134  time_gr = (grip_time_min + i_gr*grip_time_stp) *year_sec
135 
136  delta_ts = griptemp(i_kl) &
137  +(griptemp(i_gr)-griptemp(i_kl)) &
138  *(time-time_kl)/(time_gr-time_kl)
139  ! linear interpolation of the ice-core data
140 
141  end if
142 
143 else
144  delta_ts = griptemp(ndata_grip)
145 end if
146 
147 delta_ts = delta_ts * grip_temp_fact
148 ! ! modification by constant factor
149 
150 #elif (TSURFACE==5)
151 ! Enter here delta_ts scenario:
152 
153 #endif
154 
155 !-------- Sea level --------
156 
157 #if (SEA_LEVEL==1)
158 ! ------ constant sea level
159 z_sl = z_sl0
160 
161 #elif (SEA_LEVEL==2)
162 ! ------ saw-tooth-shaped palaeoclimatic sea-level forcing
163 
164 z_sl_min = -130.0_dp
165 
166 t1 = -250000.0_dp *year_sec
167 t2 = -140000.0_dp *year_sec
168 t3 = -125000.0_dp *year_sec
169 t4 = -21000.0_dp *year_sec
170 t5 = -8000.0_dp *year_sec
171 t6 = 0.0_dp *year_sec
172 
173 if (time.lt.t1) then
174  z_sl = 0.0_dp
175 else if (time.lt.t2) then
176  z_sl = z_sl_min*(time-t1)/(t2-t1)
177 else if (time.lt.t3) then
178  z_sl = -z_sl_min*(time-t3)/(t3-t2)
179 else if (time.lt.t4) then
180  z_sl = z_sl_min*(time-t3)/(t4-t3)
181 else if (time.lt.t5) then
182  z_sl = -z_sl_min*(time-t5)/(t5-t4)
183 else if (time.lt.t6) then
184  z_sl = 0.0_dp
185 else
186  z_sl = 0.0_dp
187 end if
188 
189 #elif (SEA_LEVEL==3)
190 
191 ! ------ z_sl from the SPECMAP record
192 
193 if (time/year_sec.lt.real(specmap_time_min,dp)) then
194  z_sl = specmap_zsl(0)
195 else if (time/year_sec.lt.real(specmap_time_max,dp)) then
196 
197  i_kl = floor(((time/year_sec) &
198  -real(specmap_time_min,dp))/real(specmap_time_stp,dp))
199  i_kl = max(i_kl, 0)
200 
201  i_gr = ceiling(((time/year_sec) &
202  -real(specmap_time_min,dp))/real(specmap_time_stp,dp))
203  i_gr = min(i_gr, ndata_specmap)
204 
205  if (i_kl.eq.i_gr) then
206 
207  z_sl = specmap_zsl(i_kl)
208 
209  else
210 
211  time_kl = (specmap_time_min + i_kl*specmap_time_stp) *year_sec
212  time_gr = (specmap_time_min + i_gr*specmap_time_stp) *year_sec
213 
214  z_sl = specmap_zsl(i_kl) &
215  +(specmap_zsl(i_gr)-specmap_zsl(i_kl)) &
216  *(time-time_kl)/(time_gr-time_kl)
217  ! linear interpolation of the sea-level data
218 
219  end if
220 
221 else
222  z_sl = specmap_zsl(ndata_specmap)
223 end if
224 
225 #endif
226 
227 ! ------ Time derivative of the sea level
228 
229 if ( z_sl_old > -999999.9_dp ) then
230  dzsl_dtau = (z_sl-z_sl_old)/dtime
231 else ! only dummy value for z_sl_old available
232  dzsl_dtau = 0.0_dp
233 end if
234 
235 ! ------ Minimum bedrock elevation for extent of marine ice
236 
237 #if (MARGIN==2)
238 
239 #if (MARINE_ICE_CALVING==2 || MARINE_ICE_CALVING==3)
240 z_mar = z_mar
241 #elif (MARINE_ICE_CALVING==4 || MARINE_ICE_CALVING==5)
242 z_mar = fact_z_mar*z_sl
243 #elif (MARINE_ICE_CALVING==6 || MARINE_ICE_CALVING==7)
244 if (z_sl >= -80.0_dp) then
245  z_mar = 2.5_dp*z_sl
246 else
247  z_mar = 10.25_dp*(z_sl+80.0_dp)-200.0_dp
248 end if
249 z_mar = fact_z_mar*z_mar
250 #endif
251 
252 #endif
253 
254 ! ------ Update of the mask according to the sea level
255 
256 ! ---- Check all sea and floating-ice points and their direct
257 ! neighbours
258 
259 do i=0, imax
260 do j=0, jmax
261  check_point(j,i) = .false.
262 end do
263 end do
264 
265 do i=1, imax-1
266 do j=1, jmax-1
267  if (maske(j,i) >= 2_i1b) then
268  check_point(j ,i ) = .true.
269  check_point(j ,i+1) = .true.
270  check_point(j ,i-1) = .true.
271  check_point(j+1,i ) = .true.
272  check_point(j-1,i ) = .true.
273  end if
274 end do
275 end do
276 
277 do i=1, imax-1
278 do j=1, jmax-1
279  if (check_point(j,i)) then
280  maske_neu(j,i) = mask_update_sea_level(z_sl, i, j)
281  end if
282 end do
283 end do
284 
285 ! ---- Assign new values of the mask
286 
287 do i=1, imax-1
288 do j=1, jmax-1
289  if (check_point(j,i)) then
290  maske(j,i) = maske_neu(j,i)
291  end if
292 end do
293 end do
294 
295 !-------- Surface air temperature and
296 ! accumulation-ablation function (SMB) --------
297 
298 #if (!defined(SURFACE_FORCING) || SURFACE_FORCING==1)
299 
300 errormsg = ' >>> boundary: Option SURFACE_FORCING=1 not supported any more!'
301 call error(errormsg)
302 
303 #elif (SURFACE_FORCING==2)
304 
305 do i=0, imax
306 do j=0, jmax
307 
308  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
309  temp_s(j,i) = temp_s(j,i) + delta_ts
310  ! Correction with temperature deviation delta_ts
311  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
312  ! Cut-off of positive air temperatures
313 
314  accum_present(j,i) = s_0
315  accum(j,i) = accum_present(j,i)
316 
317  ela_now = ela + dela_dts * delta_ts
318 
319  as_perp(j,i) = m_0*(zs(j,i)-ela_now)
320  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
321 
322 end do
323 end do
324 
325 #elif (SURFACE_FORCING==3)
326 
327 do i=0, imax
328 do j=0, jmax
329 
330  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
331  temp_s(j,i) = temp_s(j,i) + delta_ts
332  ! Correction with temperature deviation delta_ts
333  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
334  ! Cut-off of positive air temperatures
335 
336  accum_present(j,i) = s_0
337  accum(j,i) = accum_present(j,i)
338 
339  ela_now = ela + dela_dts * delta_ts
340 
341  if (xi(i) /= x_gip .or. eta(j) /= y_gip) then
342  ela_dist(j,i) = ela_now &
343  + ela_amp &
344  * sin(atan2((eta(j)-y_gip), (xi(i)-x_gip)) + phi_0)
345  else
346  ela_dist(j,i) = ela_now
347  end if
348 
349  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
350 
351  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
352 
353 end do
354 end do
355 
356 #elif (SURFACE_FORCING==4)
357 
358 do i=0, imax
359 do j=0, jmax
360 
361  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
362  temp_s(j,i) = temp_s(j,i) + delta_ts
363  ! Correction with temperature deviation delta_ts
364  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
365  ! Cut-off of positive air temperatures
366 
367  accum_present(j,i) = s_0
368  accum(j,i) = accum_present(j,i)
369 
370  ela_now = ela + dela_dts * delta_ts
371 
372  if (xi(i) /= x_gip .or. eta(j) /= y_gip) then
373  ela_dist(j,i) = ela_now &
374  + ela_amp &
375  * sin(atan2((eta(j)-y_gip), (xi(i)-x_gip)) + phi_0)
376  else
377  ela_dist(j,i) = ela_now
378  end if
379 
380  if (zs(j,i) <= z_gc) then
381  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
382  else ! (zs(j,i) > z_gc)
383  as_perp(j,i) = m_1*(zs(j,i)-ela_dist(j,i))+m_0*(z_gc-ela_dist(j,i))
384  end if
385 
386 end do
387 end do
388 
389 #elif (SURFACE_FORCING==5)
390 
391 asp_a = 0.0_dp
392 asp_b = 0.0_dp
393 
394 do i=0, imax
395 do j=0, jmax
396 
397  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
398  temp_s(j,i) = temp_s(j,i) + delta_ts
399  ! Correction with temperature deviation delta_ts
400  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
401  ! Cut-off of positive air temperatures
402 
403  accum_present(j,i) = s_0
404  accum(j,i) = accum_present(j,i)
405 
406  asp_a = zs(j,i+1) - zs(j,i-1)
407  asp_b = zs(j-1,i) - zs(j+1,i)
408 
409  if (asp_b < 0.0_dp) then
410  asp_dist(j,i) = mod(atan(asp_a/asp_b), 2.0_dp*pi)
411  else ! (asp_b >= 0.0_dp)
412  asp_dist(j,i) = atan(asp_a/asp_b) + pi
413  !!! singular for asp_b == 0.0_dp !!!
414  end if
415 
416  ela_now = ela + dela_dts * delta_ts
417 
418  ela_dist(j,i) = ela_now + ela_amp * sin(asp_dist(j,i) - phi_0 + pi)
419 
420  if (zs(j,i) <= z_gc) then
421  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
422  else ! (zs(j,i) > z_gc)
423  as_perp(j,i) = m_1*(zs(j,i)-ela_dist(j,i))+m_0*(z_gc-ela_dist(j,i))
424  end if
425 
426 end do
427 end do
428 
429 #elif (SURFACE_FORCING==6)
430 
431 asp_a = 0.0_dp
432 asp_b = 0.0_dp
433 
434 do i=0, imax
435 do j=0, jmax
436 
437  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
438  temp_s(j,i) = temp_s(j,i) + delta_ts
439  ! Correction with temperature deviation delta_ts
440  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
441  ! Cut-off of positive air temperatures
442 
443  accum_present(j,i) = s_0
444  accum(j,i) = accum_present(j,i)
445 
446  asp_a = zs(j,i+1) - zs(j,i-1)
447  asp_b = zs(j-1,i) - zs(j+1,i)
448 
449  if (asp_b < 0.0_dp) then
450  asp_dist(j,i) = mod(atan(asp_a/asp_b), 2.0_dp*pi)
451  else ! (asp_b >= 0.0_dp)
452  asp_dist(j,i) = atan(asp_a/asp_b) + pi
453  !!! singular for asp_b == 0.0_dp !!!
454  end if
455 
456  grad_dist(j,i) = sqrt( (asp_a/(2.0_dp*dxi))**2 &
457  + (asp_b/(2.0_dp*deta))**2 ) * 45.0_dp
458 
459  ela_now = ela + dela_dts * delta_ts
460 
461  ela_dist(j,i) = ela_now + ela_amp * sin(asp_dist(j,i) - phi_0 + pi)
462  if (grad_dist(j,i) > tgt) ela_dist(j,i) = ela_dist(j,i) + grad_dist(j,i)*2
463 
464  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
465 
466  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
467 
468 end do
469 end do
470 
471 #elif (SURFACE_FORCING==7)
472 
473 asp_a = 0.0_dp
474 asp_b = 0.0_dp
475 
476 do i=0, imax
477 do j=0, jmax
478 
479  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
480  temp_s(j,i) = temp_s(j,i) + delta_ts
481  ! Correction with temperature deviation delta_ts
482  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
483  ! Cut-off of positive air temperatures
484 
485  accum_present(j,i) = s_0
486  accum(j,i) = accum_present(j,i)
487 
488  asp_a = zs(j,i+1) - zs(j,i-1)
489  asp_b = zs(j-1,i) - zs(j+1,i)
490 
491  if (asp_b < 0.0_dp) then
492  asp_dist(j,i) = mod(atan(asp_a/asp_b), 2.0_dp*pi)
493  else ! (asp_b >= 0.0_dp)
494  asp_dist(j,i) = atan(asp_a/asp_b) + pi
495  !!! singular for asp_b == 0.0_dp !!!
496  end if
497 
498  ela_now = ela + dela_dts * delta_ts
499 
500  ela_dist(j,i) = ela_now + ela_amp * sin(asp_dist(j,i) - phi_0 + pi)
501 
502  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
503 
504  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
505 
506 end do
507 end do
508 
509 #elif (SURFACE_FORCING==8)
510 
511 do i=0, imax
512 do j=0, jmax
513  dist(j,i) = sqrt( (xi(i)-x_gip)**2 + (eta(j)-y_gip)**2 )
514  dist2(j,i) = sqrt( (xi(i)-x_gip2)**2 + (eta(j)-y_gip2)**2 )
515 end do
516 end do
517 
518 do i=0, imax
519 do j=0, jmax
520 
521  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
522  temp_s(j,i) = temp_s(j,i) + delta_ts
523  ! Correction with temperature deviation delta_ts
524  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
525  ! Cut-off of positive air temperatures
526 
527  accum_present(j,i) = s_0
528  accum(j,i) = accum_present(j,i)
529 
530  ela_now = ela + dela_dts * delta_ts
531 
532  if (dist(j,i) < 3000.0_dp) then
533  if (xi(i) /= x_gip .or. eta(j) /= y_gip) then
534  ela_dist(j,i) = ela_now &
535  + ela_amp &
536  * sin(atan2((eta(j)-y_gip), (xi(i)-x_gip)) + phi_0)
537  else
538  ela_dist(j,i) = ela_now
539  end if
540  else
541  ela_dist(j,i) = ela_now
542  end if
543 
544  if (dist2(j,i) < 1500.0_dp) then
545  if (xi(i) /= x_gip2 .or. eta(j) /= y_gip2) then
546  ela_dist(j,i) = ela_dist(j,i) &
547  + ela_amp2 &
548  * sin(atan2((eta(j)-y_gip2), (xi(i)-x_gip2)) + phi_0)
549  end if
550  end if
551 
552  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
553 
554  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
555 
556 end do
557 end do
558 
559 #endif
560 
561 runoff = accum - as_perp
562 evap = 0.0_dp
563 
564 ! ------ Prescribed SMB correction
565 
568 
569 !-------- Calving rate of grounded ice --------
570 
571 calv_grounded = 0.0_dp
572 
573 #if ((MARGIN==2) \
574  && (marine_ice_formation==2) \
575  && (marine_ice_calving==9))
576 
577 call calving_underwater_ice(z_sl)
579 
580 #endif
581 
582 if (firstcall) firstcall = .false.
583 
584 end subroutine boundary
585 
586 !-------------------------------------------------------------------------------
587 
588 end module boundary_m
589 !
integer(i4b) grip_time_stp
grip_time_stp: Time step of the data values for the surface temperature anomaly
integer(i4b) ndata_specmap
ndata_specmap: Number of data values for the sea level
subroutine error(error_message)
Main routine of error_m: Writing of error messages and stopping execution.
Definition: error_m.F90:47
integer(i4b) specmap_time_max
specmap_time_max: Maximum time of the data values for the sea level
real(dp), dimension(0:jmax, 0:imax) calv_grounded
calv_grounded(j,i): Calving rate of grounded ice
integer(i4b) specmap_time_min
specmap_time_min: Minimum time of the data values for the sea level
integer(i1b), dimension(0:jmax, 0:imax) maske_neu
maske_neu(j,i): New value of maske computed during an integration step
real(dp), dimension(0:jmax, 0:imax) as_perp
as_perp(j,i): Accumulation-ablation function at the ice surface (SMB)
real(dp), dimension(0:jmax, 0:imax) runoff
runoff(j,i): Runoff rate at the ice surface
real(dp), dimension(0:jmax, 0:imax) smb_corr_in
smb_corr_in(j,i): Prescribed SMB correction
Definition: sico_vars_m.F90:43
real(dp), dimension(0:jmax, 0:imax) temp_s
temp_s(j,i): Ice surface temperature
integer(i4b) grip_time_max
grip_time_max: Maximum time of the data values for the surface temperature anomaly ...
Declarations of global variables for SICOPOLIS (for the ANT domain).
Definition: sico_vars_m.F90:35
integer(i4b) ndata_grip
ndata_grip: Number of data values for the surface temperature anomaly
real(dp), dimension(0:jmax, 0:imax) calv_uw_ice
Calving of "underwater ice".
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), dimension(0:jmax, 0:imax) evap
evap(j,i): Evaporation rate at the ice surface
real(dp), dimension(0:jmax, 0:imax) accum
accum(j,i): Accumulation rate at the ice surface (includes liquid precipitation = rainfall!) ...
integer(i1b) function, public mask_update_sea_level(z_sl, i, j)
Main function of mask_update_m: Update of the ice-land-ocean mask due to changes of the sea level...
real(dp), dimension(:), allocatable specmap_zsl
specmap_zsl(n): Data values for the sea level
Writing of error messages and stopping execution.
Definition: error_m.F90:35
real(dp), dimension(0:jmax, 0:imax) accum_present
accum_present(j,i): Present-day accumulation rate at the ice surface (for EISMINT, ISMIP HEINO and the north and south polar caps of Mars)
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 ...
subroutine boundary(time, dtime, dxi, deta, delta_ts, glac_index, z_sl, dzsl_dtau, z_mar)
Main routine of boundary_m: Computation of the surface temperature (must be less than 0 deg C!) and o...
Definition: boundary_m.F90:61
integer(i4b) specmap_time_stp
specmap_time_stp: Time step of the data values for the sea level
real(dp), parameter pi
pi: Constant pi
integer(i4b) grip_time_min
grip_time_min: Minimum time of the data values for the surface temperature anomaly ...
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
subroutine calving_underwater_ice(z_sl)
Main routine: Calving of "underwater ice".
Computation of the surface temperature (must be less than 0 deg C!) and of the accumulation-ablation ...
Definition: boundary_m.F90:37
character(len=256) errormsg
errormsg: Error-message string
real(dp), dimension(0:jmax, 0:imax) smb_corr_prescribed
smb_corr_prescribed(j,i): Prescribed SMB correction
Update of the ice-land-ocean mask due to changes of the sea level.
integer, parameter dp
Double-precision reals.
Declarations of global variables for SICOPOLIS.
real(dp), dimension(:), allocatable griptemp
griptemp(n): Data values for the surface temperature anomaly