SICOPOLIS V3.3
calc_enhance_m.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Module : c a l c _ e n h a n c e _ m
4 !
5 !> @file
6 !!
7 !! Computation of the flow enhancement factor.
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2009-2017 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 !> Computation of the flow enhancement factor.
34 !<------------------------------------------------------------------------------
36 
37  use sico_types_m
39  use sico_vars_m
40 
41  implicit none
42 
43  private
46 
47 contains
48 
49 !-------------------------------------------------------------------------------
50 !> Computation of the flow enhancement factor.
51 !! Case ENHMOD==1:
52 !! constant for grounded ice, constant for floating ice.
53 !<------------------------------------------------------------------------------
54  subroutine calc_enhance_1()
55 
56  implicit none
57 
58  enh_t = enh_fact
59  enh_c = enh_fact
60 
61 #if (MARGIN==3) /* floating ice */
62  call calc_enhance_floating_const()
63 #endif
64 
65 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
66  call mod_enhance_dust()
67 #endif
68 
69  end subroutine calc_enhance_1
70 
71 !-------------------------------------------------------------------------------
72 !> Computation of the flow enhancement factor.
73 !! Case ENHMOD==2:
74 !! two different values depending on age for grounded ice,
75 !! constant for floating ice.
76 !<------------------------------------------------------------------------------
77  subroutine calc_enhance_2()
78 
79  implicit none
80 
81  integer(i4b) :: i, j, kc, kt
82  real(dp) :: age_trans
83 
84  age_trans = age_trans_0*year_sec
85 
86  do i=0, imax
87  do j=0, jmax
88 
89  do kt=0, ktmax
90  if (age_t(kt,j,i) < age_trans) then
91  enh_t(kt,j,i) = enh_intg ! Holocene ice
92  else
93  enh_t(kt,j,i) = enh_fact ! Pleistocene ice
94  end if
95  end do
96 
97  do kc=0, kcmax
98  if (age_c(kc,j,i) < age_trans) then
99  enh_c(kc,j,i) = enh_intg ! Holocene ice
100  else
101  enh_c(kc,j,i) = enh_fact ! Pleistocene ice
102  end if
103  end do
104 
105  end do
106  end do
107 
108 #if (MARGIN==3) /* floating ice */
109  call calc_enhance_floating_const()
110 #endif
111 
112 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
113  call mod_enhance_dust()
114 #endif
115 
116  end subroutine calc_enhance_2
117 
118 !-------------------------------------------------------------------------------
119 !> Computation of the flow enhancement factor.
120 !! Case ENHMOD==3:
121 !! two different values depending on time of deposition for grounded ice,
122 !! constant for floating ice.
123 !<------------------------------------------------------------------------------
124  subroutine calc_enhance_3(time)
126  implicit none
127 
128  real(dp), intent(in) :: time
129 
130  integer(i4b) :: i, j, kc, kt
131  real(dp) :: date_trans1, date_trans2, date_trans3
132 
133  date_trans1 = date_trans1_0*year_sec
134  date_trans2 = date_trans2_0*year_sec
135  date_trans3 = date_trans3_0*year_sec
136 
137  do i=0, imax
138  do j=0, jmax
139 
140  do kt=0, ktmax
141  if ( (time-age_t(kt,j,i)) < date_trans1 ) then
142  enh_t(kt,j,i) = enh_fact ! pre-Eemian ice
143  else if ( ((time-age_t(kt,j,i)) >= date_trans1).and. &
144  ((time-age_t(kt,j,i)) < date_trans2) ) then
145  enh_t(kt,j,i) = enh_intg ! Eemian ice
146  else if ( ((time-age_t(kt,j,i)) >= date_trans2).and. &
147  ((time-age_t(kt,j,i)) < date_trans3) ) then
148  enh_t(kt,j,i) = enh_fact ! Weichselian ice
149  else
150  enh_t(kt,j,i) = enh_intg ! Holocene ice
151  end if
152  end do
153 
154  do kc=0, kcmax
155  if ( (time-age_c(kc,j,i)) < date_trans1 ) then
156  enh_c(kc,j,i) = enh_fact ! pre-Eemian ice
157  else if ( ((time-age_c(kc,j,i)) >= date_trans1).and. &
158  ((time-age_c(kc,j,i)) < date_trans2) ) then
159  enh_c(kc,j,i) = enh_intg ! Eemian ice
160  else if ( ((time-age_c(kc,j,i)) >= date_trans2).and. &
161  ((time-age_c(kc,j,i)) < date_trans3) ) then
162  enh_c(kc,j,i) = enh_fact ! Weichselian ice
163  else
164  enh_c(kc,j,i) = enh_intg ! Holocene ice
165  end if
166  end do
167 
168  end do
169  end do
170 
171 #if (MARGIN==3) /* floating ice */
172  call calc_enhance_floating_const()
173 #endif
174 
175 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
176  call mod_enhance_dust()
177 #endif
178 
179  end subroutine calc_enhance_3
180 
181 !-------------------------------------------------------------------------------
182 !> Computation of the flow enhancement factor.
183 !! Case ENHMOD==4:
184 !! minimal anisotropic enhancement factor for grounded ice,
185 !! constant for floating ice.
186 !<------------------------------------------------------------------------------
187  subroutine calc_enhance_4()
189  implicit none
190 
191  call calc_enhance_aniso()
192 
193 #if (MARGIN==3) /* floating ice */
194  call calc_enhance_floating_const()
195 #endif
196 
197 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
198  call mod_enhance_dust()
199 #endif
200 
201  end subroutine calc_enhance_4
202 
203 !-------------------------------------------------------------------------------
204 !> Computation of the flow enhancement factor.
205 !! Case ENHMOD==5:
206 !! minimal anisotropic enhancement factor for grounded and floating ice.
207 !<------------------------------------------------------------------------------
208  subroutine calc_enhance_5()
210  implicit none
211 
212  call calc_enhance_aniso()
213 
214 #if (defined(NMARS) || defined(SMARS)) /* Martian polar caps */
215  call mod_enhance_dust()
216 #endif
217 
218  end subroutine calc_enhance_5
219 
220 !-------------------------------------------------------------------------------
221 !> Minimal anisotropic flow enhancement factor.
222 !<------------------------------------------------------------------------------
223  subroutine calc_enhance_aniso()
224 
225  implicit none
226 
227  integer(i4b) :: i, j, kc, kt
228  real(dp) :: enh_shear, enh_compr
229  real(dp) :: enh_shear_compr_diff
230 
231 #if (defined(ENH_SHEAR))
232  enh_shear = enh_shear
233 #else
234  enh_shear = 1.0_dp
235 #endif
236 
237 #if (defined(ENH_COMPR))
238  enh_compr = enh_compr
239 #else
240  enh_compr = 1.0_dp
241 #endif
242 
243  enh_shear_compr_diff = enh_shear-enh_compr
244 
245  do i=0, imax
246  do j=0, jmax
247 
248  do kt=0, ktmax
249  enh_t(kt,j,i) = enh_compr &
250  + enh_shear_compr_diff*lambda_shear_t(kt,j,i)**2
251  end do
252 
253  do kc=0, kcmax
254  enh_c(kc,j,i) = enh_compr &
255  + enh_shear_compr_diff*lambda_shear_c(kc,j,i)**2
256  end do
257 
258  end do
259  end do
260 
261  end subroutine calc_enhance_aniso
262 
263 !-------------------------------------------------------------------------------
264 !> Constant, prescribed flow enhancement factor for floating ice.
265 !<------------------------------------------------------------------------------
266  subroutine calc_enhance_floating_const()
267 
268  implicit none
269 
270  integer(i4b) :: i, j, kc, kt
271  real(dp) :: enh_shelf
272 
273 #if (defined(ENH_SHELF))
274  enh_shelf = enh_shelf
275 #else
276  enh_shelf = 1.0_dp
277 #endif
278 
279  do i=0, imax
280  do j=0, jmax
281 
282  if ( maske(j,i)==3_i2b ) then ! floating ice
283 
284  do kt=0, ktmax
285  enh_t(kt,j,i) = enh_shelf
286  end do
287 
288  do kc=0, kcmax
289  enh_c(kc,j,i) = enh_shelf
290  end do
291 
292  end if
293 
294  end do
295  end do
296 
297  end subroutine calc_enhance_floating_const
298 
299 !-------------------------------------------------------------------------------
300 !> Modification of the flow enhancement factor due to dust content
301 !! (for the Martian polar caps).
302 !<------------------------------------------------------------------------------
303  subroutine mod_enhance_dust()
304 
305  implicit none
306 
307  real(dp) :: frac_dust
308  real(dp) :: enh_mult
309 
310 #if (defined(FRAC_DUST))
311  frac_dust = frac_dust
312 #else
313  frac_dust = 0.0_dp
314 #endif
315 
316 #if (FLOW_LAW==1)
317  enh_mult = exp(-3.0_dp*2.0_dp*frac_dust)
318 #elif (FLOW_LAW==2)
319  enh_mult = exp(-1.8_dp*2.0_dp*frac_dust)
320 #elif (FLOW_LAW==3)
321  enh_mult = exp(-4.0_dp*2.0_dp*frac_dust)
322 #elif (FLOW_LAW==4)
323  stop ' >>> mod_enhance_dust: FLOW_LAW==4 has not been implemented yet!'
324 #endif
325 
326  enh_t = enh_t * enh_mult
327  enh_c = enh_c * enh_mult
328 
329  end subroutine mod_enhance_dust
330 
331 !-------------------------------------------------------------------------------
332 
333 end module calc_enhance_m
334 !
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) age_c
age_c(kc,j,i): Age in the upper (kc) ice domain
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) lambda_shear_t
lambda_shear_t(kt,j,i): Shear fraction in the lower (kt) ice domain
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) age_t
age_t(kt,j,i): Age in the lower (kt) ice domain
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) enh_c
enh_c(kc,j,i): Flow enhancement factor in the upper (kc) ice domain
subroutine, public calc_enhance_2()
Computation of the flow enhancement factor. Case ENHMOD==2: two different values depending on age for...
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 ...
subroutine, public calc_enhance_5()
Computation of the flow enhancement factor. Case ENHMOD==5: minimal anisotropic enhancement factor fo...
Declarations of global variables for SICOPOLIS (for the ANT domain).
Definition: sico_vars_m.F90:35
Computation of the flow enhancement factor.
real(dp), dimension(0:kcmax, 0:jmax, 0:imax) lambda_shear_c
lambda_shear_c(kc,j,i): Shear fraction in the upper (kc) ice domain
Declarations of kind types for SICOPOLIS.
real(dp), dimension(0:ktmax, 0:jmax, 0:imax) enh_t
enh_t(kt,j,i): Flow enhancement factor in the lower (kt) ice domain
subroutine, public calc_enhance_1()
Computation of the flow enhancement factor. Case ENHMOD==1: constant for grounded ice...
subroutine, public calc_enhance_3(time)
Computation of the flow enhancement factor. Case ENHMOD==3: two different values depending on time of...
subroutine, public calc_enhance_4()
Computation of the flow enhancement factor. Case ENHMOD==4: minimal anisotropic enhancement factor fo...
Declarations of global variables for SICOPOLIS.