SICOPOLIS V3.1
 All Classes Files Functions Variables Macros
shift_cts_downward.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Subroutine : s h i f t _ c t s _ d o w n w a r d
4 !
5 !> @file
6 !!
7 !! Downward shifting of the CTS.
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2009-2013 Ralf Greve
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 !> Downward shifting of the CTS.
34 !<------------------------------------------------------------------------------
35 subroutine shift_cts_downward(at1, at2_1, at2_2, at3_1, at3_2, &
36  at4_1, at4_2, at5, at6, at7, atr1, am1, am2, alb1, &
37  aw1, aw2, aw3, aw4, aw5, aw7, aw8, aw9, aqtld, &
38  ai1, ai2, ai3, ai4, mean_accum_inv, dzeta_t, &
39  dtime_temp, dtt_2dxi, dtt_2deta, dtime_temp_inv, &
40  i, j)
41 
42 use sico_types
43 use sico_variables
44 
45 implicit none
46 
47 integer(i4b), intent(in) :: i, j
48 real(dp), intent(in) :: 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), ai3, ai4, &
52  atr1, am1, am2, alb1
53 real(dp), intent(in) :: aw1, aw2, aw3, aw4, aw5, aw7, aw8, aw9, aqtld
54 real(dp), intent(in) :: mean_accum_inv
55 real(dp), intent(in) :: dzeta_t
56 real(dp), intent(in) :: dtt_2dxi, dtt_2deta, dtime_temp, dtime_temp_inv
57 
58 real(dp) :: zm_shift
59 real(dp) :: difftemp_a, difftemp_b, interpol
60 
61 zm_shift = 1.0_dp ! CTS shift in intervals of 1 m
62 
63 !-------- Temperature discrepancy from the computation of the main
64 ! program --------
65 
66 difftemp_a = temp_c_neu(0,j,i)-(-beta*h_c_neu(j,i))
67 if (difftemp_a >= 0.0_dp) return
68 
69 !-------- Shift CTS downward until it is too low --------
70 
71  10 continue
72 
73  zm_neu(j,i) = zm_neu(j,i) - zm_shift
74 
75 ! ------ Special case: CTS too close to the base
76 
77  if (zm_neu(j,i) <= zb(j,i)) then
78 
79  zm_shift = (zm_neu(j,i)+zm_shift)-(zb(j,i)+0.001_dp)
80  zm_neu(j,i) = zb(j,i)+0.001_dp
81  h_c_neu(j,i) = h_c_neu(j,i)+h_t_neu(j,i)-0.001_dp
82  h_t_neu(j,i) = 0.001_dp
83 ! ! CTS positioned 1 mm above ice base --------
84 
85  dh_t_dtau(j,i) = (zm_neu(j,i)-zm(j,i))*dtime_temp_inv
86  dzm_dtau(j,i) = dzb_dtau(j,i)+dh_t_dtau(j,i)
87  dh_c_dtau(j,i) = dzs_dtau(j,i)-dzm_dtau(j,i)
88 
89  am_perp(j,i) = am_perp_st(j,i) + dzm_dtau(j,i)
90 
91  call calc_temp3(at1, at2_1, at2_2, at3_1, at3_2, &
92  at4_1, at4_2, at5, at6, at7, atr1, &
93  am1, am2, alb1, &
94  aw1, aw2, aw3, aw4, aw5, aw7, aw8, aw9, aqtld, &
95  ai1, ai2, ai3, mean_accum_inv, dzeta_t, &
96  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
97 
98  difftemp_b = difftemp_a
99  difftemp_a = temp_c_neu(0,j,i)-(-beta*h_c_neu(j,i))
100 
101  if (difftemp_a >= 0.0_dp) then ! CTS remains above the base
102 
103 ! ---- Interpolate the CTS position from the last (_a) and the
104 ! last but one (_b) value, weighed with the temperature
105 ! discrepancies at the CTS --------
106 
107  interpol = difftemp_a/(difftemp_a-difftemp_b)*zm_shift
108 
109  zm_neu(j,i) = zm_neu(j,i) + interpol
110  h_c_neu(j,i) = h_c_neu(j,i) - interpol
111  h_t_neu(j,i) = h_t_neu(j,i) + interpol
112 
113  dh_t_dtau(j,i) = (zm_neu(j,i)-zm(j,i))*dtime_temp_inv
114  dzm_dtau(j,i) = dzb_dtau(j,i)+dh_t_dtau(j,i)
115  dh_c_dtau(j,i) = dzs_dtau(j,i)-dzm_dtau(j,i)
116 
117  am_perp(j,i) = am_perp_st(j,i) + dzm_dtau(j,i)
118 
119  call calc_temp3(at1, at2_1, at2_2, at3_1, at3_2, &
120  at4_1, at4_2, at5, at6, at7, atr1, &
121  am1, am2, alb1, &
122  aw1, aw2, aw3, aw4, aw5, aw7, aw8, aw9, aqtld, &
123  ai1, ai2, ai3, mean_accum_inv, dzeta_t, &
124  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
125 
126  else ! CTS disappears
127 
128  n_cts_neu(j,i) = 0
129  zm_neu(j,i) = zb(j,i)
130  h_c_neu(j,i) = h_c_neu(j,i)+h_t_neu(j,i)
131  h_t_neu(j,i) = 0.0_dp
132 
133  dh_t_dtau(j,i) = (zm_neu(j,i)-zm(j,i))*dtime_temp_inv
134  dzm_dtau(j,i) = dzb_dtau(j,i)+dh_t_dtau(j,i)
135  dh_c_dtau(j,i) = dzs_dtau(j,i)-dzm_dtau(j,i)
136 
137  am_perp(j,i) = am_perp_st(j,i) + dzm_dtau(j,i)
138 
139  call calc_temp2(at1, at2_1, at2_2, at3_1, at3_2, &
140  at4_1, at4_2, at5, at6, at7, atr1, alb1, &
141  ai1, ai2, ai4, mean_accum_inv, &
142  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
143 
144  end if
145 
146  return
147  end if
148 
149 ! ------ End of treatment of special case
150 
151  h_c_neu(j,i) = h_c_neu(j,i) + zm_shift
152  h_t_neu(j,i) = h_t_neu(j,i) - zm_shift
153 
154  dh_t_dtau(j,i) = (zm_neu(j,i)-zm(j,i))*dtime_temp_inv
155  dzm_dtau(j,i) = dzb_dtau(j,i)+dh_t_dtau(j,i)
156  dh_c_dtau(j,i) = dzs_dtau(j,i)-dzm_dtau(j,i)
157 
158  am_perp(j,i) = am_perp_st(j,i) + dzm_dtau(j,i)
159 
160  call calc_temp3(at1, at2_1, at2_2, at3_1, at3_2, &
161  at4_1, at4_2, at5, at6, at7, atr1, &
162  am1, am2, alb1, &
163  aw1, aw2, aw3, aw4, aw5, aw7, aw8, aw9, aqtld, &
164  ai1, ai2, ai3, mean_accum_inv, dzeta_t, &
165  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
166 
167  difftemp_b = difftemp_a
168  difftemp_a = temp_c_neu(0,j,i)-(-beta*h_c_neu(j,i))
169 
170 if (difftemp_a < 0.0_dp) go to 10
171 
172 !-------- Interpolate the CTS position from the last (_a) and the
173 ! last but one (_b) value, weighed with the temperature
174 ! discrepancies at the CTS --------
175 
176 interpol = difftemp_a/(difftemp_a-difftemp_b)*zm_shift
177 
178 zm_neu(j,i) = zm_neu(j,i) + interpol
179 h_c_neu(j,i) = h_c_neu(j,i) - interpol
180 h_t_neu(j,i) = h_t_neu(j,i) + interpol
181 
182 dh_t_dtau(j,i) = (zm_neu(j,i)-zm(j,i))*dtime_temp_inv
183 dzm_dtau(j,i) = dzb_dtau(j,i)+dh_t_dtau(j,i)
184 dh_c_dtau(j,i) = dzs_dtau(j,i)-dzm_dtau(j,i)
185 
186 am_perp(j,i) = am_perp_st(j,i) + dzm_dtau(j,i)
187 
188 call calc_temp3(at1, at2_1, at2_2, at3_1, at3_2, &
189  at4_1, at4_2, at5, at6, at7, atr1, &
190  am1, am2, alb1, &
191  aw1, aw2, aw3, aw4, aw5, aw7, aw8, aw9, aqtld, &
192  ai1, ai2, ai3, mean_accum_inv, dzeta_t, &
193  dtime_temp, dtt_2dxi, dtt_2deta, i, j)
194 
195 end subroutine shift_cts_downward
196 !
Generated on Thu Sep 12 2013 09:45:18 for SICOPOLIS V3.1 by   doxygen 1.8.2