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