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