SICOPOLIS V3.2
 All Classes Files Functions Variables Macros
calc_temp_enth.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Subroutine : c a l c _ t e m p _ e n t h
4 !
5 !> @file
6 !!
7 !! Computation of temperature, water content and age with the enthalpy method.
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2013-2016 Ralf Greve, Heinz Blatter
12 !!
13 !! @section License
14 !!
15 !! This file is part of SICOPOLIS.
16 !!
17 !! SICOPOLIS is free software: you can redistribute it and/or modify
18 !! it under the terms of the GNU General Public License as published by
19 !! the Free Software Foundation, either version 3 of the License, or
20 !! (at your option) any later version.
21 !!
22 !! SICOPOLIS is distributed in the hope that it will be useful,
23 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
24 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 !! GNU General Public License for more details.
26 !!
27 !! You should have received a copy of the GNU General Public License
28 !! along with SICOPOLIS. If not, see <http://www.gnu.org/licenses/>.
29 !<
30 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
31 
32 !-------------------------------------------------------------------------------
33 !> Computation of temperature, water content and age with the enthalpy method.
34 !<------------------------------------------------------------------------------
35 subroutine calc_temp_enth(dxi, deta, dzeta_c, dzeta_t, dzeta_r, &
36  dtime_temp)
37 
38 use sico_types
40 use sico_vars
41 
42 implicit none
43 
44 real(dp), intent(in) :: dxi, deta, dzeta_c, dzeta_t, dzeta_r
45 real(dp), intent(in) :: dtime_temp
46 
47 integer(i4b) :: i, j, kc, kt, kr, ii, jj
48 real(dp) :: at1(0:kcmax), at2_1(0:kcmax), at2_2(0:kcmax), &
49  at3_1(0:kcmax), at3_2(0:kcmax), at4_1(0:kcmax), &
50  at4_2(0:kcmax), at5(0:kcmax), at6(0:kcmax), at7, &
51  ai1(0:kcmax), ai2(0:kcmax), &
52  atr1, acb1, acb2, acb3, acb4, alb1, aqtlde(0:kcmax), &
53  am1, am3(0:kcmax)
54 real(dp) :: dtt_2dxi, dtt_2deta
55 
56 !-------- Term abbreviations
57 
58 at7 = 2.0_dp/rho*dtime_temp
59 
60 atr1 = kappa_r/(rho_c_r*h_r**2)*dtime_temp/(dzeta_r**2)
61 
62 if (flag_aa_nonzero) then
63  am1 = aa*beta*dzeta_c/(ea-1.0_dp)
64 else
65  am1 = beta*dzeta_c
66 end if
67 
68 if (flag_aa_nonzero) then
69  acb1 = (ea-1.0_dp)/aa/dzeta_c
70 else
71  acb1 = 1.0_dp/dzeta_c
72 end if
73 
74 acb2 = kappa_r/h_r/dzeta_r
75 acb3 = rho*g
76 acb4 = rho*g
77 
78 alb1 = h_r/kappa_r*dzeta_r
79 
80 dtt_2dxi = 0.5_dp*dtime_temp/dxi
81 dtt_2deta = 0.5_dp*dtime_temp/deta
82 
83 do kc=0, kcmax
84 
85  if (flag_aa_nonzero) then
86 
87  at1(kc) = (ea-1.0_dp)/(aa*eaz_c(kc))*dtime_temp/dzeta_c
88  at2_1(kc) = (ea-1.0_dp)/(aa*eaz_c(kc))*dtime_temp/dzeta_c
89  at2_2(kc) = (eaz_c(kc)-1.0_dp)/(aa*eaz_c(kc)) &
90  *dtime_temp/dzeta_c
91  at3_1(kc) = (ea-1.0_dp)/(aa*eaz_c(kc))*dtime_temp/dzeta_c
92  at3_2(kc) = (eaz_c(kc)-1.0_dp)/(aa*eaz_c(kc)) &
93  *dtime_temp/dzeta_c
94  at4_1(kc) = (ea-1.0_dp)/(aa*eaz_c(kc))*dtime_temp/dzeta_c
95  at4_2(kc) = (eaz_c(kc)-1.0_dp)/(aa*eaz_c(kc)) &
96  *dtime_temp/dzeta_c
97  at5(kc) = (ea-1.0_dp)/(rho*aa*eaz_c(kc)) &
98  *dtime_temp/dzeta_c
99  if (kc /= kcmax) then
100  at6(kc) = (ea-1.0_dp) &
101  /(aa*exp(aa*0.5_dp*(zeta_c(kc)+zeta_c(kc+1)))) &
102  /dzeta_c
103  else
104  at6(kc) = 0.0_dp
105  end if
106  ai1(kc) = agediff*(ea-1.0_dp)/(aa*eaz_c(kc)) &
107  *dtime_temp/dzeta_c
108  if (kc /= kcmax) then
109  ai2(kc) = (ea-1.0_dp) &
110  /(aa*exp(aa*0.5_dp*(zeta_c(kc)+zeta_c(kc+1)))) &
111  /dzeta_c
112  else
113  ai2(kc) = 0.0_dp
114  end if
115  aqtlde(kc) = (aa*eaz_c(kc))/(ea-1.0_dp)*dzeta_c/dtime_temp
116  am3(kc) = (aa*eaz_c(kc))/(ea-1.0_dp)*dzeta_c*beta
117 
118  else
119 
120  at1(kc) = dtime_temp/dzeta_c
121  at2_1(kc) = dtime_temp/dzeta_c
122  at2_2(kc) = zeta_c(kc) &
123  *dtime_temp/dzeta_c
124  at3_1(kc) = dtime_temp/dzeta_c
125  at3_2(kc) = zeta_c(kc) &
126  *dtime_temp/dzeta_c
127  at4_1(kc) = dtime_temp/dzeta_c
128  at4_2(kc) = zeta_c(kc) &
129  *dtime_temp/dzeta_c
130  at5(kc) = 1.0_dp/rho &
131  *dtime_temp/dzeta_c
132  if (kc /= kcmax) then
133  at6(kc) = 1.0_dp &
134  /dzeta_c
135  else
136  at6(kc) = 0.0_dp
137  end if
138  ai1(kc) = agediff &
139  *dtime_temp/dzeta_c
140  if (kc /= kcmax) then
141  ai2(kc) = 1.0_dp &
142  /dzeta_c
143  else
144  ai2(kc) = 0.0_dp
145  end if
146  aqtlde(kc) = dzeta_c/dtime_temp
147  am3(kc) = dzeta_c*beta
148 
149  end if
150 
151 end do
152 
153 !-------- Computation loop --------
154 
155 do i=1, imax-1 ! skipping domain margins
156 do j=1, jmax-1 ! skipping domain margins
157 
158  if (maske(j,i)==0) then ! glaciated land
159 
160 ! ------ Old vertical column cold
161 
162  if (n_cts(j,i) == -1) then
163 
164  n_cts_neu(j,i) = -1
165  kc_cts_neu(j,i) = 0
166  zm_neu(j,i) = zb(j,i)
167  h_c_neu(j,i) = h_c(j,i)
168  h_t_neu(j,i) = 0.0_dp
169 
170  call calc_temp_enth_1(at1, at2_1, at2_2, at3_1, at3_2, &
171  at4_1, at4_2, at5, at6, at7, &
172  atr1, acb1, acb2, acb3, acb4, alb1, &
173  ai1, ai2, &
174  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
175 
176 ! ---- Check whether base has become temperate
177 
178  if (temp_c_neu(0,j,i) > temp_c_m(0,j,i)-eps) then
179 
180  n_cts_neu(j,i) = 0
181  kc_cts_neu(j,i) = 0
182 
183  call calc_temp_enth_2(at1, at2_1, at2_2, at3_1, at3_2, &
184  at4_1, at4_2, at5, at6, at7, atr1, alb1, &
185  ai1, ai2, aqtlde, am3, &
186  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
187 
188  end if
189 
190 ! ------ Old vertical column with temperate base
191 
192  else if (n_cts(j,i) == 0) then
193 
194  n_cts_neu(j,i) = 0
195  kc_cts_neu(j,i) = kc_cts(j,i)
196  zm_neu(j,i) = zb(j,i)
197  h_c_neu(j,i) = h_c(j,i)
198  h_t_neu(j,i) = h_t(j,i)
199 
200  call calc_temp_enth_2(at1, at2_1, at2_2, at3_1, at3_2, &
201  at4_1, at4_2, at5, at6, at7, atr1, alb1, &
202  ai1, ai2, aqtlde, am3, &
203  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
204 
205 ! ---- Check whether temperate base becomes cold
206 
207  if ( (temp_c_neu(1,j,i)-temp_c_neu(0,j,i)) < (am1*h_c(j,i)) ) then
208 
209  n_cts_neu(j,i) = -1
210  kc_cts_neu(j,i) = 0
211 
212  call calc_temp_enth_1(at1, at2_1, at2_2, at3_1, at3_2, &
213  at4_1, at4_2, at5, at6, at7, &
214  atr1, acb1, acb2, acb3, acb4, alb1, &
215  ai1, ai2, &
216  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
217 
218  if (temp_c_neu(0,j,i) > temp_c_m(0,j,i)-eps) then
219 
220  n_cts_neu(j,i) = 0
221  kc_cts_neu(j,i) = 0
222 
223  call calc_temp_enth_2(at1, at2_1, at2_2, at3_1, at3_2, &
224  at4_1, at4_2, at5, at6, at7, atr1, alb1, &
225  ai1, ai2, aqtlde, am3, &
226  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
227 
228  end if
229 
230  end if
231 
232  end if
233 
234 #if (MARGIN==3)
235 
236  else if (maske(j,i)==3) then ! floating ice
237 
238  n_cts_neu(j,i) = -1
239  kc_cts_neu(j,i) = 0
240  zm_neu(j,i) = zb(j,i)
241  h_c_neu(j,i) = h_c(j,i)
242  h_t_neu(j,i) = 0.0_dp
243 
244  call calc_temp_enth_ssa(at1, at2_1, at2_2, at3_1, at3_2, &
245  at4_1, at4_2, at5, at6, at7, atr1, alb1, &
246  ai1, ai2, &
247  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
248 
249 ! ------ Reset temperatures above melting to the melting point
250 ! and water contents above zero to zero
251 ! (should not occur, but just in case)
252 
253  do kc=0, kcmax
254  if (temp_c_neu(kc,j,i) > temp_c_m(kc,j,i)) &
255  temp_c_neu(kc,j,i) = temp_c_m(kc,j,i)
256  if (omega_c_neu(kc,j,i) > 0.0_dp) &
257  omega_c_neu(kc,j,i) = 0.0_dp
258  end do
259 
260 #endif
261 
262  else ! maske(j,i) == 1,2 (ice-free land or sea point)
263 
264  n_cts_neu(j,i) = -1
265  kc_cts_neu(j,i) = 0
266  zm_neu(j,i) = zb(j,i)
267  h_c_neu(j,i) = h_c(j,i)
268  h_t_neu(j,i) = 0.0_dp
269 
270  call calc_temp_r(atr1, alb1, i, j)
271 
272  end if
273 
274 end do
275 end do
276 
277 !-------- Extrapolate values on margins --------
278 
279 ! ------ Lower left corner
280 
281 i=0
282 j=0
283 
284 if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
285  ! glaciated land or floating ice
286  ii=i+1
287  jj=j+1
288 
289  do kc=0, kcmax
290  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
291  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
292  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
293  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
294  end do
295 
296  do kt=0, ktmax
297  ! redundant, lower (kt) ice layer
298  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
299  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
300  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
301  end do
302 
303  do kr=0, krmax
304  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
305  end do
306 
307  n_cts_neu(j,i) = n_cts_neu(jj,ii)
308  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
309  zm_neu(j,i) = zb(j,i)
310  h_c_neu(j,i) = h_c(j,i)
311  h_t_neu(j,i) = h_t(j,i)
312 
313 else ! maske(j,i) == 1,2 (ice-free land or sea point)
314 
315  n_cts_neu(j,i) = -1
316  kc_cts_neu(j,i) = 0
317  zm_neu(j,i) = zb(j,i)
318  h_c_neu(j,i) = h_c(j,i)
319  h_t_neu(j,i) = h_t(j,i)
320 
321  call calc_temp_r(atr1, alb1, i, j)
322 
323 end if
324 
325 ! ------ Lower right corner
326 
327 i=imax
328 j=0
329 
330 if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
331  ! glaciated land or floating ice
332  ii=i-1
333  jj=j+1
334 
335  do kc=0, kcmax
336  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
337  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
338  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
339  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
340  end do
341 
342  do kt=0, ktmax
343  ! redundant, lower (kt) ice layer
344  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
345  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
346  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
347  end do
348 
349  do kr=0, krmax
350  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
351  end do
352 
353  n_cts_neu(j,i) = n_cts_neu(jj,ii)
354  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
355  zm_neu(j,i) = zb(j,i)
356  h_c_neu(j,i) = h_c(j,i)
357  h_t_neu(j,i) = h_t(j,i)
358 
359 else ! maske(j,i) == 1,2 (ice-free land or sea point)
360 
361  n_cts_neu(j,i) = -1
362  kc_cts_neu(j,i) = 0
363  zm_neu(j,i) = zb(j,i)
364  h_c_neu(j,i) = h_c(j,i)
365  h_t_neu(j,i) = h_t(j,i)
366 
367  call calc_temp_r(atr1, alb1, i, j)
368 
369 end if
370 
371 ! ------ Upper left corner
372 
373 i=0
374 j=jmax
375 
376 if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
377  ! glaciated land or floating ice
378  ii=i+1
379  jj=j-1
380 
381  do kc=0, kcmax
382  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
383  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
384  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
385  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
386  end do
387 
388  do kt=0, ktmax
389  ! redundant, lower (kt) ice layer
390  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
391  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
392  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
393  end do
394 
395  do kr=0, krmax
396  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
397  end do
398 
399  n_cts_neu(j,i) = n_cts_neu(jj,ii)
400  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
401  zm_neu(j,i) = zb(j,i)
402  h_c_neu(j,i) = h_c(j,i)
403  h_t_neu(j,i) = h_t(j,i)
404 
405 else ! maske(j,i) == 1,2 (ice-free land or sea point)
406 
407  n_cts_neu(j,i) = -1
408  kc_cts_neu(j,i) = 0
409  zm_neu(j,i) = zb(j,i)
410  h_c_neu(j,i) = h_c(j,i)
411  h_t_neu(j,i) = h_t(j,i)
412 
413  call calc_temp_r(atr1, alb1, i, j)
414 
415 end if
416 
417 ! ------ Upper right corner
418 
419 i=imax
420 j=jmax
421 
422 if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
423  ! glaciated land or floating ice
424  ii=i-1
425  jj=j-1
426 
427  do kc=0, kcmax
428  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
429  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
430  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
431  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
432  end do
433 
434  do kt=0, ktmax
435  ! redundant, lower (kt) ice layer
436  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
437  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
438  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
439  end do
440 
441  do kr=0, krmax
442  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
443  end do
444 
445  n_cts_neu(j,i) = n_cts_neu(jj,ii)
446  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
447  zm_neu(j,i) = zb(j,i)
448  h_c_neu(j,i) = h_c(j,i)
449  h_t_neu(j,i) = h_t(j,i)
450 
451 else ! maske(j,i) == 1,2 (ice-free land or sea point)
452 
453  n_cts_neu(j,i) = -1
454  kc_cts_neu(j,i) = 0
455  zm_neu(j,i) = zb(j,i)
456  h_c_neu(j,i) = h_c(j,i)
457  h_t_neu(j,i) = h_t(j,i)
458 
459  call calc_temp_r(atr1, alb1, i, j)
460 
461 end if
462 
463 ! ------ Lower and upper margins
464 
465 do i=1, imax-1
466 
467 ! ---- Lower margin
468 
469  j=0
470 
471  if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
472  ! glaciated land or floating ice
473  ii=i
474  jj=j+1
475 
476  do kc=0, kcmax
477  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
478  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
479  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
480  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
481  end do
482 
483  do kt=0, ktmax
484  ! redundant, lower (kt) ice layer
485  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
486  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
487  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
488  end do
489 
490  do kr=0, krmax
491  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
492  end do
493 
494  n_cts_neu(j,i) = n_cts_neu(jj,ii)
495  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
496  zm_neu(j,i) = zb(j,i)
497  h_c_neu(j,i) = h_c(j,i)
498  h_t_neu(j,i) = h_t(j,i)
499 
500  else ! maske(j,i) == 1,2 (ice-free land or sea point)
501 
502  n_cts_neu(j,i) = -1
503  kc_cts_neu(j,i) = 0
504  zm_neu(j,i) = zb(j,i)
505  h_c_neu(j,i) = h_c(j,i)
506  h_t_neu(j,i) = h_t(j,i)
507 
508  call calc_temp_r(atr1, alb1, i, j)
509 
510  end if
511 
512 ! ---- Upper margin
513 
514  j=jmax
515 
516  if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
517  ! glaciated land or floating ice
518  ii=i
519  jj=j-1
520 
521  do kc=0, kcmax
522  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
523  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
524  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
525  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
526  end do
527 
528  do kt=0, ktmax
529  ! redundant, lower (kt) ice layer
530  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
531  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
532  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
533  end do
534 
535  do kr=0, krmax
536  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
537  end do
538 
539  n_cts_neu(j,i) = n_cts_neu(jj,ii)
540  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
541  zm_neu(j,i) = zb(j,i)
542  h_c_neu(j,i) = h_c(j,i)
543  h_t_neu(j,i) = h_t(j,i)
544 
545  else ! maske(j,i) == 1,2 (ice-free land or sea point)
546 
547  n_cts_neu(j,i) = -1
548  kc_cts_neu(j,i) = 0
549  zm_neu(j,i) = zb(j,i)
550  h_c_neu(j,i) = h_c(j,i)
551  h_t_neu(j,i) = h_t(j,i)
552 
553  call calc_temp_r(atr1, alb1, i, j)
554 
555  end if
556 
557 end do
558 
559 ! ------ Left and right margins
560 
561 do j=1, jmax-1
562 
563 ! ---- Left margin
564 
565  i=0
566 
567  if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
568  ! glaciated land or floating ice
569  ii=i+1
570  jj=j
571 
572  do kc=0, kcmax
573  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
574  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
575  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
576  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
577  end do
578 
579  do kt=0, ktmax
580  ! redundant, lower (kt) ice layer
581  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
582  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
583  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
584  end do
585 
586  do kr=0, krmax
587  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
588  end do
589 
590  n_cts_neu(j,i) = n_cts_neu(jj,ii)
591  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
592  zm_neu(j,i) = zb(j,i)
593  h_c_neu(j,i) = h_c(j,i)
594  h_t_neu(j,i) = h_t(j,i)
595 
596  else ! maske(j,i) == 1,2 (ice-free land or sea point)
597 
598  n_cts_neu(j,i) = -1
599  kc_cts_neu(j,i) = 0
600  zm_neu(j,i) = zb(j,i)
601  h_c_neu(j,i) = h_c(j,i)
602  h_t_neu(j,i) = h_t(j,i)
603 
604  call calc_temp_r(atr1, alb1, i, j)
605 
606  end if
607 
608 ! ---- Right margin
609 
610  i=imax
611 
612  if ( (maske(j,i) == 0).or.(maske(j,i) == 3) ) then
613  ! glaciated land or floating ice
614  ii=i-1
615  jj=j
616 
617  do kc=0, kcmax
618  enth_c_neu(kc,j,i) = enth_c_neu(kc,jj,ii)
619  temp_c_neu(kc,j,i) = temp_c_neu(kc,jj,ii)
620  omega_c_neu(kc,j,i) = omega_c_neu(kc,jj,ii)
621  age_c_neu(kc,j,i) = age_c_neu(kc,jj,ii)
622  end do
623 
624  do kt=0, ktmax
625  ! redundant, lower (kt) ice layer
626  enth_t_neu(kt,j,i) = enth_t_neu(kt,jj,ii)
627  omega_t_neu(kt,j,i) = omega_t_neu(kt,jj,ii)
628  age_t_neu(kt,j,i) = age_t_neu(kt,jj,ii)
629  end do
630 
631  do kr=0, krmax
632  temp_r_neu(kr,j,i) = temp_r_neu(kr,jj,ii)
633  end do
634 
635  n_cts_neu(j,i) = n_cts_neu(jj,ii)
636  kc_cts_neu(j,i) = kc_cts_neu(jj,ii)
637  zm_neu(j,i) = zb(j,i)
638  h_c_neu(j,i) = h_c(j,i)
639  h_t_neu(j,i) = h_t(j,i)
640 
641  else ! maske(j,i) == 1,2 (ice-free land or sea point)
642 
643  n_cts_neu(j,i) = -1
644  kc_cts_neu(j,i) = 0
645  zm_neu(j,i) = zb(j,i)
646  h_c_neu(j,i) = h_c(j,i)
647  h_t_neu(j,i) = h_t(j,i)
648 
649  call calc_temp_r(atr1, alb1, i, j)
650 
651  end if
652 
653 end do
654 
655 end subroutine calc_temp_enth
656 !
subroutine calc_temp_enth(dxi, deta, dzeta_c, dzeta_t, dzeta_r, dtime_temp)
Computation of temperature, water content and age with the enthalpy method.
Declarations of kind types for SICOPOLIS.
Definition: sico_types.F90:35
Declarations of global variables for SICOPOLIS (for the ANT domain).
Definition: sico_vars.F90:35
subroutine calc_temp_enth_2(at1, at2_1, at2_2, at3_1, at3_2, at4_1, at4_2, at5, at6, at7, atr1, alb1, ai1, ai2, aqtlde, am3, dtime_temp, dtt_2dxi, dtt_2deta, i, j)
Computation of temperature and age for an ice column with a temperate base with the enthalpy method...
subroutine calc_temp_enth_1(at1, at2_1, at2_2, at3_1, at3_2, at4_1, at4_2, at5, at6, at7, atr1, acb1, acb2, acb3, acb4, alb1, ai1, ai2, dtime_temp, dtt_2dxi, dtt_2deta, i, j)
Computation of temperature and age for a cold ice column with the enthalpy method.
subroutine calc_temp_enth_ssa(at1, at2_1, at2_2, at3_1, at3_2, at4_1, at4_2, at5, at6, at7, atr1, alb1, ai1, ai2, dtime_temp, dtt_2dxi, dtt_2deta, i, j)
Computation of temperature and age for ice shelves (floating ice) with the enthalpy method...
Declarations of global variables for SICOPOLIS.
subroutine calc_temp_r(atr1, alb1, i, j)
Computation of temperature for an ice-free column.
Definition: calc_temp_r.F90:35