SICOPOLIS V5-dev  Revision 1173
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-2017 Ralf Greve
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 
43  implicit none
44 
45  public
46 
47 contains
48 
49 !-------------------------------------------------------------------------------
50 !> Main routine of boundary_m:
51 !! Computation of the surface temperature (must be less than 0 deg C!)
52 !! and of the accumulation-ablation function.
53 !<------------------------------------------------------------------------------
54 subroutine boundary(time, dtime, dxi, deta, &
55  delta_ts, glac_index, z_sl, dzsl_dtau, z_mar)
56 
57 #if ((MARGIN==2) \
58  && (marine_ice_formation==2) \
59  && (marine_ice_calving==9))
61 #endif
62 
64  use pdd_m
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, n
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
84 real(dp) :: z_sle_present, z_sle_help
85 real(dp), dimension(0:JMAX,0:IMAX,0:12) :: precip
86 real(dp), dimension(0:JMAX,0:IMAX) :: &
87  melt, melt_star
88 !real(dp), dimension(0:JMAX,0:IMAX) :: &
89 ! snowfall, rainfall, melt, melt_star
90 !snowfall, rainfall declared in sico_variables because of output5.F90
91 real(dp), dimension(0:JMAX,0:IMAX,12) :: temp_mm
92 real(dp), dimension(0:JMAX,0:IMAX) :: temp_ma
93 real(dp), dimension(12) :: temp_mm_help
94 real(dp) :: temp_jja_help
95 real(dp), dimension(0:JMAX,0:IMAX) :: et
96 real(dp) :: gamma_t, temp_diff
97 real(dp) :: gamma_p, zs_thresh, &
98  temp_rain, temp_snow, &
99  inv_delta_temp_rain_snow, coeff(0:5), inv_sqrt2_s_stat, &
100  precip_fact, frac_solid
101 real(dp) :: s_stat, beta1, beta2, pmax, mu, lambda_lti, temp_lti
102 logical, dimension(0:JMAX,0:IMAX) :: check_point
103 logical, save :: firstcall = .true.
104 
105 real(dp), parameter :: &
106  inv_twelve = 1.0_dp/12.0_dp, one_third = 1.0_dp/3.0_dp
107 
108 
109 !-------- Initialization of intent(out) variables --------
110 
111 z_sl_old = z_sl
112 
113 delta_ts = 0.0_dp
114 glac_index = 0.0_dp
115 z_sl = 0.0_dp
116 dzsl_dtau = 0.0_dp
117 z_mar = 0.0_dp
118 
119 !-------- Surface-temperature deviation from present values --------
120 
121 #if TSURFACE==1
122 delta_ts = delta_ts0
123 
124 ! ! Steady state with prescribed constant
125 ! ! air-temperature deviation
126 #elif TSURFACE==3
127 delta_ts = sine_amplit &
128  *cos(2.0_dp*pi*time/(sine_period*year_sec)) &
129  -sine_amplit
130 ! ! Sinusoidal air-temperature forcing
131 #elif TSURFACE==4
132 
133 ! ------ delta_ts from ice-core record
134 
135 if (time/year_sec.lt.real(grip_time_min,dp)) then
136  delta_ts = griptemp(0)
137 else if (time/year_sec.lt.real(grip_time_max,dp)) then
138 
139  i_kl = floor(((time/year_sec)-real(grip_time_min,dp))/real(grip_time_stp,dp))
140  i_kl = max(i_kl, 0)
141 
142  i_gr = ceiling(((time/year_sec)-real(grip_time_min,dp))/real(grip_time_stp,dp))
143  i_gr = min(i_gr, ndata_grip)
144 
145  if (i_kl.eq.i_gr) then
146 
147  delta_ts = griptemp(i_kl)
148 
149  else
150 
151  time_kl = (grip_time_min + i_kl*grip_time_stp) *year_sec
152  time_gr = (grip_time_min + i_gr*grip_time_stp) *year_sec
153 
154  delta_ts = griptemp(i_kl) &
155  +(griptemp(i_gr)-griptemp(i_kl)) &
156  *(time-time_kl)/(time_gr-time_kl)
157  ! linear interpolation of the ice-core data
158 
159  end if
160 
161 else
162  delta_ts = griptemp(ndata_grip)
163 end if
164 
165 delta_ts = delta_ts * grip_temp_fact
166 ! ! modification by constant factor
167 
168 !-------- Glacial index --------
169 
170 #elif TSURFACE==5
171 
172 if (time/year_sec < real(gi_time_min,dp)) then
173  glac_index = glacial_index(0)
174 else if (time/year_sec < real(gi_time_max,dp)) then
175 
176  i_kl = floor(((time/year_sec)-real(gi_time_min,dp))/real(gi_time_stp,dp))
177  i_kl = max(i_kl, 0)
178 
179  i_gr = ceiling(((time/year_sec)-real(gi_time_min,dp))/real(gi_time_stp,dp))
180  i_gr = min(i_gr, ndata_gi)
181 
182  if (i_kl == i_gr) then
183 
184  glac_index = glacial_index(i_kl)
185 
186  else
187 
188  time_kl = (gi_time_min + i_kl*gi_time_stp) *year_sec
189  time_gr = (gi_time_min + i_gr*gi_time_stp) *year_sec
190 
191  glac_index = glacial_index(i_kl) &
192  +(glacial_index(i_gr)-glacial_index(i_kl)) &
193  *(time-time_kl)/(time_gr-time_kl)
194  ! linear interpolation of the glacial-index data
195 
196  end if
197 
198 else
199  glac_index = glacial_index(ndata_gi)
200 end if
201 
202 #endif
203 
204 !-------- Sea level --------
205 
206 #if SEA_LEVEL==1
207 ! ------ constant sea level
208 z_sl = z_sl0
209 
210 #elif SEA_LEVEL==2
211 ! ------ saw-tooth-shaped palaeoclimatic sea-level forcing
212 
213 z_sl_min = -130.0_dp
214 
215 t1 = -250000.0_dp *year_sec
216 t2 = -140000.0_dp *year_sec
217 t3 = -125000.0_dp *year_sec
218 t4 = -21000.0_dp *year_sec
219 t5 = -8000.0_dp *year_sec
220 t6 = 0.0_dp *year_sec
221 
222 if (time.lt.t1) then
223  z_sl = 0.0_dp
224 else if (time.lt.t2) then
225  z_sl = z_sl_min*(time-t1)/(t2-t1)
226 else if (time.lt.t3) then
227  z_sl = -z_sl_min*(time-t3)/(t3-t2)
228 else if (time.lt.t4) then
229  z_sl = z_sl_min*(time-t3)/(t4-t3)
230 else if (time.lt.t5) then
231  z_sl = -z_sl_min*(time-t5)/(t5-t4)
232 else if (time.lt.t6) then
233  z_sl = 0.0_dp
234 else
235  z_sl = 0.0_dp
236 end if
237 
238 #elif SEA_LEVEL==3
239 
240 ! ------ z_sl from the SPECMAP record
241 
242 if (time/year_sec.lt.real(specmap_time_min,dp)) then
243  z_sl = specmap_zsl(0)
244 else if (time/year_sec.lt.real(specmap_time_max,dp)) then
245 
246  i_kl = floor(((time/year_sec)-real(specmap_time_min,dp))/real(specmap_time_stp,dp))
247  i_kl = max(i_kl, 0)
248 
249  i_gr = ceiling(((time/year_sec)-real(specmap_time_min,dp))/real(specmap_time_stp,dp))
250  i_gr = min(i_gr, ndata_specmap)
251 
252  if (i_kl.eq.i_gr) then
253 
254  z_sl = specmap_zsl(i_kl)
255 
256  else
257 
258  time_kl = (specmap_time_min + i_kl*specmap_time_stp) *year_sec
259  time_gr = (specmap_time_min + i_gr*specmap_time_stp) *year_sec
260 
261  z_sl = specmap_zsl(i_kl) &
262  +(specmap_zsl(i_gr)-specmap_zsl(i_kl)) &
263  *(time-time_kl)/(time_gr-time_kl)
264  ! linear interpolation of the sea-level data
265 
266  end if
267 
268 else
269  z_sl = specmap_zsl(ndata_specmap)
270 end if
271 
272 #endif
273 
274 ! ------ Time derivative of the sea level
275 
276 if ( z_sl_old > -999999.9_dp ) then
277  dzsl_dtau = (z_sl-z_sl_old)/dtime
278 else ! only dummy value for z_sl_old available
279  dzsl_dtau = 0.0_dp
280 end if
281 
282 ! ------ Minimum bedrock elevation for extent of marine ice
283 
284 #if MARGIN==2
285 
286 #if ( MARINE_ICE_CALVING==2 || MARINE_ICE_CALVING==3 )
287 z_mar = z_mar
288 #elif ( MARINE_ICE_CALVING==4 || MARINE_ICE_CALVING==5 )
289 z_mar = fact_z_mar*z_sl
290 #elif ( MARINE_ICE_CALVING==6 || MARINE_ICE_CALVING==7 )
291 if (z_sl >= -80.0_dp) then
292  z_mar = 2.5_dp*z_sl
293 else
294  z_mar = 10.25_dp*(z_sl+80.0_dp)-200.0_dp
295 end if
296 z_mar = fact_z_mar*z_mar
297 #endif
298 
299 #endif
300 
301 ! ------ Update of the mask according to the sea level
302 
303 ! ---- Check all sea and floating-ice points and their direct
304 ! neighbours
305 
306 do i=0, imax
307 do j=0, jmax
308  check_point(j,i) = .false.
309 end do
310 end do
311 
312 do i=1, imax-1
313 do j=1, jmax-1
314  if (maske(j,i).ge.2) then
315  check_point(j ,i ) = .true.
316  check_point(j ,i+1) = .true.
317  check_point(j ,i-1) = .true.
318  check_point(j+1,i ) = .true.
319  check_point(j-1,i ) = .true.
320  end if
321 end do
322 end do
323 
324 do i=1, imax-1
325 do j=1, jmax-1
326  if (check_point(j,i)) then
327  maske_neu(j,i) = mask_update_sea_level(z_sl, i, j)
328  end if
329 end do
330 end do
331 
332 ! ---- Assign new values of the mask
333 
334 do i=1, imax-1
335 do j=1, jmax-1
336  if (check_point(j,i)) then
337  maske(j,i) = maske_neu(j,i)
338  end if
339 end do
340 end do
341 
342 !-------- Surface air temperatures --------
343 
344 gamma_t = -4.5e-03_dp ! atmospheric lapse rate
345 
346 do i=0, imax
347 do j=0, jmax
348 
349 #if (TSURFACE <= 4)
350 
351 ! ------ Correction of present monthly temperatures with elevation changes
352 ! and temperature deviation delta_ts
353 
354  temp_diff = gamma_t*(zs(j,i)-zs_ref(j,i)) + delta_ts
355 
356  do n=1, 12 ! month counter
357  temp_mm(j,i,n) = temp_mm_present(j,i,n) + temp_diff
358  end do
359 
360 #elif (TSURFACE == 5)
361 
362 ! ------ Correction of present monthly temperatures with LGM anomaly and
363 ! glacial index as well as elevation changes
364 
365  temp_diff = gamma_t*(zs(j,i)-zs_ref(j,i))
366 
367  do n=1, 12 ! month counter
368  temp_mm(j,i,n) = temp_mm_present(j,i,n) &
369  + glac_index*temp_mm_lgm_anom(j,i,n) &
370  + temp_diff
371  end do
372 
373 #endif
374 
375 ! ------ Mean annual air temperature
376 
377  temp_ma(j,i) = 0.0_dp ! initialisation value
378 
379  do n=1, 12 ! month counter
380  temp_ma(j,i) = temp_ma(j,i) + temp_mm(j,i,n)*inv_twelve
381  end do
382 
383 end do
384 end do
385 
386 !-------- Accumulation-ablation function as_perp --------
387 
388 #if (ELEV_DESERT == 1)
389 
390 gamma_p = gamma_p*1.0e-03_dp ! Precipitation lapse rate
391  ! for elevation desertification, in m^(-1)
392 zs_thresh = zs_thresh ! Elevation threshold, in m
393 
394 #endif
395 
396 #if (SOLID_PRECIP == 1) /* Marsiat (1994) */
397 
398 temp_rain = 7.0_dp ! Threshold monthly mean temperature for
399  ! precipitation = 100% rain, in deg C
400 temp_snow = -10.0_dp ! Threshold monthly mean temperature for &
401  ! precipitation = 100% snow, in deg C
402 
403 inv_delta_temp_rain_snow = 1.0_dp/(temp_rain-temp_snow)
404 
405 #elif (SOLID_PRECIP == 2) /* Bales et al. (2009) */
406 
407 temp_rain = 7.2_dp ! Threshold monthly mean temperature for &
408  ! precipitation = 100% rain, in deg C
409 temp_snow = -11.6_dp ! Threshold monthly mean temperature for &
410  ! precipitation = 100% snow, in deg C
411 
412 coeff(0) = 5.4714e-01_dp ! Coefficients
413 coeff(1) = -9.1603e-02_dp ! of
414 coeff(2) = -3.314e-03_dp ! the
415 coeff(3) = 4.66e-04_dp ! fifth-order
416 coeff(4) = 3.8e-05_dp ! polynomial
417 coeff(5) = 6.0e-07_dp ! fit
418 
419 #elif (SOLID_PRECIP == 3) /* Huybrechts and de Wolde (1999) */
420 
421 temp_rain = 2.0_dp ! Threshold instantaneous temperature for &
422  ! precipitation = 100% rain, in deg C
423 temp_snow = temp_rain ! Threshold instantaneous temperature for &
424  ! precipitation = 100% snow, in deg C
425 
426 s_stat = s_stat_0 ! Standard deviation of the air termperature
427  ! (same parameter as in the PDD model)
428 
429 inv_sqrt2_s_stat = 1.0_dp/(sqrt(2.0_dp)*s_stat)
430 
431 #endif
432 
433 #if (ABLSURFACE==1 || ABLSURFACE==2)
434 
435 s_stat = s_stat_0
436 beta1 = beta1_0 *(0.001_dp/86400.0_dp)*(rho_w/rho)
437  ! (mm WE)/(d*deg C) --> (m IE)/(s*deg C)
438 beta2 = beta2_0 *(0.001_dp/86400.0_dp)*(rho_w/rho)
439  ! (mm WE)/(d*deg C) --> (m IE)/(s*deg C)
440 
441 pmax = pmax_0
442 mu = mu_0 *(1000.0_dp*86400.0_dp)*(rho/rho_w)
443  ! (d*deg C)/(mm WE) --> (s*deg C)/(m IE)
444 
445 #elif (ABLSURFACE==3)
446 
447 lambda_lti = lambda_lti *(0.001_dp/year_sec)*(rho_w/rho)
448  ! (mm WE)/(a*deg C) --> (m IE)/(s*deg C)
449 temp_lti = temp_lti
450 mu = 0.0_dp ! no superimposed ice considered
451 
452 #endif
453 
454 do i=0, imax
455 do j=0, jmax
456 
457 ! ------ Accumulation
458 
459 #if (ACCSURFACE <= 3)
460 
461 ! ---- Elevation desertification of precipitation
462 
463 #if (ELEV_DESERT == 0)
464 
465  precip_fact = 1.0_dp ! no elevation desertification
466 
467 #elif (ELEV_DESERT == 1)
468 
469  if (zs_ref(j,i) < zs_thresh) then
470  precip_fact &
471  = exp(gamma_p*(max(zs(j,i),zs_thresh)-zs_thresh))
472  else
473  precip_fact &
474  = exp(gamma_p*(max(zs(j,i),zs_thresh)-zs_ref(j,i)))
475  end if
476 
477 #else
478  stop ' >>> boundary: Parameter ELEV_DESERT must be either 0 or 1!'
479 #endif
480 
481  do n=1, 12 ! month counter
482  precip(j,i,n) = precip_present(j,i,n)*precip_fact
483  end do
484 
485 #endif
486 
487 ! ---- Precipitation change related to changing climate
488 
489 #if ACCSURFACE==1
490  precip_fact = accfact
491 #elif ACCSURFACE==2
492  precip_fact = 1.0_dp + gamma_s*delta_ts
493 #elif ACCSURFACE==3
494  precip_fact = exp(gamma_s*delta_ts)
495 #endif
496 
497 #if (ACCSURFACE <= 3)
498 
499  precip(j,i,0) = 0.0_dp ! initialisation value for mean annual precip
500 
501  do n=1, 12 ! month counter
502  precip(j,i,n) = precip(j,i,n)*precip_fact ! monthly precip
503  precip(j,i,0) = precip(j,i,0) + precip(j,i,n)*inv_twelve
504  ! mean annual precip
505  end do
506 
507 #elif (ACCSURFACE == 5)
508 
509  precip(j,i,0) = 0.0_dp ! initialisation value for mean annual precip
510 
511  do n=1, 12 ! month counter
512 
513 #if (PRECIP_ANOM_INTERPOL==1)
514  precip_fact = 1.0_dp-glac_index+glac_index*precip_lgm_anom(j,i,n)
515  ! interpolation with a linear function
516 #elif (PRECIP_ANOM_INTERPOL==2)
517  precip_fact = exp(-glac_index*gamma_precip_lgm_anom(j,i,n))
518  ! interpolation with an exponential function
519 #endif
520 
521  precip(j,i,n) = precip_present(j,i,n)*precip_fact ! monthly precip
522  precip(j,i,0) = precip(j,i,0) + precip(j,i,n)*inv_twelve
523  ! mean annual precip
524  end do
525 
526 #endif
527 
528 ! ---- Annual accumulation, snowfall and rainfall rates
529 
530  accum(j,i) = precip(j,i,0)
531 
532  snowfall(j,i) = 0.0_dp ! initialisation value
533 
534  do n=1, 12 ! month counter
535 
536 #if (SOLID_PRECIP == 1) /* Marsiat (1994) */
537 
538  if (temp_mm(j,i,n) >= temp_rain) then
539  frac_solid = 0.0_dp
540  else if (temp_mm(j,i,n) <= temp_snow) then
541  frac_solid = 1.0_dp
542  else
543  frac_solid = (temp_rain-temp_mm(j,i,n))*inv_delta_temp_rain_snow
544  end if
545 
546 #elif (SOLID_PRECIP == 2) /* Bales et al. (2009) */
547 
548  if (temp_mm(j,i,n) >= temp_rain) then
549  frac_solid = 0.0_dp
550  else if (temp_mm(j,i,n) <= temp_snow) then
551  frac_solid = 1.0_dp
552  else
553  frac_solid = coeff(0) + temp_mm(j,i,n) * ( coeff(1) &
554  + temp_mm(j,i,n) * ( coeff(2) &
555  + temp_mm(j,i,n) * ( coeff(3) &
556  + temp_mm(j,i,n) * ( coeff(4) &
557  + temp_mm(j,i,n) * coeff(5) ) ) ) )
558  ! evaluation of 5th-order polynomial by Horner scheme
559  end if
560 
561 #elif (SOLID_PRECIP == 3) /* Huybrechts and de Wolde (1999) */
562 
563  frac_solid = 1.0_dp &
564  - 0.5_dp*erfc((temp_rain-temp_mm(j,i,n))*inv_sqrt2_s_stat)
565 
566 #endif
567 
568  snowfall(j,i) = snowfall(j,i) + precip(j,i,n)*frac_solid*inv_twelve
569 
570  end do
571 
572  rainfall(j,i) = precip(j,i,0) - snowfall(j,i)
573 
574  if (snowfall(j,i) < 0.0_dp) snowfall(j,i) = 0.0_dp ! correction of
575  if (rainfall(j,i) < 0.0_dp) rainfall(j,i) = 0.0_dp ! negative values
576 
577 ! ------ Ablation
578 
579 ! ---- Runoff
580 
581 #if (ABLSURFACE==1 || ABLSURFACE==2)
582 
583 ! ---- Temperature excess ET
584 
585  do n=1, 12 ! month counter
586  temp_mm_help(n) = temp_mm(j,i,n)
587  end do
588 
589  call pdd(temp_mm_help, s_stat, et(j,i))
590 
591 
592 ! ---- Formation rate of superimposed ice (melt_star), melt rate (melt)
593 ! and runoff rate (runoff)
594 
595 #if (ABLSURFACE==1)
596 
597  if ((beta1*et(j,i)) <= (pmax*snowfall(j,i))) then
598  melt_star(j,i) = beta1*et(j,i)
599  melt(j,i) = 0.0_dp
600  runoff(j,i) = melt(j,i)+rainfall(j,i)
601  else
602  melt_star(j,i) = pmax*snowfall(j,i)
603  melt(j,i) = beta2*(et(j,i)-melt_star(j,i)/beta1)
604  runoff(j,i) = melt(j,i)+rainfall(j,i)
605  end if
606 
607 #elif (ABLSURFACE==2)
608 
609  if ( rainfall(j,i) <= (pmax*snowfall(j,i)) ) then
610 
611  if ( (rainfall(j,i)+beta1*et(j,i)) <= (pmax*snowfall(j,i)) ) then
612  melt_star(j,i) = rainfall(j,i)+beta1*et(j,i)
613  melt(j,i) = 0.0_dp
614  runoff(j,i) = melt(j,i)
615  else
616  melt_star(j,i) = pmax*snowfall(j,i)
617  melt(j,i) = beta2 &
618  *(et(j,i)-(melt_star(j,i)-rainfall(j,i))/beta1)
619  runoff(j,i) = melt(j,i)
620  end if
621 
622  else
623 
624  melt_star(j,i) = pmax*snowfall(j,i)
625  melt(j,i) = beta2*et(j,i)
626  runoff(j,i) = melt(j,i) + rainfall(j,i)-pmax*snowfall(j,i)
627 
628  end if
629 
630 #endif
631 
632 #elif (ABLSURFACE==3)
633 
634  temp_jja_help = one_third*(temp_mm(j,i,6)+temp_mm(j,i,7)+temp_mm(j,i,8))
635 
636  melt_star(j,i) = 0.0_dp ! no superimposed ice considered
637  melt(j,i) = lambda_lti*max((temp_jja_help-temp_lti), 0.0_dp)
638  runoff(j,i) = melt(j,i) + rainfall(j,i)
639 
640 #endif
641 
642 ! ---- Evaporation
643 
644  evap(j,i) = 0.0_dp
645 
646 end do
647 end do
648 
649 ! ------ Accumulation minus runoff
650 
651 as_perp = accum - evap - runoff
652 
653  !hack for initial present-day conditions:
654  !set as_perp to zero, hold bedrock fixed and let zs evolve freely (initial zs smoothing)
655  !as_perp = 0.0_dp
656 
657 ! ------ Ice-surface temperature (10-m firn temperature) temp_s,
658 ! including empirical firn-warming correction due to
659 ! refreezing meltwater when superimposed ice is formed
660 
661 where (melt_star >= melt)
662  temp_s = temp_ma + mu*(melt_star-melt)
663 elsewhere
664  temp_s = temp_ma
665 end where
666 
667 where (temp_s > -0.001_dp) temp_s = -0.001_dp
668  ! Cut-off of positive air temperatures
669 
670 !-------- Calving rate of grounded ice --------
671 
672 calv_grounded = 0.0_dp
673 
674 #if ((MARGIN==2) \
675  && (marine_ice_formation==2) \
676  && (marine_ice_calving==9))
677 
678 call calving_underwater_ice(z_sl)
680 
681 #endif
682 
683 if (firstcall) firstcall = .false.
684 
685 end subroutine boundary
686 
687 !-------------------------------------------------------------------------------
688 
689 end module boundary_m
690 !
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
real(dp), dimension(0:jmax, 0:imax, 12) temp_mm_present
temp_mm_present(j,i,n): Present-day mean monthly surface temperature
real(dp) rho_w
RHO_W: Density of pure water.
integer(i4b) specmap_time_max
specmap_time_max: Maximum time of the data values for the sea level
real(dp) mu_0
MU_0: Firn-warming correction.
real(dp), dimension(0:jmax, 0:imax) calv_grounded
calv_grounded(j,i): Calving rate of grounded ice
real(dp), dimension(0:jmax, 0:imax) snowfall
snowfall(j,i): Snowfall rate at the ice surface
Definition: sico_vars_m.F90:53
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 ...
integer(i4b) specmap_time_min
specmap_time_min: Minimum time of the data values for the sea level
Computation of the positive degree days (PDD) with statistical temperature fluctuations; based on sem...
Definition: pdd_m.F90:37
real(dp), dimension(0:jmax, 0:imax) as_perp
as_perp(j,i): Accumulation-ablation function at the ice surface (SMB)
real(dp), dimension(:), allocatable glacial_index
glacial_index(n): Data values for the glacial index
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) runoff
runoff(j,i): Runoff rate at the ice surface
real(dp), dimension(0:jmax, 0:imax) zs_ref
zs_ref(j,i): Reference elevation for precip_present, temp_ma_present and temp_mj_present ...
integer(i2b) 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(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_gi
ndata_gi: Number of data values for the glacial index
integer(i4b) ndata_grip
ndata_grip: Number of data values for the surface temperature anomaly
real(dp) beta1_0
BETA1_0: Degree-day factor for snow.
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) pmax_0
PMAX_0: Saturation factor for the formation of superimposed ice.
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(i4b) gi_time_min
gi_time_min: Minimum time of the data values for the glacial index
real(dp), dimension(:), allocatable specmap_zsl
specmap_zsl(n): Data values for the sea level
real(dp) beta2_0
BETA2_0: Degree-day factor for ice.
integer(i4b) gi_time_stp
gi_time_stp: Time step of the data values for the glacial index
real(dp), dimension(0:jmax, 0:imax) rainfall
rainfall(j,i): Rainfall rate at the ice surface
Definition: sico_vars_m.F90:55
integer(i4b) gi_time_max
gi_time_max: Maximum time of the data values for the glacial index
real(dp), dimension(0:jmax, 0:imax, 12) gamma_precip_lgm_anom
gamma_precip_lgm_anom(j,i,n): negative natural logarithm of precip_lgm_anom(j,i,n) ...
integer(i2b), dimension(0:jmax, 0:imax) maske_neu
maske_neu(j,i): New value of maske computed during an integration step
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:56
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 ...
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
subroutine, public pdd(temp_mm, s_stat, ET)
Main subroutine of pdd_m: Computation of the positive degree days (PDD) with statistical temperature ...
Definition: pdd_m.F90:57
real(dp) rho
RHO: Density of ice.
Update of the ice-land-ocean mask due to changes of the sea level.
integer, parameter dp
Double-precision reals.
real(dp), dimension(0:jmax, 0:imax, 12) precip_present
precip_present(j,i,n): Present-day mean monthly precipitation rate at the ice surface ...
real(dp), dimension(0:jmax, 0:imax, 12) precip_lgm_anom
precip_lgm_anom(j,i,n): LGM anomaly (ratio LGM/present) of the mean monthly precipitation rate at the...
real(dp), dimension(0:jmax, 0:imax, 12) temp_mm_lgm_anom
temp_mm_lgm_anom(j,i,n): LGM anomaly (difference LGM - present) of the mean monthly surface temperatu...
Declarations of global variables for SICOPOLIS.
real(dp), dimension(:), allocatable griptemp
griptemp(n): Data values for the surface temperature anomaly