New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
icethd_pnd.F90 in NEMO/trunk/src/ICE – NEMO

source: NEMO/trunk/src/ICE/icethd_pnd.F90 @ 12489

Last change on this file since 12489 was 12489, checked in by davestorkey, 4 years ago

Preparation for new timestepping scheme #2390.
Main changes:

  1. Initial euler timestep now handled in stp and not in TRA/DYN routines.
  2. Renaming of all timestep parameters. In summary, the namelist parameter is now rn_Dt and the current timestep is rDt (and rDt_ice, rDt_trc etc).
  3. Renaming of a few miscellaneous parameters, eg. atfp -> rn_atfp (namelist parameter used everywhere) and rau0 -> rho0.

This version gives bit-comparable results to the previous version of the trunk.
  • Property svn:keywords set to Id
File size: 12.0 KB
Line 
1MODULE icethd_pnd 
2   !!======================================================================
3   !!                     ***  MODULE  icethd_pnd   ***
4   !!   sea-ice: Melt ponds on top of sea ice
5   !!======================================================================
6   !! history :       !  2012     (O. Lecomte)       Adaptation from Flocco and Turner
7   !!                 !  2017     (M. Vancoppenolle, O. Lecomte, C. Rousset) Implementation
8   !!            4.0  !  2018     (many people)      SI3 [aka Sea Ice cube]
9   !!----------------------------------------------------------------------
10#if defined key_si3
11   !!----------------------------------------------------------------------
12   !!   'key_si3' :                                     SI3 sea-ice model
13   !!----------------------------------------------------------------------
14   !!   ice_thd_pnd_init : some initialization and namelist read
15   !!   ice_thd_pnd      : main calling routine
16   !!----------------------------------------------------------------------
17   USE phycst         ! physical constants
18   USE dom_oce        ! ocean space and time domain
19   USE ice            ! sea-ice: variables
20   USE ice1D          ! sea-ice: thermodynamics variables
21   USE icetab         ! sea-ice: 1D <==> 2D transformation
22   !
23   USE in_out_manager ! I/O manager
24   USE lib_mpp        ! MPP library
25   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
26   USE timing         ! Timing
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC   ice_thd_pnd_init    ! routine called by icestp.F90
32   PUBLIC   ice_thd_pnd         ! routine called by icestp.F90
33
34   INTEGER ::              nice_pnd    ! choice of the type of pond scheme
35   !                                   ! associated indices:
36   INTEGER, PARAMETER ::   np_pndNO  = 0   ! No pond scheme
37   INTEGER, PARAMETER ::   np_pndCST = 1   ! Constant pond scheme
38   INTEGER, PARAMETER ::   np_pndH12 = 2   ! Evolutive pond scheme (Holland et al. 2012)
39
40   !!----------------------------------------------------------------------
41   !! NEMO/ICE 4.0 , NEMO Consortium (2018)
42   !! $Id$
43   !! Software governed by the CeCILL license (see ./LICENSE)
44   !!----------------------------------------------------------------------
45CONTAINS
46
47   SUBROUTINE ice_thd_pnd
48      !!-------------------------------------------------------------------
49      !!               ***  ROUTINE ice_thd_pnd   ***
50      !!               
51      !! ** Purpose :   change melt pond fraction
52      !!               
53      !! ** Method  :   brut force
54      !!-------------------------------------------------------------------
55      !
56      SELECT CASE ( nice_pnd )
57      !
58      CASE (np_pndCST)   ;   CALL pnd_CST    !==  Constant melt ponds  ==!
59         !
60      CASE (np_pndH12)   ;   CALL pnd_H12    !==  Holland et al 2012 melt ponds  ==!
61         !
62      END SELECT
63      !
64   END SUBROUTINE ice_thd_pnd 
65
66
67   SUBROUTINE pnd_CST 
68      !!-------------------------------------------------------------------
69      !!                ***  ROUTINE pnd_CST  ***
70      !!
71      !! ** Purpose :   Compute melt pond evolution
72      !!
73      !! ** Method  :   Melt pond fraction and thickness are prescribed
74      !!                to non-zero values when t_su = 0C
75      !!
76      !! ** Tunable parameters : pond fraction (rn_apnd), pond depth (rn_hpnd)
77      !!               
78      !! ** Note   : Coupling with such melt ponds is only radiative
79      !!             Advection, ridging, rafting... are bypassed
80      !!
81      !! ** References : Bush, G.W., and Trump, D.J. (2017)
82      !!-------------------------------------------------------------------
83      INTEGER  ::   ji        ! loop indices
84      !!-------------------------------------------------------------------
85      DO ji = 1, npti
86         !
87         IF( a_i_1d(ji) > 0._wp .AND. t_su_1d(ji) >= rt0 ) THEN
88            a_ip_frac_1d(ji) = rn_apnd
89            h_ip_1d(ji)      = rn_hpnd   
90            a_ip_1d(ji)      = a_ip_frac_1d(ji) * a_i_1d(ji)
91         ELSE
92            a_ip_frac_1d(ji) = 0._wp
93            h_ip_1d(ji)      = 0._wp   
94            a_ip_1d(ji)      = 0._wp
95         ENDIF
96         !
97      END DO
98      !
99   END SUBROUTINE pnd_CST
100
101
102   SUBROUTINE pnd_H12
103      !!-------------------------------------------------------------------
104      !!                ***  ROUTINE pnd_H12  ***
105      !!
106      !! ** Purpose    : Compute melt pond evolution
107      !!
108      !! ** Method     : Empirical method. A fraction of meltwater is accumulated in ponds
109      !!                 and sent to ocean when surface is freezing
110      !!
111      !!                 pond growth:      Vp = Vp + dVmelt
112      !!                    with dVmelt = R/rhow * ( rhoi*dh_i + rhos*dh_s ) * a_i
113      !!                 pond contraction: Vp = Vp * exp(0.01*MAX(Tp-Tsu,0)/Tp)
114      !!                    with Tp = -2degC
115      !! 
116      !! ** Tunable parameters : (no real expertise yet, ideas?)
117      !!
118      !! ** Note       : Stolen from CICE for quick test of the melt pond
119      !!                 radiation and freshwater interfaces
120      !!                 Coupling can be radiative AND freshwater
121      !!                 Advection, ridging, rafting are called
122      !!
123      !! ** References : Holland, M. M. et al (J Clim 2012)
124      !!-------------------------------------------------------------------
125      REAL(wp), PARAMETER ::   zrmin       = 0.15_wp  ! minimum fraction of available meltwater retained for melt ponding
126      REAL(wp), PARAMETER ::   zrmax       = 0.70_wp  ! maximum     -           -         -         -            -
127      REAL(wp), PARAMETER ::   zpnd_aspect = 0.8_wp   ! pond aspect ratio
128      REAL(wp), PARAMETER ::   zTp         = -2._wp   ! reference temperature
129      !
130      REAL(wp) ::   zfr_mlt          ! fraction of available meltwater retained for melt ponding
131      REAL(wp) ::   zdv_mlt          ! available meltwater for melt ponding
132      REAL(wp) ::   z1_Tp            ! inverse reference temperature
133      REAL(wp) ::   z1_rhow          ! inverse freshwater density
134      REAL(wp) ::   z1_zpnd_aspect   ! inverse pond aspect ratio
135      REAL(wp) ::   zfac, zdum
136      !
137      INTEGER  ::   ji   ! loop indices
138      !!-------------------------------------------------------------------
139      z1_rhow        = 1._wp / rhow 
140      z1_zpnd_aspect = 1._wp / zpnd_aspect
141      z1_Tp          = 1._wp / zTp 
142
143      DO ji = 1, npti
144         !                                                        !--------------------------------!
145         IF( h_i_1d(ji) < rn_himin) THEN                          ! Case ice thickness < rn_himin  !
146            !                                                     !--------------------------------!
147            !--- Remove ponds on thin ice
148            a_ip_1d(ji)      = 0._wp
149            a_ip_frac_1d(ji) = 0._wp
150            h_ip_1d(ji)      = 0._wp
151            !                                                     !--------------------------------!
152         ELSE                                                     ! Case ice thickness >= rn_himin !
153            !                                                     !--------------------------------!
154            v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji)   ! record pond volume at previous time step
155            !
156            ! available meltwater for melt ponding [m, >0] and fraction
157            zdv_mlt = -( dh_i_sum(ji)*rhoi + dh_s_mlt(ji)*rhos ) * z1_rhow * a_i_1d(ji)
158            zfr_mlt = zrmin + ( zrmax - zrmin ) * a_i_1d(ji)  ! from CICE doc
159            !zfr_mlt = zrmin + zrmax * a_i_1d(ji)             ! from Holland paper
160            !
161            !--- Pond gowth ---!
162            ! v_ip should never be negative, otherwise code crashes
163            v_ip_1d(ji) = MAX( 0._wp, v_ip_1d(ji) + zfr_mlt * zdv_mlt )
164            !
165            ! melt pond mass flux (<0)
166            IF( zdv_mlt > 0._wp ) THEN
167               zfac = zfr_mlt * zdv_mlt * rhow * r1_Dt_ice
168               wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac
169               !
170               ! adjust ice/snow melting flux to balance melt pond flux (>0)
171               zdum = zfac / ( wfx_snw_sum_1d(ji) + wfx_sum_1d(ji) )
172               wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) * (1._wp + zdum)
173               wfx_sum_1d(ji)     = wfx_sum_1d(ji)     * (1._wp + zdum)
174            ENDIF
175            !
176            !--- Pond contraction (due to refreezing) ---!
177            v_ip_1d(ji) = v_ip_1d(ji) * EXP( 0.01_wp * MAX( zTp+rt0 - t_su_1d(ji), 0._wp ) * z1_Tp )
178            !
179            ! Set new pond area and depth assuming linear relation between h_ip and a_ip_frac
180            !    h_ip = zpnd_aspect * a_ip_frac = zpnd_aspect * a_ip/a_i
181            a_ip_1d(ji)      = SQRT( v_ip_1d(ji) * z1_zpnd_aspect * a_i_1d(ji) )
182            a_ip_frac_1d(ji) = a_ip_1d(ji) / a_i_1d(ji)
183            h_ip_1d(ji)      = zpnd_aspect * a_ip_frac_1d(ji)
184            !
185         ENDIF
186      END DO
187      !
188   END SUBROUTINE pnd_H12
189
190
191   SUBROUTINE ice_thd_pnd_init 
192      !!-------------------------------------------------------------------
193      !!                  ***  ROUTINE ice_thd_pnd_init   ***
194      !!
195      !! ** Purpose : Physical constants and parameters linked to melt ponds
196      !!              over sea ice
197      !!
198      !! ** Method  :  Read the namthd_pnd  namelist and check the melt pond 
199      !!               parameter values called at the first timestep (nit000)
200      !!
201      !! ** input   :   Namelist namthd_pnd 
202      !!-------------------------------------------------------------------
203      INTEGER  ::   ios, ioptio   ! Local integer
204      !!
205      NAMELIST/namthd_pnd/  ln_pnd, ln_pnd_H12, ln_pnd_CST, rn_apnd, rn_hpnd, ln_pnd_alb
206      !!-------------------------------------------------------------------
207      !
208      READ  ( numnam_ice_ref, namthd_pnd, IOSTAT = ios, ERR = 901)
209901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namthd_pnd  in reference namelist' )
210      READ  ( numnam_ice_cfg, namthd_pnd, IOSTAT = ios, ERR = 902 )
211902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namthd_pnd in configuration namelist' )
212      IF(lwm) WRITE ( numoni, namthd_pnd )
213      !
214      IF(lwp) THEN                        ! control print
215         WRITE(numout,*)
216         WRITE(numout,*) 'ice_thd_pnd_init: ice parameters for melt ponds'
217         WRITE(numout,*) '~~~~~~~~~~~~~~~~'
218         WRITE(numout,*) '   Namelist namicethd_pnd:'
219         WRITE(numout,*) '      Melt ponds activated or not                                     ln_pnd     = ', ln_pnd
220         WRITE(numout,*) '         Evolutive  melt pond fraction and depth (Holland et al 2012) ln_pnd_H12 = ', ln_pnd_H12
221         WRITE(numout,*) '         Prescribed melt pond fraction and depth                      ln_pnd_CST = ', ln_pnd_CST
222         WRITE(numout,*) '            Prescribed pond fraction                                  rn_apnd    = ', rn_apnd
223         WRITE(numout,*) '            Prescribed pond depth                                     rn_hpnd    = ', rn_hpnd
224         WRITE(numout,*) '         Melt ponds affect albedo or not                              ln_pnd_alb = ', ln_pnd_alb
225      ENDIF
226      !
227      !                             !== set the choice of ice pond scheme ==!
228      ioptio = 0
229      IF( .NOT.ln_pnd ) THEN   ;   ioptio = ioptio + 1   ;   nice_pnd = np_pndNO     ;   ENDIF
230      IF( ln_pnd_CST  ) THEN   ;   ioptio = ioptio + 1   ;   nice_pnd = np_pndCST    ;   ENDIF
231      IF( ln_pnd_H12  ) THEN   ;   ioptio = ioptio + 1   ;   nice_pnd = np_pndH12    ;   ENDIF
232      IF( ioptio /= 1 )   &
233         & CALL ctl_stop( 'ice_thd_pnd_init: choose either none (ln_pnd=F) or only one pond scheme (ln_pnd_H12 or ln_pnd_CST)' )
234      !
235      SELECT CASE( nice_pnd )
236      CASE( np_pndNO )         
237         IF( ln_pnd_alb ) THEN ; ln_pnd_alb = .FALSE. ; CALL ctl_warn( 'ln_pnd_alb=false when no ponds' ) ; ENDIF
238      END SELECT
239      !
240   END SUBROUTINE ice_thd_pnd_init
241   
242#else
243   !!----------------------------------------------------------------------
244   !!   Default option          Empty module          NO SI3 sea-ice model
245   !!----------------------------------------------------------------------
246#endif 
247
248   !!======================================================================
249END MODULE icethd_pnd 
Note: See TracBrowser for help on using the repository browser.