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 branches/UKMO/CO6_KD490/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

source: branches/UKMO/CO6_KD490/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssr.F90 @ 6610

Last change on this file since 6610 was 6610, checked in by deazer, 8 years ago

Adding Direct MO Fluxes to CO6
File size: 14.3 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 iom            ! I/O manager
22   USE in_out_manager ! I/O manager
23   USE lib_mpp        ! distribued memory computing library
24   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
25   USE timing         ! Timing
26   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC   sbc_ssr        ! routine called in sbcmod
32   PUBLIC   sbc_ssr_init   ! 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
37   !                                   !!* Namelist namsbc_ssr *
38   INTEGER, PUBLIC ::   nn_sstr         ! SST/SSS restoring indicator
39   INTEGER, PUBLIC ::   nn_sssr         ! SST/SSS restoring indicator
40   REAL(wp)        ::   rn_dqdt         ! restoring factor on SST and SSS
41   REAL(wp)        ::   rn_deds         ! restoring factor on SST and SSS
42   LOGICAL         ::   ln_sssr_bnd     ! flag to bound erp term
43   REAL(wp)        ::   rn_sssr_bnd     ! ABS(Max./Min.) value of erp term [mm/day]
44   LOGICAL         ::   ln_UKMO_haney   ! UKMO specific flag to calculate Haney forcing 
45
46   REAL(wp) , ALLOCATABLE, DIMENSION(:) ::   buffer   ! Temporary buffer for exchange
47   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sst   ! structure of input SST (file informations, fields read)
48   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sss   ! structure of input SSS (file informations, fields read)
49
50   !! * Substitutions
51#  include "domzgr_substitute.h90"
52   !!----------------------------------------------------------------------
53   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
54   !! $Id$
55   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
56   !!----------------------------------------------------------------------
57CONTAINS
58
59   SUBROUTINE sbc_ssr( kt )
60      !!---------------------------------------------------------------------
61      !!                     ***  ROUTINE sbc_ssr  ***
62      !!
63      !! ** Purpose :   Add to heat and/or freshwater fluxes a damping term
64      !!                toward observed SST and/or SSS.
65      !!
66      !! ** Method  : - Read namelist namsbc_ssr
67      !!              - Read observed SST and/or SSS
68      !!              - at each nscb time step
69      !!                   add a retroaction term on qns    (nn_sstr = 1)
70      !!                   add a damping term on sfx        (nn_sssr = 1)
71      !!                   add a damping term on emp        (nn_sssr = 2)
72      !!---------------------------------------------------------------------
73      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
74      !!
75      INTEGER  ::   ji, jj   ! dummy loop indices
76      REAL(wp) ::   zerp     ! local scalar for evaporation damping
77      REAL(wp) ::   zqrp     ! local scalar for heat flux damping
78      REAL(wp) ::   zsrp     ! local scalar for unit conversion of rn_deds factor
79      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor
80      INTEGER  ::   ierror   ! return error code
81      !!
82      REAL(wp) ::   sst1,sst2                      ! sea surface temperature
83      REAL(wp) ::   e_sst1, e_sst2                 ! saturation vapour pressure
84      REAL(wp) ::   qs1,qs2                        ! specific humidity
85      REAL(wp) ::   pr_tmp                         ! temporary variable for pressure
86 
87      REAL(wp), DIMENSION(jpi,jpj) ::  hny_frc1    ! Haney forcing for sensible heat, correction for latent heat   
88      REAL(wp), DIMENSION(jpi,jpj) ::  hny_frc2    ! Haney forcing for sensible heat, correction for latent heat   
89      !!
90      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
91      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read
92      !!----------------------------------------------------------------------
93      !
94      IF( nn_timing == 1 )  CALL timing_start('sbc_ssr')
95      !
96      IF( nn_sstr + nn_sssr /= 0 ) THEN
97         !
98         IF( nn_sstr == 1)   CALL fld_read( kt, nn_fsbc, sf_sst )   ! Read SST data and provides it at kt
99         IF( nn_sssr >= 1)   CALL fld_read( kt, nn_fsbc, sf_sss )   ! Read SSS data and provides it at kt
100         !
101         !                                         ! ========================= !
102         IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN      !    Add restoring term     !
103            !                                      ! ========================= !
104            !
105            IF( nn_sstr == 1 ) THEN                                   !* Temperature restoring term
106                  IF( ln_UKMO_haney ) THEN
107                     DO jj = 1, jpj
108                        DO ji = 1, jpi
109                           sst1   =  sst_m(ji,jj)
110                           sst2   =  sf_sst(1)%fnow(ji,jj,1)   
111                           e_sst1 = 10**((0.7859+0.03477*sst1)/(1.+0.00412*sst1))
112                           e_sst2 = 10**((0.7859+0.03477*sst2)/(1.+0.00412*sst2))         
113                           pr_tmp = 0.01*pressnow(ji,jj)  !pr_tmp = 1012.0
114                           qs1    = (0.62197*e_sst1)/(pr_tmp-0.378*e_sst1)
115                           qs2    = (0.62197*e_sst2)/(pr_tmp-0.378*e_sst2)
116                           hny_frc1(ji,jj) = sst1-sst2                   
117                           hny_frc2(ji,jj) = qs1-qs2                     
118                          !Might need to mask off land points.
119                           hny_frc1(ji,jj)=-hny_frc1(ji,jj)*wndm(ji,jj)*1.42
120                           hny_frc2(ji,jj)=-hny_frc2(ji,jj)*wndm(ji,jj)*4688.0
121                           qns(ji,jj)=qns(ji,jj)+hny_frc1(ji,jj)+hny_frc2(ji,jj)   
122                           qrp(ji,jj) = 0.e0
123                        END DO
124                     END DO
125                  ELSE
126                     DO jj = 1, jpj
127                        DO ji = 1, jpi
128                           zqrp = rn_dqdt * ( sst_m(ji,jj) - sf_sst(1)%fnow(ji,jj,1) )
129                           qns(ji,jj) = qns(ji,jj) + zqrp
130                           qrp(ji,jj) = zqrp
131                        END DO
132                     END DO
133                  ENDIF
134               CALL iom_put( "qrp", qrp )                             ! heat flux damping
135            ENDIF
136            !
137            IF( nn_sssr == 1 ) THEN                                   !* Salinity damping term (salt flux only (sfx))
138               zsrp = rn_deds / rday                                  ! from [mm/day] to [kg/m2/s]
139!CDIR COLLAPSE
140               DO jj = 1, jpj
141                  DO ji = 1, jpi
142                     zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths
143                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) ) 
144                     sfx(ji,jj) = sfx(ji,jj) + zerp                 ! salt flux
145                     erp(ji,jj) = zerp / MAX( sss_m(ji,jj), 1.e-20 ) ! converted into an equivalent volume flux (diagnostic only)
146                  END DO
147               END DO
148               CALL iom_put( "erp", erp )                             ! freshwater flux damping
149               !
150            ELSEIF( nn_sssr == 2 ) THEN                               !* Salinity damping term (volume flux (emp) and associated heat flux (qns)
151               zsrp = rn_deds / rday                                  ! from [mm/day] to [kg/m2/s]
152               zerp_bnd = rn_sssr_bnd / rday                          !       -              -   
153!CDIR COLLAPSE
154               DO jj = 1, jpj
155                  DO ji = 1, jpi                           
156                     zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths
157                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) )   &
158                        &        / MAX(  sss_m(ji,jj), 1.e-20   )
159                     IF( ln_sssr_bnd )   zerp = SIGN( 1., zerp ) * MIN( zerp_bnd, ABS(zerp) )
160                     emp(ji,jj) = emp (ji,jj) + zerp
161                     qns(ji,jj) = qns(ji,jj) - zerp * rcp * sst_m(ji,jj)
162                     erp(ji,jj) = zerp
163                  END DO
164               END DO
165               CALL iom_put( "erp", erp )                             ! freshwater flux damping
166            ENDIF
167            !
168         ENDIF
169         !
170      ENDIF
171      !
172      IF( nn_timing == 1 )  CALL timing_stop('sbc_ssr')
173      !
174   END SUBROUTINE sbc_ssr
175
176 
177   SUBROUTINE sbc_ssr_init
178      !!---------------------------------------------------------------------
179      !!                  ***  ROUTINE sbc_ssr_init  ***
180      !!
181      !! ** Purpose :   initialisation of surface damping term
182      !!
183      !! ** Method  : - Read namelist namsbc_ssr
184      !!              - Read observed SST and/or SSS if required
185      !!---------------------------------------------------------------------
186      INTEGER  ::   ji, jj   ! dummy loop indices
187      REAL(wp) ::   zerp     ! local scalar for evaporation damping
188      REAL(wp) ::   zqrp     ! local scalar for heat flux damping
189      REAL(wp) ::   zsrp     ! local scalar for unit conversion of rn_deds factor
190      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor
191      INTEGER  ::   ierror   ! return error code
192      !!
193      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
194      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read
195      NAMELIST/namsbc_ssr/ cn_dir, nn_sstr, nn_sssr, rn_dqdt, rn_deds, sn_sst, sn_sss, ln_sssr_bnd, rn_sssr_bnd, ln_UKMO_haney
196      INTEGER     ::  ios
197      !!----------------------------------------------------------------------
198      !
199 
200      REWIND( numnam_ref )              ! Namelist namsbc_ssr in reference namelist :
201      READ  ( numnam_ref, namsbc_ssr, IOSTAT = ios, ERR = 901)
202901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_ssr in reference namelist', lwp )
203
204      REWIND( numnam_cfg )              ! Namelist namsbc_ssr in configuration namelist :
205      READ  ( numnam_cfg, namsbc_ssr, IOSTAT = ios, ERR = 902 )
206902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_ssr in configuration namelist', lwp )
207      IF(lwm) WRITE ( numond, namsbc_ssr )
208
209      IF(lwp) THEN                 !* control print
210         WRITE(numout,*)
211         WRITE(numout,*) 'sbc_ssr : SST and/or SSS damping term '
212         WRITE(numout,*) '~~~~~~~ '
213         WRITE(numout,*) '   Namelist namsbc_ssr :'
214         WRITE(numout,*) '      SST restoring term (Yes=1)             nn_sstr     = ', nn_sstr
215         WRITE(numout,*) '      SSS damping term (Yes=1, salt flux)    nn_sssr     = ', nn_sssr
216         WRITE(numout,*) '                       (Yes=2, volume flux) '
217         WRITE(numout,*) '      dQ/dT (restoring magnitude on SST)     rn_dqdt     = ', rn_dqdt, ' W/m2/K'
218         WRITE(numout,*) '      dE/dS (restoring magnitude on SST)     rn_deds     = ', rn_deds, ' mm/day'
219         WRITE(numout,*) '      flag to bound erp term                 ln_sssr_bnd = ', ln_sssr_bnd
220         WRITE(numout,*) '      ABS(Max./Min.) erp threshold           rn_sssr_bnd = ', rn_sssr_bnd, ' mm/day'
221         WRITE(numout,*) '      Haney forcing                          ln_UKMO_haney = ', ln_UKMO_haney
222      ENDIF
223      !
224      !                            !* Allocate erp and qrp array
225      ALLOCATE( qrp(jpi,jpj), erp(jpi,jpj), STAT=ierror )
226      IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate erp and qrp array' )
227      !
228      IF( nn_sstr == 1 ) THEN      !* set sf_sst structure & allocate arrays
229         !
230         ALLOCATE( sf_sst(1), STAT=ierror )
231         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst structure' )
232         ALLOCATE( sf_sst(1)%fnow(jpi,jpj,1), STAT=ierror )
233         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst now array' )
234         !
235         ! fill sf_sst with sn_sst and control print
236         CALL fld_fill( sf_sst, (/ sn_sst /), cn_dir, 'sbc_ssr', 'SST restoring term toward SST data', 'namsbc_ssr' )
237         IF( sf_sst(1)%ln_tint )   ALLOCATE( sf_sst(1)%fdta(jpi,jpj,1,2), STAT=ierror )
238         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst data array' )
239         !
240      ENDIF
241      !
242      IF( nn_sssr >= 1 ) THEN      !* set sf_sss structure & allocate arrays
243         !
244         ALLOCATE( sf_sss(1), STAT=ierror )
245         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss structure' )
246         ALLOCATE( sf_sss(1)%fnow(jpi,jpj,1), STAT=ierror )
247         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss now array' )
248         !
249         ! fill sf_sss with sn_sss and control print
250         CALL fld_fill( sf_sss, (/ sn_sss /), cn_dir, 'sbc_ssr', 'SSS restoring term toward SSS data', 'namsbc_ssr' )
251         IF( sf_sss(1)%ln_tint )   ALLOCATE( sf_sss(1)%fdta(jpi,jpj,1,2), STAT=ierror )
252         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss data array' )
253         !
254      ENDIF
255      !
256      !                            !* Initialize qrp and erp if no restoring
257      IF( nn_sstr /= 1                   )   qrp(:,:) = 0._wp
258      IF( nn_sssr /= 1 .OR. nn_sssr /= 2 )   erp(:,:) = 0._wp
259      !
260   END SUBROUTINE sbc_ssr_init
261     
262   !!======================================================================
263END MODULE sbcssr
Note: See TracBrowser for help on using the repository browser.