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.
sbcssr.F90 in NEMO/branches/UKMO/NEMO_4.0.4_CO9_shelf_climate/src/OCE/SBC – NEMO

source: NEMO/branches/UKMO/NEMO_4.0.4_CO9_shelf_climate/src/OCE/SBC/sbcssr.F90 @ 15312

Last change on this file since 15312 was 15312, checked in by hadjt, 2 years ago

Haney Correction for radiative fluxes

SBC/sbc_ssr.F90 edited to add the Haney Correction to the radiative fluxes
namelist keyword added to namsbc_ssr: ln_UKMO_haney

The version is committed in restrospect, so may not include all the changes to work.
File size: 16.2 KB
Line 
1MODULE sbcssr
2   !!======================================================================
3   !!                       ***  MODULE  sbcssr  ***
4   !! Surface module :  heat and fresh water fluxes a restoring term toward observed SST/SSS
5   !!======================================================================
6   !! History :  3.0  !  2006-06  (G. Madec)  Original code
7   !!            3.2  !  2009-04  (B. Lemaire)  Introduce iom_put
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   sbc_ssr       : add to sbc a restoring term toward SST/SSS climatology
12   !!   sbc_ssr_init  : initialisation of surface restoring
13   !!----------------------------------------------------------------------
14   USE oce            ! ocean dynamics and tracers
15   USE dom_oce        ! ocean space and time domain
16   USE sbc_oce        ! surface boundary condition
17   USE phycst         ! physical constants
18   USE sbcrnf         ! surface boundary condition : runoffs
19   !
20   USE fldread        ! read input fields
21   USE in_out_manager ! I/O manager
22   USE iom            ! I/O manager
23   USE lib_mpp        ! distribued memory computing library
24   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
25   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC   sbc_ssr        ! routine called in sbcmod
31   PUBLIC   sbc_ssr_init   ! routine called in sbcmod
32   PUBLIC   sbc_ssr_alloc  ! routine called in sbcmod
33
34   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   erp   !: evaporation damping   [kg/m2/s]
35   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   qrp   !: heat flux damping        [w/m2]
36   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   coefice   !: under ice relaxation coefficient
37
38   !                                   !!* Namelist namsbc_ssr *
39   INTEGER, PUBLIC ::   nn_sstr         ! SST/SSS restoring indicator
40   INTEGER, PUBLIC ::   nn_sssr         ! SST/SSS restoring indicator
41   REAL(wp)        ::   rn_dqdt         ! restoring factor on SST and SSS
42   REAL(wp)        ::   rn_deds         ! restoring factor on SST and SSS
43   LOGICAL         ::   ln_sssr_bnd     ! flag to bound erp term
44   REAL(wp)        ::   rn_sssr_bnd     ! ABS(Max./Min.) value of erp term [mm/day]
45   INTEGER         ::   nn_sssr_ice     ! Control of restoring under ice
46   ! JT
47   LOGICAL         ::   ln_UKMO_haney   ! UKMO specific flag to calculate Haney forcing 
48   ! JT
49
50   REAL(wp) , ALLOCATABLE, DIMENSION(:) ::   buffer   ! Temporary buffer for exchange
51   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sst   ! structure of input SST (file informations, fields read)
52   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sss   ! structure of input SSS (file informations, fields read)
53
54   !!----------------------------------------------------------------------
55   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
56   !! $Id$
57   !! Software governed by the CeCILL license (see ./LICENSE)
58   !!----------------------------------------------------------------------
59CONTAINS
60
61   SUBROUTINE sbc_ssr( kt )
62      !!---------------------------------------------------------------------
63      !!                     ***  ROUTINE sbc_ssr  ***
64      !!
65      !! ** Purpose :   Add to heat and/or freshwater fluxes a damping term
66      !!                toward observed SST and/or SSS.
67      !!
68      !! ** Method  : - Read namelist namsbc_ssr
69      !!              - Read observed SST and/or SSS
70      !!              - at each nscb time step
71      !!                   add a retroaction term on qns    (nn_sstr = 1)
72      !!                   add a damping term on sfx        (nn_sssr = 1)
73      !!                   add a damping term on emp        (nn_sssr = 2)
74      !!---------------------------------------------------------------------
75      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
76      !!
77      INTEGER  ::   ji, jj   ! dummy loop indices
78      REAL(wp) ::   zerp     ! local scalar for evaporation damping
79      REAL(wp) ::   zqrp     ! local scalar for heat flux damping
80      REAL(wp) ::   zsrp     ! local scalar for unit conversion of rn_deds factor
81      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor
82      INTEGER  ::   ierror   ! return error code
83      !!
84      ! JT
85      REAL(wp) ::   sst1,sst2                      ! sea surface temperature
86      REAL(wp) ::   e_sst1, e_sst2                 ! saturation vapour pressure
87      REAL(wp) ::   qs1,qs2                        ! specific humidity
88      REAL(wp) ::   pr_tmp                         ! temporary variable for pressure
89 
90      REAL(wp), DIMENSION(jpi,jpj) ::  hny_frc1    ! Haney forcing for sensible heat, correction for latent heat   
91      REAL(wp), DIMENSION(jpi,jpj) ::  hny_frc2    ! Haney forcing for sensible heat, correction for latent heat   
92      ! JT
93      !!
94      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
95      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read
96      !!----------------------------------------------------------------------
97      !
98      IF( nn_sstr + nn_sssr /= 0 ) THEN
99         !
100         IF( nn_sstr == 1)   CALL fld_read( kt, nn_fsbc, sf_sst )   ! Read SST data and provides it at kt
101         IF( nn_sssr >= 1)   CALL fld_read( kt, nn_fsbc, sf_sss )   ! Read SSS data and provides it at kt
102         !
103         !                                         ! ========================= !
104         IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN      !    Add restoring term     !
105            !                                      ! ========================= !
106            !
107            IF( nn_sstr == 1 ) THEN                                   !* Temperature restoring term
108               ! JT
109               IF( ln_UKMO_haney ) THEN
110                  DO jj = 1, jpj
111                     DO ji = 1, jpi
112                        sst1   =  sst_m(ji,jj)
113                        sst2   =  sf_sst(1)%fnow(ji,jj,1)   
114                        e_sst1 = 10**((0.7859+0.03477*sst1)/(1.+0.00412*sst1))
115                        e_sst2 = 10**((0.7859+0.03477*sst2)/(1.+0.00412*sst2))         
116                        pr_tmp = 0.01*pressnow(ji,jj)  !pr_tmp = 1012.0
117                        qs1    = (0.62197*e_sst1)/(pr_tmp-0.378*e_sst1)
118                        qs2    = (0.62197*e_sst2)/(pr_tmp-0.378*e_sst2)
119                        hny_frc1(ji,jj) = sst1-sst2                   
120                        hny_frc2(ji,jj) = qs1-qs2                     
121                       !Might need to mask off land points.
122                        hny_frc1(ji,jj)=-hny_frc1(ji,jj)*wndm(ji,jj)*1.42
123                        hny_frc2(ji,jj)=-hny_frc2(ji,jj)*wndm(ji,jj)*4688.0
124                        ! JT Have masked out the land points
125                        qns(ji,jj)=qns(ji,jj)+(hny_frc1(ji,jj)+hny_frc2(ji,jj))*tmask(ji,jj,1)
126                        qrp(ji,jj) = 0.e0
127                     END DO
128                  END DO
129               ELSE
130              ! JT
131                  DO jj = 1, jpj
132                     DO ji = 1, jpi
133                        zqrp = rn_dqdt * ( sst_m(ji,jj) - sf_sst(1)%fnow(ji,jj,1) ) * tmask(ji,jj,1)
134                        qns(ji,jj) = qns(ji,jj) + zqrp
135                        qrp(ji,jj) = zqrp
136                     END DO
137                  END DO
138               ENDIF
139            ENDIF
140            ! JT
141            ! JT CALL iom_put( "qrp", qrp )                             ! heat flux damping
142            ! JT
143            !
144            IF( nn_sssr /= 0 .AND. nn_sssr_ice /= 1 ) THEN
145              ! use fraction of ice ( fr_i ) to adjust relaxation under ice if nn_sssr_ice .ne. 1
146              ! n.b. coefice is initialised and fixed to 1._wp if nn_sssr_ice = 1
147               DO jj = 1, jpj
148                  DO ji = 1, jpi
149                     SELECT CASE ( nn_sssr_ice )
150                       CASE ( 0 )    ;  coefice(ji,jj) = 1._wp - fr_i(ji,jj)              ! no/reduced damping under ice
151                       CASE  DEFAULT ;  coefice(ji,jj) = 1._wp + ( nn_sssr_ice - 1 ) * fr_i(ji,jj) ! reinforced damping (x nn_sssr_ice) under ice )
152                     END SELECT
153                  END DO
154               END DO
155            ENDIF
156            !
157            IF( nn_sssr == 1 ) THEN                                   !* Salinity damping term (salt flux only (sfx))
158               zsrp = rn_deds / rday                                  ! from [mm/day] to [kg/m2/s]
159               DO jj = 1, jpj
160                  DO ji = 1, jpi
161                     zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths
162                        &        *   coefice(ji,jj)            &      ! Optional control of damping under sea-ice
163                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) ) * tmask(ji,jj,1)
164                     sfx(ji,jj) = sfx(ji,jj) + zerp                 ! salt flux
165                     erp(ji,jj) = zerp / MAX( sss_m(ji,jj), 1.e-20 ) ! converted into an equivalent volume flux (diagnostic only)
166                  END DO
167               END DO
168
169               !
170            ELSEIF( nn_sssr == 2 ) THEN                               !* Salinity damping term (volume flux (emp) and associated heat flux (qns)
171               zsrp = rn_deds / rday                                  ! from [mm/day] to [kg/m2/s]
172               zerp_bnd = rn_sssr_bnd / rday                          !       -              -   
173               DO jj = 1, jpj
174                  DO ji = 1, jpi                           
175                     zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths
176                        &        *   coefice(ji,jj)            &      ! Optional control of damping under sea-ice
177                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) )   &
178                        &        / MAX(  sss_m(ji,jj), 1.e-20   ) * tmask(ji,jj,1)
179                     IF( ln_sssr_bnd )   zerp = SIGN( 1., zerp ) * MIN( zerp_bnd, ABS(zerp) )
180                     emp(ji,jj) = emp (ji,jj) + zerp
181                     qns(ji,jj) = qns(ji,jj) - zerp * rcp * sst_m(ji,jj)
182                     erp(ji,jj) = zerp
183                  END DO
184               END DO
185               ! JT CALL iom_put( "erp", erp )                             ! freshwater flux damping
186            ENDIF
187            !
188         ENDIF
189         !
190      ENDIF
191      !
192   END SUBROUTINE sbc_ssr
193
194 
195   SUBROUTINE sbc_ssr_init
196      !!---------------------------------------------------------------------
197      !!                  ***  ROUTINE sbc_ssr_init  ***
198      !!
199      !! ** Purpose :   initialisation of surface damping term
200      !!
201      !! ** Method  : - Read namelist namsbc_ssr
202      !!              - Read observed SST and/or SSS if required
203      !!---------------------------------------------------------------------
204      INTEGER  ::   ji, jj   ! dummy loop indices
205      REAL(wp) ::   zerp     ! local scalar for evaporation damping
206      REAL(wp) ::   zqrp     ! local scalar for heat flux damping
207      REAL(wp) ::   zsrp     ! local scalar for unit conversion of rn_deds factor
208      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor
209      INTEGER  ::   ierror   ! return error code
210      !!
211      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
212      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read
213      NAMELIST/namsbc_ssr/ cn_dir, nn_sstr, nn_sssr, rn_dqdt, rn_deds, sn_sst, &
214              & sn_sss, ln_sssr_bnd, rn_sssr_bnd, nn_sssr_ice, &
215              & ln_UKMO_haney    ! JT
216      INTEGER     ::  ios
217      !!----------------------------------------------------------------------
218      !
219      IF(lwp) THEN
220         WRITE(numout,*)
221         WRITE(numout,*) 'sbc_ssr : SST and/or SSS damping term '
222         WRITE(numout,*) '~~~~~~~ '
223      ENDIF
224      !
225      REWIND( numnam_ref )              ! Namelist namsbc_ssr in reference namelist :
226      READ  ( numnam_ref, namsbc_ssr, IOSTAT = ios, ERR = 901)
227901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc_ssr in reference namelist' )
228
229      REWIND( numnam_cfg )              ! Namelist namsbc_ssr in configuration namelist :
230      READ  ( numnam_cfg, namsbc_ssr, IOSTAT = ios, ERR = 902 )
231902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namsbc_ssr in configuration namelist' )
232      IF(lwm) WRITE ( numond, namsbc_ssr )
233
234      IF(lwp) THEN                 !* control print
235         WRITE(numout,*) '   Namelist namsbc_ssr :'
236         WRITE(numout,*) '      SST restoring term (Yes=1)             nn_sstr        = ', nn_sstr
237         WRITE(numout,*) '         dQ/dT (restoring magnitude on SST)     rn_dqdt     = ', rn_dqdt, ' W/m2/K'
238         WRITE(numout,*) '      SSS damping term (Yes=1, salt   flux)  nn_sssr        = ', nn_sssr
239         WRITE(numout,*) '                       (Yes=2, volume flux) '
240         WRITE(numout,*) '         dE/dS (restoring magnitude on SST)     rn_deds     = ', rn_deds, ' mm/day'
241         WRITE(numout,*) '         flag to bound erp term                 ln_sssr_bnd = ', ln_sssr_bnd
242         WRITE(numout,*) '         ABS(Max./Min.) erp threshold           rn_sssr_bnd = ', rn_sssr_bnd, ' mm/day'
243         WRITE(numout,*) '      Haney forcing                          ln_UKMO_haney = ', ln_UKMO_haney
244         WRITE(numout,*) '      Cntrl of surface restoration under ice nn_sssr_ice    = ', nn_sssr_ice
245         WRITE(numout,*) '          ( 0 = no restoration under ice)'
246         WRITE(numout,*) '          ( 1 = restoration everywhere  )'
247         WRITE(numout,*) '          (>1 = enhanced restoration under ice  )'
248      ENDIF
249      !
250      IF( nn_sstr == 1 ) THEN      !* set sf_sst structure & allocate arrays
251         !
252         ALLOCATE( sf_sst(1), STAT=ierror )
253         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst structure' )
254         ALLOCATE( sf_sst(1)%fnow(jpi,jpj,1), STAT=ierror )
255         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst now array' )
256         !
257         ! fill sf_sst with sn_sst and control print
258         CALL fld_fill( sf_sst, (/ sn_sst /), cn_dir, 'sbc_ssr', 'SST restoring term toward SST data', 'namsbc_ssr', no_print )
259         IF( sf_sst(1)%ln_tint )   ALLOCATE( sf_sst(1)%fdta(jpi,jpj,1,2), STAT=ierror )
260         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst data array' )
261         !
262      ENDIF
263      !
264      IF( nn_sssr >= 1 ) THEN      !* set sf_sss structure & allocate arrays
265         !
266         ALLOCATE( sf_sss(1), STAT=ierror )
267         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss structure' )
268         ALLOCATE( sf_sss(1)%fnow(jpi,jpj,1), STAT=ierror )
269         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss now array' )
270         !
271         ! fill sf_sss with sn_sss and control print
272         CALL fld_fill( sf_sss, (/ sn_sss /), cn_dir, 'sbc_ssr', 'SSS restoring term toward SSS data', 'namsbc_ssr', no_print )
273         IF( sf_sss(1)%ln_tint )   ALLOCATE( sf_sss(1)%fdta(jpi,jpj,1,2), STAT=ierror )
274         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss data array' )
275         !
276      ENDIF
277      !
278      coefice(:,:) = 1._wp         !  Initialise coefice to 1._wp ; will not need to be changed if nn_sssr_ice=1
279      !                            !* Initialize qrp and erp if no restoring
280      IF( nn_sstr /= 1                   )   qrp(:,:) = 0._wp
281      IF( nn_sssr /= 1 .OR. nn_sssr /= 2 )   erp(:,:) = 0._wp
282      !
283   END SUBROUTINE sbc_ssr_init
284         
285   INTEGER FUNCTION sbc_ssr_alloc()
286      !!----------------------------------------------------------------------
287      !!               ***  FUNCTION sbc_ssr_alloc  ***
288      !!----------------------------------------------------------------------
289      sbc_ssr_alloc = 0       ! set to zero if no array to be allocated
290      IF( .NOT. ALLOCATED( erp ) ) THEN
291         ALLOCATE( qrp(jpi,jpj), erp(jpi,jpj), coefice(jpi,jpj), STAT= sbc_ssr_alloc )
292         !
293         IF( lk_mpp                  )   CALL mpp_sum ( 'sbcssr', sbc_ssr_alloc )
294         IF( sbc_ssr_alloc /= 0 )   CALL ctl_warn('sbc_ssr_alloc: failed to allocate arrays.')
295         !
296      ENDIF
297   END FUNCTION
298     
299   !!======================================================================
300END MODULE sbcssr
Note: See TracBrowser for help on using the repository browser.