SICOPOLIS V5-dev  Revision 1348
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  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), dimension(0:JMAX,0:IMAX) :: ela_dist,asp_dist,grad_dist
85 real(dp), dimension(0:JMAX,0:IMAX) :: dist,dist2
86 logical, dimension(0:JMAX,0:IMAX) :: check_point
87 logical, save :: firstcall = .true.
88 
89 !-------- Initialization of intent(out) variables --------
90 
91 z_sl_old = z_sl
92 
93 delta_ts = 0.0_dp
94 glac_index = 0.0_dp
95 z_sl = 0.0_dp
96 dzsl_dtau = 0.0_dp
97 z_mar = 0.0_dp
98 
99 !-------- Surface-temperature deviation from present values --------
100 
101 #if (TSURFACE==1)
102 delta_ts = delta_ts0
103 ! ! Steady state with prescribed constant
104 ! ! air-temperature deviation
105 #elif (TSURFACE==3)
106 delta_ts = sine_amplit &
107  *cos(2.0_dp*pi*time/(sine_period*year_sec)) &
108  -sine_amplit
109 ! ! Sinusoidal air-temperature forcing
110 #elif (TSURFACE==4)
111 
112 ! ------ delta_ts from the GRIP record
113 
114 if (time/year_sec.lt.real(grip_time_min,dp)) then
115  delta_ts = griptemp(0)
116 else if (time/year_sec.lt.real(grip_time_max,dp)) then
117 
118  i_kl = floor(((time/year_sec) &
119  -real(grip_time_min,dp))/real(grip_time_stp,dp))
120  i_kl = max(i_kl, 0)
121 
122  i_gr = ceiling(((time/year_sec) &
123  -real(grip_time_min,dp))/real(grip_time_stp,dp))
124  i_gr = min(i_gr, ndata_grip)
125 
126  if (i_kl.eq.i_gr) then
127 
128  delta_ts = griptemp(i_kl)
129 
130  else
131 
132  time_kl = (grip_time_min + i_kl*grip_time_stp) *year_sec
133  time_gr = (grip_time_min + i_gr*grip_time_stp) *year_sec
134 
135  delta_ts = griptemp(i_kl) &
136  +(griptemp(i_gr)-griptemp(i_kl)) &
137  *(time-time_kl)/(time_gr-time_kl)
138  ! linear interpolation of the ice-core data
139 
140  end if
141 
142 else
143  delta_ts = griptemp(ndata_grip)
144 end if
145 
146 delta_ts = delta_ts * grip_temp_fact
147 ! ! modification by constant factor
148 
149 #elif (TSURFACE==5)
150 ! Enter here delta_ts scenario:
151 
152 #endif
153 
154 !-------- Sea level --------
155 
156 #if (SEA_LEVEL==1)
157 ! ------ constant sea level
158 z_sl = z_sl0
159 
160 #elif (SEA_LEVEL==2)
161 ! ------ saw-tooth-shaped palaeoclimatic sea-level forcing
162 
163 z_sl_min = -130.0_dp
164 
165 t1 = -250000.0_dp *year_sec
166 t2 = -140000.0_dp *year_sec
167 t3 = -125000.0_dp *year_sec
168 t4 = -21000.0_dp *year_sec
169 t5 = -8000.0_dp *year_sec
170 t6 = 0.0_dp *year_sec
171 
172 if (time.lt.t1) then
173  z_sl = 0.0_dp
174 else if (time.lt.t2) then
175  z_sl = z_sl_min*(time-t1)/(t2-t1)
176 else if (time.lt.t3) then
177  z_sl = -z_sl_min*(time-t3)/(t3-t2)
178 else if (time.lt.t4) then
179  z_sl = z_sl_min*(time-t3)/(t4-t3)
180 else if (time.lt.t5) then
181  z_sl = -z_sl_min*(time-t5)/(t5-t4)
182 else if (time.lt.t6) then
183  z_sl = 0.0_dp
184 else
185  z_sl = 0.0_dp
186 end if
187 
188 #elif (SEA_LEVEL==3)
189 
190 ! ------ z_sl from the SPECMAP record
191 
192 if (time/year_sec.lt.real(specmap_time_min,dp)) then
193  z_sl = specmap_zsl(0)
194 else if (time/year_sec.lt.real(specmap_time_max,dp)) then
195 
196  i_kl = floor(((time/year_sec) &
197  -real(specmap_time_min,dp))/real(specmap_time_stp,dp))
198  i_kl = max(i_kl, 0)
199 
200  i_gr = ceiling(((time/year_sec) &
201  -real(specmap_time_min,dp))/real(specmap_time_stp,dp))
202  i_gr = min(i_gr, ndata_specmap)
203 
204  if (i_kl.eq.i_gr) then
205 
206  z_sl = specmap_zsl(i_kl)
207 
208  else
209 
210  time_kl = (specmap_time_min + i_kl*specmap_time_stp) *year_sec
211  time_gr = (specmap_time_min + i_gr*specmap_time_stp) *year_sec
212 
213  z_sl = specmap_zsl(i_kl) &
214  +(specmap_zsl(i_gr)-specmap_zsl(i_kl)) &
215  *(time-time_kl)/(time_gr-time_kl)
216  ! linear interpolation of the sea-level data
217 
218  end if
219 
220 else
221  z_sl = specmap_zsl(ndata_specmap)
222 end if
223 
224 #endif
225 
226 ! ------ Time derivative of the sea level
227 
228 if ( z_sl_old > -999999.9_dp ) then
229  dzsl_dtau = (z_sl-z_sl_old)/dtime
230 else ! only dummy value for z_sl_old available
231  dzsl_dtau = 0.0_dp
232 end if
233 
234 ! ------ Minimum bedrock elevation for extent of marine ice
235 
236 #if (MARGIN==2)
237 
238 #if (MARINE_ICE_CALVING==2 || MARINE_ICE_CALVING==3)
239 z_mar = z_mar
240 #elif (MARINE_ICE_CALVING==4 || MARINE_ICE_CALVING==5)
241 z_mar = fact_z_mar*z_sl
242 #elif (MARINE_ICE_CALVING==6 || MARINE_ICE_CALVING==7)
243 if (z_sl >= -80.0_dp) then
244  z_mar = 2.5_dp*z_sl
245 else
246  z_mar = 10.25_dp*(z_sl+80.0_dp)-200.0_dp
247 end if
248 z_mar = fact_z_mar*z_mar
249 #endif
250 
251 #endif
252 
253 ! ------ Update of the mask according to the sea level
254 
255 ! ---- Check all sea and floating-ice points and their direct
256 ! neighbours
257 
258 do i=0, imax
259 do j=0, jmax
260  check_point(j,i) = .false.
261 end do
262 end do
263 
264 do i=1, imax-1
265 do j=1, jmax-1
266  if (maske(j,i).ge.2_i1b) then
267  check_point(j ,i ) = .true.
268  check_point(j ,i+1) = .true.
269  check_point(j ,i-1) = .true.
270  check_point(j+1,i ) = .true.
271  check_point(j-1,i ) = .true.
272  end if
273 end do
274 end do
275 
276 do i=1, imax-1
277 do j=1, jmax-1
278  if (check_point(j,i)) then
279  maske_neu(j,i) = mask_update_sea_level(z_sl, i, j)
280  end if
281 end do
282 end do
283 
284 ! ---- Assign new values of the mask
285 
286 do i=1, imax-1
287 do j=1, jmax-1
288  if (check_point(j,i)) then
289  maske(j,i) = maske_neu(j,i)
290  end if
291 end do
292 end do
293 
294 !-------- Surface air temperature and accumulation-ablation function --------
295 
296 #if (!defined(SURFACE_FORCING) || SURFACE_FORCING==1)
297 
298 errormsg = ' >>> boundary: Option SURFACE_FORCING=1 not supported any more!'
299 call error(errormsg)
300 
301 #elif (SURFACE_FORCING==2)
302 
303 do i=0, imax
304 do j=0, jmax
305 
306  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
307  temp_s(j,i) = temp_s(j,i) + delta_ts
308  ! Correction with temperature deviation delta_ts
309  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
310  ! Cut-off of positive air temperatures
311 
312  accum_present(j,i) = s_0
313  accum(j,i) = accum_present(j,i)
314 
315  as_perp(j,i) = m_0*(zs(j,i)-ela)
316  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
317 
318 end do
319 end do
320 
321 #elif (SURFACE_FORCING==3)
322 
323 do i=0, imax
324 do j=0, jmax
325 
326  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
327  temp_s(j,i) = temp_s(j,i) + delta_ts
328  ! Correction with temperature deviation delta_ts
329  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
330  ! Cut-off of positive air temperatures
331 
332  accum_present(j,i) = s_0
333  accum(j,i) = accum_present(j,i)
334 
335  if (xi(i) >= x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan((eta(j)-y_gip)/(xi(i)-x_gip)) + phi_0)
336  if (xi(i) < x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan(-(eta(j)-y_gip)/(xi(i)-x_gip)) - phi_0)
337  if (xi(i) == x_gip .AND. eta(j) == y_gip) ela_dist(j,i) = ela
338 
339  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
340  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = accum(j,i)
341 
342 end do
343 end do
344 
345 #elif (SURFACE_FORCING==4)
346 
347 do i=0, imax
348 do j=0, jmax
349 
350  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
351  temp_s(j,i) = temp_s(j,i) + delta_ts
352  ! Correction with temperature deviation delta_ts
353  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
354  ! Cut-off of positive air temperatures
355 
356  if (xi(i) >= x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan((eta(j)-y_gip)/(xi(i)-x_gip)) + phi_0)
357  if (xi(i) < x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan(-(eta(j)-y_gip)/(xi(i)-x_gip)) - phi_0)
358  if (xi(i) == x_gip .AND. eta(j) == y_gip) ela_dist(j,i) = ela
359 
360  if (zs(j,i) <= z_gc) as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
361  if (zs(j,i) > z_gc) as_perp(j,i) = m_1*(zs(j,i)-ela_dist(j,i))+m_0*(z_gc-ela_dist(j,i))
362 
363 end do
364 end do
365 
366 !open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/SMB/smb_sico.dat')
367 !do j = 0, JMAX
368 ! write(13,*) as_perp(j,:)
369 !end do
370 !close(13)
371 
372 #elif (SURFACE_FORCING==5)
373 
374 asp_a = 0
375 asp_b = 0
376 
377 do i=0, imax
378 do j=0, jmax
379 
380  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
381  temp_s(j,i) = temp_s(j,i) + delta_ts
382  ! Correction with temperature deviation delta_ts
383  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
384  ! Cut-off of positive air temperatures
385 
386  asp_a = zs(j,i+1) - zs(j,i-1)
387  asp_b = zs(j-1,i) - zs(j+1,i)
388  if (asp_b < 0) asp_dist(j,i) = mod(atan(asp_a/asp_b),6.2832)
389  if (asp_b >= 0) asp_dist(j,i) = atan(asp_a/asp_b) + 3.1416
390 
391  ela_dist(j,i) = ela + ela_amp * sin(asp_dist(j,i) - phi_0 + 3.1416)
392 
393 ! if (xi(i) >= x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan((eta(j)-y_gip)/(xi(i)-x_gip)) + phi_0)
394 ! if (xi(i) < x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan(-(eta(j)-y_gip)/(xi(i)-x_gip)) - phi_0)
395 ! if (xi(i) == x_gip .AND. eta(j) == y_gip) ela_dist(j,i) = ela
396 
397  if (zs(j,i) <= z_gc) as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
398  if (zs(j,i) > z_gc) as_perp(j,i) = m_1*(zs(j,i)-ela_dist(j,i))+m_0*(z_gc-ela_dist(j,i))
399 
400 end do
401 end do
402 
403 !open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/Aspect/asp_fortran.dat')
404 !do j = 0, JMAX
405 ! write(13,*) asp_dist(j,:)
406 !end do
407 !close(13)
408 
409 !open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/Aspect/zs_fortran.dat')
410 !do j = 0, JMAX
411 ! write(13,*) zs(j,:)
412 !end do
413 !close(13)
414 
415 !open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/Aspect/ela_fortran.dat')
416 !do j = 0, JMAX
417 ! write(13,*) ela_dist(j,:)
418 !end do
419 !close(13)
420 
421 #elif (SURFACE_FORCING==6)
422 
423 asp_a = 0
424 asp_b = 0
425 
426 do i=0, imax
427 do j=0, jmax
428 
429  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
430  temp_s(j,i) = temp_s(j,i) + delta_ts
431  ! Correction with temperature deviation delta_ts
432  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
433  ! Cut-off of positive air temperatures
434 
435  asp_a = zs(j,i+1) - zs(j,i-1)
436  asp_b = zs(j-1,i) - zs(j+1,i)
437  if (asp_b < 0) asp_dist(j,i) = mod(atan(asp_a/asp_b),6.2832)
438  if (asp_b >= 0) asp_dist(j,i) = atan(asp_a/asp_b) + 3.1416
439 
440  grad_dist(j,i) = sqrt((asp_a/(2*dxi))**2 + (asp_b/(2*deta))**2)*45
441 
442  ela_dist(j,i) = ela + ela_amp * sin(asp_dist(j,i) - phi_0 + 3.1416)
443  if (grad_dist(j,i) > tgt) ela_dist(j,i) = ela_dist(j,i) + grad_dist(j,i)*2
444 
445  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
446  if (as_perp(j,i) > s_0) as_perp(j,i) = s_0
447 
448 end do
449 end do
450 
451 open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/TopGrad/smb_fortran.dat')
452 do j = 0, jmax
453  write(13,*) as_perp(j,:)
454 end do
455 close(13)
456 
457 !open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/TopGrad/vars.dat')
458 ! write(13,*) dxi,deta
459 !close(13)
460 
461 #elif (SURFACE_FORCING==7)
462 
463 asp_a = 0
464 asp_b = 0
465 
466 do i=0, imax
467 do j=0, jmax
468 
469  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
470  temp_s(j,i) = temp_s(j,i) + delta_ts
471  ! Correction with temperature deviation delta_ts
472  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
473  ! Cut-off of positive air temperatures
474 
475  asp_a = zs(j,i+1) - zs(j,i-1)
476  asp_b = zs(j-1,i) - zs(j+1,i)
477  if (asp_b < 0) asp_dist(j,i) = mod(atan(asp_a/asp_b),6.2832)
478  if (asp_b >= 0) asp_dist(j,i) = atan(asp_a/asp_b) + 3.1416
479 
480  ela_dist(j,i) = ela + ela_amp * sin(asp_dist(j,i) - phi_0 + 3.1416)
481 
482  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
483  if (as_perp(j,i) > s_0) as_perp(j,i) = s_0
484 
485 end do
486 end do
487 
488 open(unit = 13, file = '/home/matthiasscheiter/Dropbox/sicopolis_erg/PureAspect/ela_fortran.dat')
489 do j = 0, jmax
490  write(13,*) ela_dist(j,:)
491 end do
492 close(13)
493 
494 #elif (SURFACE_FORCING==8)
495 
496 do i=0, imax
497 do j=0, jmax
498  dist(j,i) = sqrt( (xi(i)-x_gip)**2 + (eta(j)-y_gip)**2 )
499  dist2(j,i) = sqrt( (xi(i)-x_gip2)**2 + (eta(j)-y_gip2)**2 )
500 end do
501 end do
502 
503 
504 do i=0, imax
505 do j=0, jmax
506 
507  temp_s(j,i) = temp_0 - gamma_t*zs(j,i)
508  temp_s(j,i) = temp_s(j,i) + delta_ts
509  ! Correction with temperature deviation delta_ts
510  if (temp_s(j,i) > -0.001_dp) temp_s(j,i) = -0.001_dp
511  ! Cut-off of positive air temperatures
512 
513  if (dist(j,i) .lt. 3000.) then
514  if (xi(i) >= x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan((eta(j)-y_gip)/(xi(i)-x_gip)) + phi_0)
515  if (xi(i) < x_gip) ela_dist(j,i) = ela + ela_amp * sin(atan(-(eta(j)-y_gip)/(xi(i)-x_gip)) - phi_0)
516  if (xi(i) == x_gip .AND. eta(j) == y_gip) ela_dist(j,i) = ela
517  else
518  ela_dist(j,i) = ela
519  end if
520 
521  if (dist2(j,i) .lt. 1500.) then
522  if (xi(i) >= x_gip2) ela_dist(j,i) = ela_dist(j,i) + ela_amp2 * sin(atan((eta(j)-y_gip2)/(xi(i)-x_gip2)) + phi_0)
523  if (xi(i) < x_gip2) ela_dist(j,i) = ela_dist(j,i) + ela_amp2 * sin(atan(-(eta(j)-y_gip2)/(xi(i)-x_gip2)) - phi_0)
524  if (xi(i) == x_gip2 .AND. eta(j) == y_gip2) ela_dist(j,i) = ela_dist(j,i)
525  end if
526 
527  as_perp(j,i) = m_0*(zs(j,i)-ela_dist(j,i))
528  if (as_perp(j,i) > accum(j,i)) as_perp(j,i) = s_0
529 
530 end do
531 end do
532 
533 #endif
534 
535 runoff = accum - as_perp
536 evap = 0.0_dp
537 
538 !-------- Calving rate of grounded ice --------
539 
540 calv_grounded = 0.0_dp
541 
542 #if ((MARGIN==2) \
543  && (marine_ice_formation==2) \
544  && (marine_ice_calving==9))
545 
546 call calving_underwater_ice(z_sl)
548 
549 #endif
550 
551 if (firstcall) firstcall = .false.
552 
553 end subroutine boundary
554 
555 !-------------------------------------------------------------------------------
556 
557 end module boundary_m
558 !
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) 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
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