SICOPOLIS V3.2
 All Classes Files Functions Variables Macros
calc_temp_ssa.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Subroutine : c a l c _ t e m p _ s s a
4 !
5 !> @file
6 !!
7 !! Computation of temperature and age for ice shelves (floating ice).
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2009-2016 Ralf Greve, Tatsuru Sato
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 and age for ice shelves (floating ice).
34 !<------------------------------------------------------------------------------
35 subroutine calc_temp_ssa(at1, at2_1, at2_2, at3_1, at3_2, &
36  at4_1, at4_2, at5, at6, at7, atr1, alb1, &
37  ai1, ai2, &
38  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
39 
40 use sico_types
41 use sico_variables
42 use sico_vars
43 use ice_material_quantities, only : kappa_val, c_val
44 use sico_sle_solvers, only : tri_sle
45 
46 implicit none
47 
48 integer(i4b), intent(in) :: i, j
49 real(dp), intent(in) :: at1(0:kcmax), at2_1(0:kcmax), at2_2(0:kcmax), &
50  at3_1(0:kcmax), at3_2(0:kcmax), at4_1(0:kcmax), &
51  at4_2(0:kcmax), at5(0:kcmax), at6(0:kcmax), at7, &
52  ai1(0:kcmax), ai2(0:kcmax), &
53  atr1, alb1
54 real(dp), intent(in) :: dtime_temp, dtt_2dxi, dtt_2deta
55 
56 integer(i4b) :: kc, kt, kr
57 real(dp) :: ct1(0:kcmax), ct2(0:kcmax), ct3(0:kcmax), ct4(0:kcmax), &
58  ct5(0:kcmax), ct6(0:kcmax), ct7(0:kcmax), ctr1, clb1
59 real(dp) :: ct1_sg(0:kcmax), ct2_sg(0:kcmax), ct3_sg(0:kcmax), &
60  ct4_sg(0:kcmax), adv_vert_sg(0:kcmax), abs_adv_vert_sg(0:kcmax)
61 real(dp) :: ci1(0:kcmax), ci2(0:kcmax)
62 real(dp) :: temp_c_help(0:kcmax)
63 real(dp) :: vx_c_help, vy_c_help
64 real(dp) :: adv_vert_help
65 real(dp) :: dtt_dxi, dtt_deta
66 real(dp) :: lgs_a0(0:kcmax+ktmax+krmax+imax+jmax), &
67  lgs_a1(0:kcmax+ktmax+krmax+imax+jmax), &
68  lgs_a2(0:kcmax+ktmax+krmax+imax+jmax), &
69  lgs_x(0:kcmax+ktmax+krmax+imax+jmax), &
70  lgs_b(0:kcmax+ktmax+krmax+imax+jmax)
71 real(dp), parameter :: zero=0.0_dp
72 
73 !-------- Check for boundary points --------
74 
75 if ((i == 0).or.(i == imax).or.(j == 0).or.(j == jmax)) &
76  stop ' calc_temp_ssa: Boundary points not allowed.'
77 
78 !-------- Abbreviations --------
79 
80 ctr1 = atr1
81 clb1 = alb1*q_geo(j,i)
82 
83 #if ADV_VERT==1
84 
85 do kc=1, kcmax-1
86  ct1(kc) = at1(kc)/h_c(j,i)*0.5_dp*(vz_c(kc,j,i)+vz_c(kc-1,j,i))
87 end do
88 
89 kc=0
90 ct1_sg(kc) = 0.5_dp*(at1(kc)+at1(kc+1))/h_c(j,i)*vz_c(kc,j,i)
91  ! only needed for kc=0 ...
92 kc=kcmax-1
93 ct1_sg(kc) = 0.5_dp*(at1(kc)+at1(kc+1))/h_c(j,i)*vz_c(kc,j,i)
94  ! ... and kc=KCMAX-1
95 
96 #elif ( ADV_VERT==2 || ADV_VERT==3 )
97 
98 do kc=0, kcmax-1
99  ct1_sg(kc) = 0.5_dp*(at1(kc)+at1(kc+1))/h_c(j,i)*vz_c(kc,j,i)
100 end do
101 
102 #endif
103 
104 do kc=0, kcmax
105 
106  ct2(kc) = ( at2_1(kc)*dzm_dtau(j,i) &
107  +at2_2(kc)*dh_c_dtau(j,i) )/h_c(j,i)
108  ct3(kc) = ( at3_1(kc)*dzm_dxi_g(j,i) &
109  +at3_2(kc)*dh_c_dxi_g(j,i) )/h_c(j,i) &
110  *0.5_dp*(vx_c(kc,j,i)+vx_c(kc,j,i-1)) *insq_g11_g(j,i)
111  ct4(kc) = ( at4_1(kc)*dzm_deta_g(j,i) &
112  +at4_2(kc)*dh_c_deta_g(j,i) )/h_c(j,i) &
113  *0.5_dp*(vy_c(kc,j,i)+vy_c(kc,j-1,i)) *insq_g22_g(j,i)
114  ct5(kc) = at5(kc) &
115  /c_val(temp_c(kc,j,i)) &
116  /h_c(j,i)
117  ct7(kc) = 2.0_dp*at7 &
118  /c_val(temp_c(kc,j,i)) &
119  *viscosity(de_ssa(j,i), &
120  temp_c(kc,j,i), temp_c_m(kc,j,i), 0.0_dp, &
121  enh_c(kc,j,i), 0_i2b) &
122  *de_ssa(j,i)**2
123  ci1(kc) = ai1(kc)/h_c(j,i)
124 
125 end do
126 
127 #if ADV_VERT==1
128 
129 kc=0
130 ct2_sg(kc) = 0.5_dp*(ct2(kc)+ct2(kc+1))
131 ct3_sg(kc) = 0.5_dp*(ct3(kc)+ct3(kc+1))
132 ct4_sg(kc) = 0.5_dp*(ct4(kc)+ct4(kc+1))
133 adv_vert_sg(kc) = ct1_sg(kc)-ct2_sg(kc)-ct3_sg(kc)-ct4_sg(kc)
134 abs_adv_vert_sg(kc) = abs(adv_vert_sg(kc)) ! only needed for kc=0 ...
135 kc=kcmax-1
136 ct2_sg(kc) = 0.5_dp*(ct2(kc)+ct2(kc+1))
137 ct3_sg(kc) = 0.5_dp*(ct3(kc)+ct3(kc+1))
138 ct4_sg(kc) = 0.5_dp*(ct4(kc)+ct4(kc+1))
139 adv_vert_sg(kc) = ct1_sg(kc)-ct2_sg(kc)-ct3_sg(kc)-ct4_sg(kc)
140 abs_adv_vert_sg(kc) = abs(adv_vert_sg(kc)) ! ... and kc=KCMAX-1
141 
142 #elif ( ADV_VERT==2 || ADV_VERT==3 )
143 
144 do kc=0, kcmax-1
145  ct2_sg(kc) = 0.5_dp*(ct2(kc)+ct2(kc+1))
146  ct3_sg(kc) = 0.5_dp*(ct3(kc)+ct3(kc+1))
147  ct4_sg(kc) = 0.5_dp*(ct4(kc)+ct4(kc+1))
148  adv_vert_sg(kc) = ct1_sg(kc)-ct2_sg(kc)-ct3_sg(kc)-ct4_sg(kc)
149  abs_adv_vert_sg(kc) = abs(adv_vert_sg(kc))
150 end do
151 
152 #endif
153 
154 do kc=0, kcmax-1
155  temp_c_help(kc) = 0.5_dp*(temp_c(kc,j,i)+temp_c(kc+1,j,i))
156  ct6(kc) = at6(kc) &
157  *kappa_val(temp_c_help(kc)) &
158  /h_c(j,i)
159  ci2(kc) = ai2(kc)/h_c(j,i)
160 end do
161 
162 #if (ADV_HOR==3)
163 dtt_dxi = 2.0_dp*dtt_2dxi
164 dtt_deta = 2.0_dp*dtt_2deta
165 #endif
166 
167 !-------- Set up the equations for the bedrock temperature --------
168 
169 kr=0
170 lgs_a1(kr) = 1.0_dp
171 lgs_a2(kr) = -1.0_dp
172 lgs_b(kr) = clb1
173 
174 #if Q_LITHO==1
175 ! (coupled heat-conducting bedrock)
176 
177 do kr=1, krmax-1
178  lgs_a0(kr) = - ctr1
179  lgs_a1(kr) = 1.0_dp + 2.0_dp*ctr1
180  lgs_a2(kr) = - ctr1
181  lgs_b(kr) = temp_r(kr,j,i)
182 end do
183 
184 #elif Q_LITHO==0
185 ! (no coupled heat-conducting bedrock)
186 
187 do kr=1, krmax-1
188  lgs_a0(kr) = 1.0_dp
189  lgs_a1(kr) = 0.0_dp
190  lgs_a2(kr) = -1.0_dp
191  lgs_b(kr) = 2.0_dp*clb1
192 end do
193 
194 #endif
195 
196 kr=krmax
197 lgs_a0(kr) = 0.0_dp
198 lgs_a1(kr) = 1.0_dp
199 lgs_b(kr) = temp_c_m(0,j,i)-delta_tm_sw
200 
201 !-------- Solve system of linear equations --------
202 
203 call tri_sle(lgs_a0, lgs_a1, lgs_a2, lgs_x, lgs_b, krmax)
204 
205 !-------- Assign the result --------
206 
207 do kr=0, krmax
208  temp_r_neu(kr,j,i) = lgs_x(kr)
209 end do
210 
211 !-------- Set up the equations for the ice temperature --------
212 
213 kc=0
214 lgs_a1(kc) = 1.0_dp
215 lgs_a2(kc) = 0.0_dp
216 lgs_b(kc) = temp_c_m(0,j,i)-delta_tm_sw
217 
218 do kc=1, kcmax-1
219 
220 #if ADV_VERT==1
221 
222  lgs_a0(kc) = -0.5_dp*(ct1(kc)-ct2(kc)-ct3(kc)-ct4(kc)) &
223  -ct5(kc)*ct6(kc-1)
224  lgs_a1(kc) = 1.0_dp+ct5(kc)*(ct6(kc)+ct6(kc-1))
225  lgs_a2(kc) = 0.5_dp*(ct1(kc)-ct2(kc)-ct3(kc)-ct4(kc)) &
226  -ct5(kc)*ct6(kc)
227 
228 #elif ADV_VERT==2
229 
230  lgs_a0(kc) &
231  = -0.5_dp*(adv_vert_sg(kc-1)+abs_adv_vert_sg(kc-1)) &
232  -ct5(kc)*ct6(kc-1)
233  lgs_a1(kc) &
234  = 1.0_dp &
235  +0.5_dp*(adv_vert_sg(kc-1)+abs_adv_vert_sg(kc-1)) &
236  -0.5_dp*(adv_vert_sg(kc) -abs_adv_vert_sg(kc) ) &
237  +ct5(kc)*(ct6(kc)+ct6(kc-1))
238  lgs_a2(kc) &
239  = 0.5_dp*(adv_vert_sg(kc) -abs_adv_vert_sg(kc) ) &
240  -ct5(kc)*ct6(kc)
241 
242 #elif ADV_VERT==3
243 
244  adv_vert_help = 0.5_dp*(adv_vert_sg(kc)+adv_vert_sg(kc-1))
245 
246  lgs_a0(kc) &
247  = -max(adv_vert_help, 0.0_dp) &
248  -ct5(kc)*ct6(kc-1)
249  lgs_a1(kc) &
250  = 1.0_dp &
251  +max(adv_vert_help, 0.0_dp)-min(adv_vert_help, 0.0_dp) &
252  +ct5(kc)*(ct6(kc)+ct6(kc-1))
253  lgs_a2(kc) &
254  = min(adv_vert_help, 0.0_dp) &
255  -ct5(kc)*ct6(kc)
256 
257 #endif
258 
259 #if ADV_HOR==2
260 
261  lgs_b(kc) = temp_c(kc,j,i) + ct7(kc) &
262  -dtt_2dxi* &
263  ( (vx_c(kc,j,i)-abs(vx_c(kc,j,i))) &
264  *(temp_c(kc,j,i+1)-temp_c(kc,j,i)) &
265  *insq_g11_sgx(j,i) &
266  +(vx_c(kc,j,i-1)+abs(vx_c(kc,j,i-1))) &
267  *(temp_c(kc,j,i)-temp_c(kc,j,i-1)) &
268  *insq_g11_sgx(j,i-1) ) &
269  -dtt_2deta* &
270  ( (vy_c(kc,j,i)-abs(vy_c(kc,j,i))) &
271  *(temp_c(kc,j+1,i)-temp_c(kc,j,i)) &
272  *insq_g22_sgy(j,i) &
273  +(vy_c(kc,j-1,i)+abs(vy_c(kc,j-1,i))) &
274  *(temp_c(kc,j,i)-temp_c(kc,j-1,i)) &
275  *insq_g22_sgy(j-1,i) )
276 
277 #elif ADV_HOR==3
278 
279  vx_c_help = 0.5_dp*(vx_c(kc,j,i)+vx_c(kc,j,i-1))
280  vy_c_help = 0.5_dp*(vy_c(kc,j,i)+vy_c(kc,j-1,i))
281 
282  lgs_b(kc) = temp_c(kc,j,i) + ct7(kc) &
283  -dtt_dxi* &
284  ( min(vx_c_help, 0.0_dp) &
285  *(temp_c(kc,j,i+1)-temp_c(kc,j,i)) &
286  *insq_g11_sgx(j,i) &
287  +max(vx_c_help, 0.0_dp) &
288  *(temp_c(kc,j,i)-temp_c(kc,j,i-1)) &
289  *insq_g11_sgx(j,i-1) ) &
290  -dtt_deta* &
291  ( min(vy_c_help, 0.0_dp) &
292  *(temp_c(kc,j+1,i)-temp_c(kc,j,i)) &
293  *insq_g22_sgy(j,i) &
294  +max(vy_c_help, 0.0_dp) &
295  *(temp_c(kc,j,i)-temp_c(kc,j-1,i)) &
296  *insq_g22_sgy(j-1,i) )
297 
298 #endif
299 
300 end do
301 
302 kc=kcmax
303 lgs_a0(kc) = 0.0_dp
304 lgs_a1(kc) = 1.0_dp
305 lgs_b(kc) = temp_s(j,i)
306 
307 !-------- Solve system of linear equations --------
308 
309 call tri_sle(lgs_a0, lgs_a1, lgs_a2, lgs_x, lgs_b, kcmax)
310 
311 !-------- Assign the result --------
312 
313 do kc=0, kcmax
314  temp_c_neu(kc,j,i) = lgs_x(kc)
315 end do
316 
317 !-------- Set water content in the non-existing temperate layer
318 ! to zero --------
319 
320 do kt=0, ktmax
321  omega_t_neu(kt,j,i) = 0.0_dp
322 end do
323 
324 !-------- Water drainage from the non-existing temperate layer --------
325 
326 q_tld(j,i) = 0.0_dp
327 
328 !-------- Set up the equations for the age of cold ice --------
329 
330 kc=0 ! adv_vert_sg(0) <= 0
331 lgs_a1(kc) = 1.0_dp - min(adv_vert_sg(kc), 0.0_dp) ! (directed downward)
332 lgs_a2(kc) = min(adv_vert_sg(kc), 0.0_dp) ! assumed/enforced
333 
334 #if ADV_HOR==2
335 
336 lgs_b(kc) = age_c(kc,j,i) + dtime_temp &
337  -dtt_2dxi* &
338  ( (vx_c(kc,j,i)-abs(vx_c(kc,j,i))) &
339  *(age_c(kc,j,i+1)-age_c(kc,j,i)) &
340  *insq_g11_sgx(j,i) &
341  +(vx_c(kc,j,i-1)+abs(vx_c(kc,j,i-1))) &
342  *(age_c(kc,j,i)-age_c(kc,j,i-1)) &
343  *insq_g11_sgx(j,i-1) ) &
344  -dtt_2deta* &
345  ( (vy_c(kc,j,i)-abs(vy_c(kc,j,i))) &
346  *(age_c(kc,j+1,i)-age_c(kc,j,i)) &
347  *insq_g22_sgy(j,i) &
348  +(vy_c(kc,j-1,i)+abs(vy_c(kc,j-1,i))) &
349  *(age_c(kc,j,i)-age_c(kc,j-1,i)) &
350  *insq_g22_sgy(j-1,i) )
351 
352 #elif ADV_HOR==3
353 
354 vx_c_help = 0.5_dp*(vx_c(kc,j,i)+vx_c(kc,j,i-1))
355 vy_c_help = 0.5_dp*(vy_c(kc,j,i)+vy_c(kc,j-1,i))
356 
357 lgs_b(kc) = age_c(kc,j,i) + dtime_temp &
358  -dtt_dxi* &
359  ( min(vx_c_help, 0.0_dp) &
360  *(age_c(kc,j,i+1)-age_c(kc,j,i)) &
361  *insq_g11_sgx(j,i) &
362  +max(vx_c_help, 0.0_dp) &
363  *(age_c(kc,j,i)-age_c(kc,j,i-1)) &
364  *insq_g11_sgx(j,i-1) ) &
365  -dtt_deta* &
366  ( min(vy_c_help, 0.0_dp) &
367  *(age_c(kc,j+1,i)-age_c(kc,j,i)) &
368  *insq_g22_sgy(j,i) &
369  +max(vy_c_help, 0.0_dp) &
370  *(age_c(kc,j,i)-age_c(kc,j-1,i)) &
371  *insq_g22_sgy(j-1,i) )
372 
373 #endif
374 
375 do kc=1, kcmax-1
376 
377 #if ADV_VERT==1
378 
379  lgs_a0(kc) = -0.5_dp*(ct1(kc)-ct2(kc)-ct3(kc)-ct4(kc)) &
380  -ci1(kc)*ci2(kc-1)
381  lgs_a1(kc) = 1.0_dp+ci1(kc)*(ci2(kc)+ci2(kc-1))
382  lgs_a2(kc) = 0.5_dp*(ct1(kc)-ct2(kc)-ct3(kc)-ct4(kc)) &
383  -ci1(kc)*ci2(kc)
384 
385 #elif ADV_VERT==2
386 
387  lgs_a0(kc) = -0.5_dp*(adv_vert_sg(kc-1)+abs_adv_vert_sg(kc-1))
388  lgs_a1(kc) = 1.0_dp &
389  +0.5_dp*(adv_vert_sg(kc-1)+abs_adv_vert_sg(kc-1)) &
390  -0.5_dp*(adv_vert_sg(kc) -abs_adv_vert_sg(kc) )
391  lgs_a2(kc) = 0.5_dp*(adv_vert_sg(kc) -abs_adv_vert_sg(kc) )
392 
393 #elif ADV_VERT==3
394 
395  adv_vert_help = 0.5_dp*(adv_vert_sg(kc)+adv_vert_sg(kc-1))
396 
397  lgs_a0(kc) = -max(adv_vert_help, 0.0_dp)
398  lgs_a1(kc) = 1.0_dp &
399  +max(adv_vert_help, 0.0_dp)-min(adv_vert_help, 0.0_dp)
400  lgs_a2(kc) = min(adv_vert_help, 0.0_dp)
401 
402 #endif
403 
404 #if ADV_HOR==2
405 
406  lgs_b(kc) = age_c(kc,j,i) + dtime_temp &
407  -dtt_2dxi* &
408  ( (vx_c(kc,j,i)-abs(vx_c(kc,j,i))) &
409  *(age_c(kc,j,i+1)-age_c(kc,j,i)) &
410  *insq_g11_sgx(j,i) &
411  +(vx_c(kc,j,i-1)+abs(vx_c(kc,j,i-1))) &
412  *(age_c(kc,j,i)-age_c(kc,j,i-1)) &
413  *insq_g11_sgx(j,i-1) ) &
414  -dtt_2deta* &
415  ( (vy_c(kc,j,i)-abs(vy_c(kc,j,i))) &
416  *(age_c(kc,j+1,i)-age_c(kc,j,i)) &
417  *insq_g22_sgy(j,i) &
418  +(vy_c(kc,j-1,i)+abs(vy_c(kc,j-1,i))) &
419  *(age_c(kc,j,i)-age_c(kc,j-1,i)) &
420  *insq_g22_sgy(j-1,i) )
421 
422 #elif ADV_HOR==3
423 
424  vx_c_help = 0.5_dp*(vx_c(kc,j,i)+vx_c(kc,j,i-1))
425  vy_c_help = 0.5_dp*(vy_c(kc,j,i)+vy_c(kc,j-1,i))
426 
427  lgs_b(kc) = age_c(kc,j,i) + dtime_temp &
428  -dtt_dxi* &
429  ( min(vx_c_help, 0.0_dp) &
430  *(age_c(kc,j,i+1)-age_c(kc,j,i)) &
431  *insq_g11_sgx(j,i) &
432  +max(vx_c_help, 0.0_dp) &
433  *(age_c(kc,j,i)-age_c(kc,j,i-1)) &
434  *insq_g11_sgx(j,i-1) ) &
435  -dtt_deta* &
436  ( min(vy_c_help, 0.0_dp) &
437  *(age_c(kc,j+1,i)-age_c(kc,j,i)) &
438  *insq_g22_sgy(j,i) &
439  +max(vy_c_help, 0.0_dp) &
440  *(age_c(kc,j,i)-age_c(kc,j-1,i)) &
441  *insq_g22_sgy(j-1,i) )
442 
443 #endif
444 
445 end do
446 
447 kc=kcmax
448 if (as_perp(j,i) >= zero) then
449  lgs_a0(kc) = 0.0_dp
450  lgs_a1(kc) = 1.0_dp
451  lgs_b(kc) = 0.0_dp
452 else
453  lgs_a0(kc) = -max(adv_vert_sg(kc-1), 0.0_dp)
454  lgs_a1(kc) = 1.0_dp + max(adv_vert_sg(kc-1), 0.0_dp)
455  ! adv_vert_sg(KCMAX-1) >= 0 (directed upward)
456  ! assumed/enforced
457 #if ADV_HOR==2
458 
459  lgs_b(kc) = age_c(kc,j,i) + dtime_temp &
460  -dtt_2dxi* &
461  ( (vx_c(kc,j,i)-abs(vx_c(kc,j,i))) &
462  *(age_c(kc,j,i+1)-age_c(kc,j,i)) &
463  *insq_g11_sgx(j,i) &
464  +(vx_c(kc,j,i-1)+abs(vx_c(kc,j,i-1))) &
465  *(age_c(kc,j,i)-age_c(kc,j,i-1)) &
466  *insq_g11_sgx(j,i-1) ) &
467  -dtt_2deta* &
468  ( (vy_c(kc,j,i)-abs(vy_c(kc,j,i))) &
469  *(age_c(kc,j+1,i)-age_c(kc,j,i)) &
470  *insq_g22_sgy(j,i) &
471  +(vy_c(kc,j-1,i)+abs(vy_c(kc,j-1,i))) &
472  *(age_c(kc,j,i)-age_c(kc,j-1,i)) &
473  *insq_g22_sgy(j-1,i) )
474 
475 #elif ADV_HOR==3
476 
477  vx_c_help = 0.5_dp*(vx_c(kc,j,i)+vx_c(kc,j,i-1))
478  vy_c_help = 0.5_dp*(vy_c(kc,j,i)+vy_c(kc,j-1,i))
479 
480  lgs_b(kc) = age_c(kc,j,i) + dtime_temp &
481  -dtt_dxi* &
482  ( min(vx_c_help, 0.0_dp) &
483  *(age_c(kc,j,i+1)-age_c(kc,j,i)) &
484  *insq_g11_sgx(j,i) &
485  +max(vx_c_help, 0.0_dp) &
486  *(age_c(kc,j,i)-age_c(kc,j,i-1)) &
487  *insq_g11_sgx(j,i-1) ) &
488  -dtt_deta* &
489  ( min(vy_c_help, 0.0_dp) &
490  *(age_c(kc,j+1,i)-age_c(kc,j,i)) &
491  *insq_g22_sgy(j,i) &
492  +max(vy_c_help, 0.0_dp) &
493  *(age_c(kc,j,i)-age_c(kc,j-1,i)) &
494  *insq_g22_sgy(j-1,i) )
495 
496 #endif
497 
498 end if
499 
500 !-------- Solve system of linear equations --------
501 
502 call tri_sle(lgs_a0, lgs_a1, lgs_a2, lgs_x, lgs_b, kcmax)
503 
504 !-------- Assign the result,
505 ! restriction to interval [0, AGE_MAX yr] --------
506 
507 do kc=0, kcmax
508 
509  age_c_neu(kc,j,i) = lgs_x(kc)
510 
511  if (age_c_neu(kc,j,i) < (age_min*year_sec)) &
512  age_c_neu(kc,j,i) = 0.0_dp
513  if (age_c_neu(kc,j,i) > (age_max*year_sec)) &
514  age_c_neu(kc,j,i) = age_max*year_sec
515 
516 end do
517 
518 !-------- Age of the ice in the non-existing temperate layer --------
519 
520 do kt=0, ktmax
521  age_t_neu(kt,j,i) = age_c_neu(0,j,i)
522 end do
523 
524 end subroutine calc_temp_ssa
525 !
sico_sle_solvers::tri_sle
subroutine tri_sle(a0, a1, a2, x, b, nrows)
Solution of a system of linear equations Ax=b with tridiagonal matrix A.
Definition: sico_sle_solvers.F90:122
sico_types
Declarations of kind types for SICOPOLIS.
Definition: sico_types.F90:35
viscosity
real(dp) function viscosity(de_val, temp_val, temp_m_val, omega_val, enh_val, i_flag_cold_temp)
Ice viscosity as a function of the effective strain rate and the temperature (in cold ice) or the wat...
Definition: viscosity.F90:39
sico_vars
Declarations of global variables for SICOPOLIS (for the ANT domain).
Definition: sico_vars.F90:35
sico_sle_solvers
Solvers for systems of linear equations used by SICOPOLIS.
Definition: sico_sle_solvers.F90:35
ice_material_quantities::kappa_val
real(dp) function, public kappa_val(temp_val)
Heat conductivity of ice: Linear interpolation of tabulated values in KAPPA(.).
Definition: ice_material_quantities.F90:217
ice_material_quantities::c_val
real(dp) function, public c_val(temp_val)
Specific heat of ice: Linear interpolation of tabulated values in C(.).
Definition: ice_material_quantities.F90:254
ice_material_quantities
Material quantities of ice: Rate factor, heat conductivity, specific heat (heat capacity).
Definition: ice_material_quantities.F90:37
calc_temp_ssa
subroutine calc_temp_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).
Definition: calc_temp_ssa.F90:35
sico_variables
Declarations of global variables for SICOPOLIS.
Definition: sico_variables.F90:35
Generated on Wed Mar 2 2016 11:26:55 for SICOPOLIS V3.2 by   doxygen 1.8.5