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 branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icethd_pnd.F90 @ 9119

Last change on this file since 9119 was 9019, checked in by timgraham, 6 years ago

Merge of dev_CNRS_2017 into branch
File size: 12.9 KB
Line 
1MODULE icethd_pnd 
2   !!======================================================================
3   !!                     ***  MODULE  icethd_pnd   ***
4   !!   Melt ponds
5   !!======================================================================
6   !! history :       ! Original code by Daniela Flocco and Adrian Turner
7   !!            1.0  ! 2012    (O. Lecomte) Adaptation for scientific tests (NEMO3.1)
8   !!            2.0  ! 2017    (M. Vancoppenolle, O. Lecomte, C. Rousset) Implementation in NEMO4
9   !!----------------------------------------------------------------------
10#if defined key_lim3
11   !!----------------------------------------------------------------------
12   !!   'key_lim3' :                                     ESIM 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   !! * Substitutions
41#  include "vectopt_loop_substitute.h90"
42   !!----------------------------------------------------------------------
43   !! NEMO/ICE 4.0 , NEMO Consortium (2017)
44   !! $Id: icethd_pnd.F90 8420 2017-10-05 13:07:10Z clem $
45   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
46   !!----------------------------------------------------------------------
47CONTAINS
48
49   SUBROUTINE ice_thd_pnd
50      !!-------------------------------------------------------------------
51      !!               ***  ROUTINE ice_thd_pnd   ***
52      !!               
53      !! ** Purpose :   change melt pond fraction
54      !!               
55      !! ** Method  :   brut force
56      !!
57      !! ** Action  : -
58      !!              -
59      !!-------------------------------------------------------------------
60
61      SELECT CASE ( nice_pnd )
62
63      CASE (np_pndCST)
64         !                             !-------------------------------!
65         CALL pnd_CST                  ! Constant melt ponds           !
66         !                             !-------------------------------!
67      CASE (np_pndH12)
68         !                             !-------------------------------!
69         CALL pnd_H12                  ! Holland et al 2012 melt ponds !
70         !                             !-------------------------------!
71      END SELECT
72
73   END SUBROUTINE ice_thd_pnd 
74
75   SUBROUTINE pnd_CST 
76      !!-------------------------------------------------------------------
77      !!                ***  ROUTINE pnd_CST  ***
78      !!
79      !! ** Purpose    : Compute melt pond evolution
80      !!
81      !! ** Method     : Melt pond fraction and thickness are prescribed
82      !!                 to non-zero values when t_su = 0C
83      !!
84      !! ** Tunable parameters : pond fraction (rn_apnd), pond depth (rn_hpnd)
85      !!               
86      !! ** Note       : Coupling with such melt ponds is only radiative
87      !!                 Advection, ridging, rafting... are bypassed
88      !!
89      !! ** References : Bush, G.W., and Trump, D.J. (2017)
90      !!   
91      !!-------------------------------------------------------------------
92      INTEGER  ::   ji        ! loop indices
93      !!-------------------------------------------------------------------
94      DO ji = 1, npti
95         
96         IF( a_i_1d(ji) > 0._wp .AND. t_su_1d(ji) >= rt0 ) THEN
97            a_ip_frac_1d(ji) = rn_apnd
98            h_ip_1d(ji)      = rn_hpnd   
99            a_ip_1d(ji)      = a_ip_frac_1d(ji) * a_i_1d(ji)
100         ELSE
101            a_ip_frac_1d(ji) = 0._wp
102            h_ip_1d(ji)      = 0._wp   
103            a_ip_1d(ji)      = 0._wp
104         ENDIF
105         
106      END DO
107     
108   END SUBROUTINE pnd_CST
109
110   SUBROUTINE pnd_H12
111      !!-------------------------------------------------------------------
112      !!                ***  ROUTINE pnd_H12  ***
113      !!
114      !! ** Purpose    : Compute melt pond evolution
115      !!
116      !! ** Method     : Empirical method. A fraction of meltwater is accumulated in ponds
117      !!                 and sent to ocean when surface is freezing
118      !!
119      !!                 pond growth:      Vp = Vp + dVmelt
120      !!                    with dVmelt = R/rhow * ( rhoi*dh_i + rhos*dh_s ) * a_i
121      !!                 pond contraction: Vp = Vp * exp(0.01*MAX(Tp-Tsu,0)/Tp)
122      !!                    with Tp = -2degC
123      !! 
124      !! ** Tunable parameters : (no real expertise yet, ideas?)
125      !!
126      !! ** Note       : Stolen from CICE for quick test of the melt pond
127      !!                 radiation and freshwater interfaces
128      !!                 Coupling can be radiative AND freshwater
129      !!                 Advection, ridging, rafting are called
130      !!
131      !! ** References : Holland, M. M. et al (J Clim 2012)
132      !!   
133      !!-------------------------------------------------------------------
134      REAL(wp), PARAMETER ::   zrmin       = 0.15_wp  ! minimum fraction of available meltwater retained for melt ponding
135      REAL(wp), PARAMETER ::   zrmax       = 0.70_wp  ! maximum   ''           ''       ''        ''            ''
136      REAL(wp), PARAMETER ::   zpnd_aspect = 0.8_wp   ! pond aspect ratio
137      REAL(wp), PARAMETER ::   zTp         = -2._wp   ! reference temperature
138
139      REAL(wp) ::   zfr_mlt          ! fraction of available meltwater retained for melt ponding
140      REAL(wp) ::   zdv_mlt          ! available meltwater for melt ponding
141      REAL(wp) ::   z1_Tp            ! inverse reference temperature
142      REAL(wp) ::   z1_rhofw         ! inverse freshwater density
143      REAL(wp) ::   z1_zpnd_aspect   ! inverse pond aspect ratio
144      REAL(wp) ::   zfac, zdum
145
146      INTEGER  ::   ji   ! loop indices
147      !!-------------------------------------------------------------------
148      z1_rhofw       = 1. / rhofw 
149      z1_zpnd_aspect = 1. / zpnd_aspect
150      z1_Tp          = 1._wp / zTp 
151
152      DO ji = 1, npti
153         !                                                        !--------------------------------!
154         IF( h_i_1d(ji) < rn_himin) THEN                          ! Case ice thickness < rn_himin  !
155            !                                                     !--------------------------------!
156            !--- Remove ponds on thin ice
157            a_ip_1d(ji)      = 0._wp
158            a_ip_frac_1d(ji) = 0._wp
159            h_ip_1d(ji)      = 0._wp
160            !                                                     !--------------------------------!
161         ELSE                                                     ! Case ice thickness >= rn_himin !
162            !                                                     !--------------------------------!
163            v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji)   ! record pond volume at previous time step
164
165            ! available meltwater for melt ponding [m, >0] and fraction
166            zdv_mlt = -( dh_i_surf(ji)*rhoic + dh_s_mlt(ji)*rhosn ) * z1_rhofw * a_i_1d(ji)
167            zfr_mlt = zrmin + ( zrmax - zrmin ) * a_i_1d(ji)  ! from CICE doc
168            !zfr_mlt = zrmin + zrmax * a_i_1d(ji)             ! from Holland paper
169
170            !--- Pond gowth ---!
171            ! v_ip should never be negative, otherwise code crashes
172            ! MV: as far as I saw, UM5 can create very small negative v_ip values (not Prather)
173            v_ip_1d(ji) = MAX( 0._wp, v_ip_1d(ji) + zfr_mlt * zdv_mlt )
174
175            ! melt pond mass flux (<0)
176            IF( ln_pnd_fwb .AND. zdv_mlt > 0._wp ) THEN
177               zfac = zfr_mlt * zdv_mlt * rhofw * r1_rdtice
178               wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac
179
180               ! adjust ice/snow melting flux to balance melt pond flux (>0)
181               zdum = zfac / ( wfx_snw_sum_1d(ji) + wfx_sum_1d(ji) )
182               wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) * (1._wp + zdum)
183               wfx_sum_1d(ji)     = wfx_sum_1d(ji)     * (1._wp + zdum)
184            ENDIF
185
186            !--- Pond contraction (due to refreezing) ---!
187            v_ip_1d(ji) = v_ip_1d(ji) * EXP( 0.01_wp * MAX( zTp+rt0 - t_su_1d(ji), 0._wp ) * z1_Tp )
188
189            ! Set new pond area and depth assuming linear relation between h_ip and a_ip_frac
190            !    h_ip = zpnd_aspect * a_ip_frac = zpnd_aspect * a_ip/a_i
191            a_ip_1d(ji)      = SQRT( v_ip_1d(ji) * z1_zpnd_aspect * a_i_1d(ji) )
192            a_ip_frac_1d(ji) = a_ip_1d(ji) / a_i_1d(ji)
193            h_ip_1d(ji)      = zpnd_aspect * a_ip_frac_1d(ji)
194
195         ENDIF
196      END DO
197
198   END SUBROUTINE pnd_H12
199
200   SUBROUTINE ice_thd_pnd_init 
201      !!-------------------------------------------------------------------
202      !!                  ***  ROUTINE ice_thd_pnd_init   ***
203      !!
204      !! ** Purpose : Physical constants and parameters linked to melt ponds
205      !!              over sea ice
206      !!
207      !! ** Method  :  Read the namthd_pnd  namelist and check the melt pond 
208      !!               parameter values called at the first timestep (nit000)
209      !!
210      !! ** input   :   Namelist namthd_pnd 
211      !!-------------------------------------------------------------------
212      INTEGER  ::   ios, ioptio                 ! Local integer output status for namelist read
213      NAMELIST/namthd_pnd/  ln_pnd_H12, ln_pnd_fwb, ln_pnd_CST, rn_apnd, rn_hpnd, ln_pnd_alb
214      !!-------------------------------------------------------------------
215
216      REWIND( numnam_ice_ref )              ! Namelist namthd_pnd  in reference namelist : Melt Ponds 
217      READ  ( numnam_ice_ref, namthd_pnd, IOSTAT = ios, ERR = 901)
218901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namthd_pnd  in reference namelist', lwp )
219
220      REWIND( numnam_ice_cfg )              ! Namelist namthd_pnd  in configuration namelist : Melt Ponds
221      READ  ( numnam_ice_cfg, namthd_pnd, IOSTAT = ios, ERR = 902 )
222902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namthd_pnd in configuration namelist', lwp )
223      IF(lwm) WRITE ( numoni, namthd_pnd )
224     
225      IF(lwp) THEN                        ! control print
226         WRITE(numout,*)
227         WRITE(numout,*) 'ice_thd_pnd_init: ice parameters for melt ponds'
228         WRITE(numout,*) '~~~~~~~~~~~~~~~~'
229         WRITE(numout,*) '   Namelist namicethd_pnd:'
230         WRITE(numout,*) '      Evolutive melt pond fraction and depth (Holland et al 2012)  ln_pnd_H12 = ', ln_pnd_H12
231         WRITE(numout,*) '         Melt ponds store fresh water or not                       ln_pnd_fwb = ', ln_pnd_fwb
232         WRITE(numout,*) '      Prescribed melt pond fraction and depth                      ln_pnd_Cst = ', ln_pnd_CST
233         WRITE(numout,*) '         Prescribed pond fraction                                  rn_apnd    = ', rn_apnd
234         WRITE(numout,*) '         Prescribed pond depth                                     rn_hpnd    = ', rn_hpnd
235         WRITE(numout,*) '      Melt ponds affect albedo or not                              ln_pnd_alb = ', ln_pnd_alb
236      ENDIF
237      !
238      !                             !== set the choice of ice pond scheme ==!
239      ioptio = 0
240                                                            nice_pnd = np_pndNO
241      IF( ln_pnd_CST ) THEN   ;   ioptio = ioptio + 1   ;   nice_pnd = np_pndCST    ;   ENDIF
242      IF( ln_pnd_H12 ) THEN   ;   ioptio = ioptio + 1   ;   nice_pnd = np_pndH12    ;   ENDIF
243      IF( ioptio > 1 )   CALL ctl_stop( 'ice_thd_pnd_init: choose one and only one pond scheme (ln_pnd_H12 or ln_pnd_CST)' )
244
245      SELECT CASE( nice_pnd )
246      CASE( np_pndNO )         
247         IF(ln_pnd_fwb) THEN ; ln_pnd_fwb = .FALSE. ; CALL ctl_warn( 'ln_pnd_fwb=false when no ponds' ) ; ENDIF
248         IF(ln_pnd_alb) THEN ; ln_pnd_alb = .FALSE. ; CALL ctl_warn( 'ln_pnd_alb=false when no ponds' ) ; ENDIF
249      CASE( np_pndCST)
250         IF(ln_pnd_fwb) THEN ; ln_pnd_fwb = .FALSE. ; CALL ctl_warn( 'ln_pnd_fwb=false when ln_pnd_CST=true' ) ; ENDIF
251      END SELECT
252      !
253   END SUBROUTINE ice_thd_pnd_init
254   
255#else
256   !!----------------------------------------------------------------------
257   !!   Default option          Empty module          NO ESIM sea-ice model
258   !!----------------------------------------------------------------------
259#endif 
260
261   !!======================================================================
262END MODULE icethd_pnd 
Note: See TracBrowser for help on using the repository browser.