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.
albedo.F90 in branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

source: branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/OPA_SRC/SBC/albedo.F90 @ 13191

Last change on this file since 13191 was 13191, checked in by jwhile, 4 years ago

Updates for 1d runnig

File size: 15.3 KB
RevLine 
[152]1MODULE albedo
2   !!======================================================================
3   !!                       ***  MODULE  albedo  ***
4   !! Ocean forcing:  bulk thermohaline forcing of the ocean (or ice)
5   !!=====================================================================
[1601]6   !! History :  8.0  ! 2001-04  (LIM 1.0)
7   !!   NEMO     1.0  ! 2003-07  (C. Ethe, G. Madec)  Optimization (old name:shine)
8   !!             -   ! 2004-11  (C. Talandier)  add albedo_init
9   !!             -   ! 2001-06  (M. Vancoppenolle) LIM 3.0
10   !!             -   ! 2006-08  (G. Madec)  cleaning for surface module
[6498]11   !!            3.6  ! 2016-01  (C. Rousset) new parameterization for sea ice albedo
[152]12   !!----------------------------------------------------------------------
[1601]13
14   !!----------------------------------------------------------------------
[3625]15   !!   albedo_ice    : albedo for   ice (clear and overcast skies)
16   !!   albedo_oce    : albedo for ocean (clear and overcast skies)
17   !!   albedo_init   : initialisation of albedo computation
[152]18   !!----------------------------------------------------------------------
[3625]19   USE phycst         ! physical constants
20   USE in_out_manager ! I/O manager
21   USE lib_mpp        ! MPP library
22   USE wrk_nemo       ! work arrays
[13191]23   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
[11442]24   USE stopack
[152]25
26   IMPLICIT NONE
27   PRIVATE
28
[1601]29   PUBLIC   albedo_ice   ! routine called sbcice_lim.F90
30   PUBLIC   albedo_oce   ! routine called by ???
[165]31
[888]32   INTEGER  ::   albd_init = 0      !: control flag for initialization
[13191]33
[6498]34   REAL(wp) ::   rmue     = 0.40    !  cosine of local solar altitude
35   REAL(wp) ::   ralb_oce = 0.066   ! ocean or lead albedo (Pegau and Paulson, Ann. Glac. 2001)
36   REAL(wp) ::   c1       = 0.05    ! snow thickness (only for nn_ice_alb=0)
37   REAL(wp) ::   c2       = 0.10    !  "        "
38   REAL(wp) ::   rcloud   = 0.06    ! cloud effect on albedo (only-for nn_ice_alb=0)
[13191]39
[4147]40   !                             !!* namelist namsbc_alb
[6498]41   INTEGER  ::   nn_ice_alb
42   REAL(wp) ::   rn_albice
[152]43
44   !!----------------------------------------------------------------------
[2528]45   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
[888]46   !! $Id$
[2715]47   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
[152]48   !!----------------------------------------------------------------------
49CONTAINS
50
[888]51   SUBROUTINE albedo_ice( pt_ice, ph_ice, ph_snw, pa_ice_cs, pa_ice_os )
[152]52      !!----------------------------------------------------------------------
[888]53      !!               ***  ROUTINE albedo_ice  ***
[13191]54      !!
55      !! ** Purpose :   Computation of the albedo of the snow/ice system
56      !!
[6498]57      !! ** Method  :   Two schemes are available (from namelist parameter nn_ice_alb)
58      !!                  0: the scheme is that of Shine & Henderson-Sellers (JGR 1985) for clear-skies
59      !!                  1: the scheme is "home made" (for cloudy skies) and based on Brandt et al. (J. Climate 2005)
60      !!                                                                           and Grenfell & Perovich (JGR 2004)
61      !!                Description of scheme 1:
62      !!                  1) Albedo dependency on ice thickness follows the findings from Brandt et al (2005)
63      !!                     which are an update of Allison et al. (JGR 1993) ; Brandt et al. 1999
64      !!                     0-5cm  : linear function of ice thickness
65      !!                     5-150cm: log    function of ice thickness
66      !!                     > 150cm: constant
67      !!                  2) Albedo dependency on snow thickness follows the findings from Grenfell & Perovich (2004)
68      !!                     i.e. it increases as -EXP(-snw_thick/0.02) during freezing and -EXP(-snw_thick/0.03) during melting
69      !!                  3) Albedo dependency on clouds is speculated from measurements of Grenfell and Perovich (2004)
70      !!                     i.e. cloudy-clear albedo depend on cloudy albedo following a 2d order polynomial law
71      !!                  4) The needed 4 parameters are: dry and melting snow, freezing ice and bare puddled ice
[152]72      !!
[6498]73      !! ** Note    :   The parameterization from Shine & Henderson-Sellers presents several misconstructions:
74      !!                  1) ice albedo when ice thick. tends to 0 is different than ocean albedo
[13191]75      !!                  2) for small ice thick. covered with some snow (<3cm?), albedo is larger
[6498]76      !!                     under melting conditions than under freezing conditions
[13191]77      !!                  3) the evolution of ice albedo as a function of ice thickness shows
[6498]78      !!                     3 sharp inflexion points (at 5cm, 100cm and 150cm) that look highly unrealistic
79      !!
80      !! References :   Shine & Henderson-Sellers 1985, JGR, 90(D1), 2243-2250.
81      !!                Brandt et al. 2005, J. Climate, vol 18
[13191]82      !!                Grenfell & Perovich 2004, JGR, vol 109
[888]83      !!----------------------------------------------------------------------
[1463]84      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   pt_ice      !  ice surface temperature (Kelvin)
[888]85      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_ice      !  sea-ice thickness
86      REAL(wp), INTENT(in   ), DIMENSION(:,:,:) ::   ph_snw      !  snow thickness
87      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   pa_ice_cs   !  albedo of ice under clear    sky
88      REAL(wp), INTENT(  out), DIMENSION(:,:,:) ::   pa_ice_os   !  albedo of ice under overcast sky
[719]89      !!
[6498]90      INTEGER  ::   ji, jj, jl         ! dummy loop indices
91      INTEGER  ::   ijpl               ! number of ice categories (3rd dim of ice input arrays)
92      REAL(wp)            ::   ralb_im, ralb_sf, ralb_sm, ralb_if
93      REAL(wp)            ::   zswitch, z1_c1, z1_c2
94      REAL(wp)                            ::   zalb_sm, zalb_sf, zalb_st ! albedo of snow melting, freezing, total
95      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zalb, zalb_it             ! intermediate variable & albedo of ice (snow free)
[152]96      !!---------------------------------------------------------------------
[6498]97
98      ijpl = SIZE( pt_ice, 3 )                     ! number of ice categories
[13191]99
[6498]100      CALL wrk_alloc( jpi,jpj,ijpl, zalb, zalb_it )
[165]101
[13191]102      IF( albd_init == 0 )   CALL albedo_init      ! initialization
[888]103
[13191]104
[6498]105      SELECT CASE ( nn_ice_alb )
[2715]106
[6498]107      !------------------------------------------
108      !  Shine and Henderson-Sellers (1985)
109      !------------------------------------------
110      CASE( 0 )
[13191]111
[6498]112         ralb_sf = 0.80       ! dry snow
113         ralb_sm = 0.65       ! melting snow
114         ralb_if = 0.72       ! bare frozen ice
[13191]115         ralb_im = rn_albice  ! bare puddled ice
116
[6498]117         !  Computation of ice albedo (free of snow)
118         WHERE     ( ph_snw == 0._wp .AND. pt_ice >= rt0_ice )   ;   zalb(:,:,:) = ralb_im
119         ELSE WHERE                                              ;   zalb(:,:,:) = ralb_if
120         END  WHERE
[13191]121
[6498]122         WHERE     ( 1.5  < ph_ice                     )  ;  zalb_it = zalb
123         ELSE WHERE( 1.0  < ph_ice .AND. ph_ice <= 1.5 )  ;  zalb_it = 0.472  + 2.0 * ( zalb - 0.472 ) * ( ph_ice - 1.0 )
124         ELSE WHERE( 0.05 < ph_ice .AND. ph_ice <= 1.0 )  ;  zalb_it = 0.2467 + 0.7049 * ph_ice              &
125            &                                                                 - 0.8608 * ph_ice * ph_ice     &
126            &                                                                 + 0.3812 * ph_ice * ph_ice * ph_ice
127         ELSE WHERE                                       ;  zalb_it = 0.1    + 3.6    * ph_ice
128         END WHERE
[13191]129
[6498]130         DO jl = 1, ijpl
131            DO jj = 1, jpj
132               DO ji = 1, jpi
133                  ! freezing snow
134                  ! no effect of underlying ice layer IF snow thickness > c1. Albedo does not depend on snow thick if > c2
[13191]135                  !                                        !  freezing snow
[6498]136                  zswitch   = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ( ph_snw(ji,jj,jl) - c1 ) ) )
137                  zalb_sf   = ( 1._wp - zswitch ) * (  zalb_it(ji,jj,jl)  &
138                     &                           + ph_snw(ji,jj,jl) * ( ralb_sf - zalb_it(ji,jj,jl) ) / c1  )   &
[13191]139                     &        +         zswitch   * ralb_sf
[888]140
[6498]141                  ! melting snow
142                  ! no effect of underlying ice layer. Albedo does not depend on snow thick IF > c2
143                  zswitch   = MAX( 0._wp , SIGN( 1._wp , ph_snw(ji,jj,jl) - c2 ) )
144                  zalb_sm = ( 1._wp - zswitch ) * ( ralb_im + ph_snw(ji,jj,jl) * ( ralb_sm - ralb_im ) / c2 )   &
[13191]145                      &     +         zswitch   *   ralb_sm
[6498]146                  !
147                  ! snow albedo
[13191]148                  zswitch  =  MAX( 0._wp , SIGN( 1._wp , pt_ice(ji,jj,jl) - rt0_snow ) )
[6498]149                  zalb_st  =  zswitch * zalb_sm + ( 1._wp - zswitch ) * zalb_sf
[13191]150
[6498]151                  ! Ice/snow albedo
152                  zswitch   = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ph_snw(ji,jj,jl) ) )
153                  pa_ice_cs(ji,jj,jl) =  zswitch * zalb_st + ( 1._wp - zswitch ) * zalb_it(ji,jj,jl)
154                  !
155               END DO
[833]156            END DO
[13191]157
158#if defined key_traldf_c2d || key_traldf_c3d
[11442]159            IF ( ln_stopack .AND. nn_spp_icealb > 0 ) &
160               & CALL spp_gen( 1, pa_ice_cs(:,:,jl), nn_spp_icealb, rn_icealb_sd, jk_spp_alb, jl )
[13191]161#else
162            IF ( ln_stopack .AND. nn_spp_icealb > 0 ) &
163               & CALL ctl_stop( 'albedo_ice: parameter perturbation will only work with '// &
164                                'key_traldf_c2d or key_traldf_c3d')
165#endif
[833]166         END DO
[6498]167
168         pa_ice_os(:,:,:) = pa_ice_cs(:,:,:) + rcloud       ! Oberhuber correction for overcast sky
169
170      !------------------------------------------
171      !  New parameterization (2016)
172      !------------------------------------------
[13191]173      CASE( 1 )
[6498]174
175         ralb_im = rn_albice  ! bare puddled ice
176! compilation of values from literature
177         ralb_sf = 0.85      ! dry snow
178         ralb_sm = 0.75      ! melting snow
179         ralb_if = 0.60      ! bare frozen ice
180! Perovich et al 2002 (Sheba) => the only dataset for which all types of ice/snow were retrieved
181!         ralb_sf = 0.85       ! dry snow
182!         ralb_sm = 0.72       ! melting snow
183!         ralb_if = 0.65       ! bare frozen ice
184! Brandt et al 2005 (East Antarctica)
185!         ralb_sf = 0.87      ! dry snow
186!         ralb_sm = 0.82      ! melting snow
187!         ralb_if = 0.54      ! bare frozen ice
[13191]188!
[6498]189         !  Computation of ice albedo (free of snow)
[13191]190         z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) )
[6498]191         z1_c2 = 1. / 0.05
192         WHERE     ( ph_snw == 0._wp .AND. pt_ice >= rt0_ice )   ;   zalb = ralb_im
193         ELSE WHERE                                              ;   zalb = ralb_if
194         END  WHERE
[13191]195
[6498]196         WHERE     ( 1.5  < ph_ice                     )  ;  zalb_it = zalb
197         ELSE WHERE( 0.05 < ph_ice .AND. ph_ice <= 1.5 )  ;  zalb_it = zalb     + ( 0.18 - zalb     ) * z1_c1 *  &
198            &                                                                     ( LOG(1.5) - LOG(ph_ice) )
199         ELSE WHERE                                       ;  zalb_it = ralb_oce + ( 0.18 - ralb_oce ) * z1_c2 * ph_ice
200         END WHERE
201
202         z1_c1 = 1. / 0.02
203         z1_c2 = 1. / 0.03
204         !  Computation of the snow/ice albedo
205         DO jl = 1, ijpl
206            DO jj = 1, jpj
207               DO ji = 1, jpi
208                  zalb_sf = ralb_sf - ( ralb_sf - zalb_it(ji,jj,jl)) * EXP( - ph_snw(ji,jj,jl) * z1_c1 );
209                  zalb_sm = ralb_sm - ( ralb_sm - zalb_it(ji,jj,jl)) * EXP( - ph_snw(ji,jj,jl) * z1_c2 );
210
211                   ! snow albedo
[13191]212                  zswitch = MAX( 0._wp , SIGN( 1._wp , pt_ice(ji,jj,jl) - rt0_snow ) )
[6498]213                  zalb_st = zswitch * zalb_sm + ( 1._wp - zswitch ) * zalb_sf
214
[13191]215                  ! Ice/snow albedo
[6498]216                  zswitch             = MAX( 0._wp , SIGN( 1._wp , - ph_snw(ji,jj,jl) ) )
217                  pa_ice_os(ji,jj,jl) = ( 1._wp - zswitch ) * zalb_st + zswitch *  zalb_it(ji,jj,jl)
218
219              END DO
220            END DO
[13191]221
222#if defined key_traldf_c2d || key_traldf_c3d
[11442]223            IF ( ln_stopack .AND. nn_spp_icealb > 0 ) &
224               & CALL spp_gen( 1, pa_ice_os(:,:,jl), nn_spp_icealb, rn_icealb_sd, jk_spp_alb, jl )
[13191]225#else
226            IF ( ln_stopack .AND. nn_spp_icealb > 0 ) &
227               & CALL ctl_stop( 'albedo_ice: parameter perturbation will only work with '// &
228                                'key_traldf_c2d or key_traldf_c3d')
229#endif
[6498]230         END DO
231         ! Effect of the clouds (2d order polynomial)
[13191]232         pa_ice_cs = pa_ice_os - ( - 0.1010 * pa_ice_os * pa_ice_os + 0.1933 * pa_ice_os - 0.0148 );
[6498]233
234      END SELECT
[13191]235
[6498]236      CALL wrk_dealloc( jpi,jpj,ijpl, zalb, zalb_it )
[888]237      !
238   END SUBROUTINE albedo_ice
[833]239
240
[888]241   SUBROUTINE albedo_oce( pa_oce_os , pa_oce_cs )
[152]242      !!----------------------------------------------------------------------
[888]243      !!               ***  ROUTINE albedo_oce  ***
[13191]244      !!
[888]245      !! ** Purpose :   Computation of the albedo of the ocean
246      !!----------------------------------------------------------------------
[2715]247      REAL(wp), DIMENSION(:,:), INTENT(out) ::   pa_oce_os   !  albedo of ocean under overcast sky
248      REAL(wp), DIMENSION(:,:), INTENT(out) ::   pa_oce_cs   !  albedo of ocean under clear sky
[152]249      !!
[13191]250      REAL(wp) :: zcoef
[152]251      !!----------------------------------------------------------------------
[888]252      !
[6498]253      zcoef = 0.05 / ( 1.1 * rmue**1.4 + 0.15 )   ! Parameterization of Briegled and Ramanathan, 1982
[13191]254      pa_oce_cs(:,:) = zcoef
[6498]255      pa_oce_os(:,:) = 0.06                       ! Parameterization of Kondratyev, 1969 and Payne, 1972
[888]256      !
257   END SUBROUTINE albedo_oce
[152]258
259
[165]260   SUBROUTINE albedo_init
261      !!----------------------------------------------------------------------
262      !!                 ***  ROUTINE albedo_init  ***
263      !!
264      !! ** Purpose :   initializations for the albedo parameters
265      !!
[1601]266      !! ** Method  :   Read the namelist namsbc_alb
[165]267      !!----------------------------------------------------------------------
[4147]268      INTEGER  ::   ios                 ! Local integer output status for namelist read
[13191]269      NAMELIST/namsbc_alb/ nn_ice_alb, rn_albice
[719]270      !!----------------------------------------------------------------------
[1601]271      !
272      albd_init = 1                     ! indicate that the initialization has been done
273      !
[4147]274      REWIND( numnam_ref )              ! Namelist namsbc_alb in reference namelist : Albedo parameters
275      READ  ( numnam_ref, namsbc_alb, IOSTAT = ios, ERR = 901)
276901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_alb in reference namelist', lwp )
277
278      REWIND( numnam_cfg )              ! Namelist namsbc_alb in configuration namelist : Albedo parameters
279      READ  ( numnam_cfg, namsbc_alb, IOSTAT = ios, ERR = 902 )
280902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_alb in configuration namelist', lwp )
[4624]281      IF(lwm) WRITE ( numond, namsbc_alb )
[1601]282      !
283      IF(lwp) THEN                      ! Control print
[165]284         WRITE(numout,*)
[1601]285         WRITE(numout,*) 'albedo : set albedo parameters'
286         WRITE(numout,*) '~~~~~~~'
287         WRITE(numout,*) '   Namelist namsbc_alb : albedo '
[6498]288         WRITE(numout,*) '      choose the albedo parameterization                  nn_ice_alb = ', nn_ice_alb
289         WRITE(numout,*) '      albedo of bare puddled ice                          rn_albice  = ', rn_albice
[165]290      ENDIF
[888]291      !
292   END SUBROUTINE albedo_init
[719]293
[152]294   !!======================================================================
295END MODULE albedo
Note: See TracBrowser for help on using the repository browser.