source: branches/UKMO/dev_r5518_GO6_package_FOAMv14/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssr.F90 @ 11442

Last change on this file since 11442 was 11442, checked in by mattmartin, 21 months ago

Introduction of stochastic physics in NEMO, based on Andrea Storto's code.
For details, see ticket https://code.metoffice.gov.uk/trac/utils/ticket/251.

File size: 12.9 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   USE stopack
28   USE wrk_nemo       ! Memory Allocation
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   sbc_ssr        ! routine called in sbcmod
34   PUBLIC   sbc_ssr_init   ! routine called in sbcmod
35
36   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   erp   !: evaporation damping   [kg/m2/s]
37   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   qrp   !: heat flux damping        [w/m2]
38
39   !                                   !!* Namelist namsbc_ssr *
40   INTEGER, PUBLIC ::   nn_sstr         ! SST/SSS restoring indicator
41   INTEGER, PUBLIC ::   nn_sssr         ! SST/SSS restoring indicator
42   REAL(wp)        ::   rn_dqdt         ! restoring factor on SST and SSS
43   REAL(wp)        ::   rn_deds         ! restoring factor on SST and SSS
44   LOGICAL         ::   ln_sssr_bnd     ! flag to bound erp term
45   REAL(wp)        ::   rn_sssr_bnd     ! ABS(Max./Min.) value of erp term [mm/day]
46
47   REAL(wp) , ALLOCATABLE, DIMENSION(:) ::   buffer   ! Temporary buffer for exchange
48   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sst   ! structure of input SST (file informations, fields read)
49   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sss   ! structure of input SSS (file informations, fields read)
50
51   !! * Substitutions
52#  include "domzgr_substitute.h90"
53   !!----------------------------------------------------------------------
54   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
55   !! $Id$
56   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
57   !!----------------------------------------------------------------------
58CONTAINS
59
60   SUBROUTINE sbc_ssr( kt )
61      !!---------------------------------------------------------------------
62      !!                     ***  ROUTINE sbc_ssr  ***
63      !!
64      !! ** Purpose :   Add to heat and/or freshwater fluxes a damping term
65      !!                toward observed SST and/or SSS.
66      !!
67      !! ** Method  : - Read namelist namsbc_ssr
68      !!              - Read observed SST and/or SSS
69      !!              - at each nscb time step
70      !!                   add a retroaction term on qns    (nn_sstr = 1)
71      !!                   add a damping term on sfx        (nn_sssr = 1)
72      !!                   add a damping term on emp        (nn_sssr = 2)
73      !!---------------------------------------------------------------------
74      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
75      !!
76      INTEGER  ::   ji, jj   ! dummy loop indices
77      REAL(wp) ::   zerp     ! local scalar for evaporation damping
78      REAL(wp) ::   zqrp     ! local scalar for heat flux damping
79      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor
80      REAL(wp), POINTER, DIMENSION(:,:) :: rn_dqdt_s, zsrp
81      INTEGER  ::   ierror   ! return error code
82      !!
83      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
84      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read
85      !!----------------------------------------------------------------------
86      !
87      IF( nn_timing == 1 )  CALL timing_start('sbc_ssr')
88      !
89      IF( nn_sstr + nn_sssr /= 0 ) THEN
90         !
91         IF( nn_sstr == 1)   CALL fld_read( kt, nn_fsbc, sf_sst )   ! Read SST data and provides it at kt
92         IF( nn_sssr >= 1)   CALL fld_read( kt, nn_fsbc, sf_sss )   ! Read SSS data and provides it at kt
93         !
94         !                                         ! ========================= !
95         IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN      !    Add restoring term     !
96            !                                      ! ========================= !
97            !
98            IF( nn_sstr == 1 ) THEN                                   !* Temperature restoring term
99
100               CALL wrk_alloc( jpi, jpj, rn_dqdt_s )
101               rn_dqdt_s(:,:) = rn_dqdt
102
103               IF( ln_stopack .AND. nn_spp_dqdt > 0 ) &
104                  & CALL spp_gen( kt, rn_dqdt_s, nn_spp_dqdt, rn_dqdt_sd, jk_spp_dqdt )
105               DO jj = 1, jpj
106                  DO ji = 1, jpi
107                     zqrp = rn_dqdt_s(ji,jj) * ( sst_m(ji,jj) - sf_sst(1)%fnow(ji,jj,1) )
108                     qns(ji,jj) = qns(ji,jj) + zqrp
109                     qrp(ji,jj) = zqrp
110                  END DO
111               END DO
112               CALL iom_put( "qrp", qrp )                             ! heat flux damping
113               CALL wrk_dealloc( jpi, jpj, rn_dqdt_s )
114            ENDIF
115            !
116            IF( nn_sssr == 1 ) THEN                                   !* Salinity damping term (salt flux only (sfx))
117               CALL wrk_alloc( jpi, jpj, zsrp)
118               zsrp(:,:) = rn_deds
119               IF( ln_stopack .AND. nn_spp_dedt > 0 ) &
120                  & CALL spp_gen(kt, zsrp, nn_spp_dedt, rn_dedt_sd, jk_spp_deds )
121!CDIR COLLAPSE
122               DO jj = 1, jpj
123                  DO ji = 1, jpi
124                     zerp = (zsrp(ji,jj)/rday) * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths
125                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) ) 
126                     sfx(ji,jj) = sfx(ji,jj) + zerp                 ! salt flux
127                     erp(ji,jj) = zerp / MAX( sss_m(ji,jj), 1.e-20 ) ! converted into an equivalent volume flux (diagnostic only)
128                  END DO
129               END DO
130               CALL iom_put( "erp", erp )                             ! freshwater flux damping
131               CALL wrk_dealloc( jpi,jpj, zsrp )
132               !
133            ELSEIF( nn_sssr == 2 ) THEN                               !* Salinity damping term (volume flux (emp) and associated heat flux (qns)
134               CALL wrk_alloc( jpi, jpj, zsrp)
135               zsrp(:,:) = rn_deds
136               IF( ln_stopack .AND. nn_spp_dedt > 0 ) &
137                  & CALL spp_gen( kt, zsrp, nn_spp_dedt, rn_dedt_sd, jk_spp_deds )
138               zerp_bnd = rn_sssr_bnd / rday                          !       -              -   
139!CDIR COLLAPSE
140               DO jj = 1, jpj
141                  DO ji = 1, jpi                           
142                     zerp = (zsrp(ji,jj)/rday) * ( 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                        &        / MAX(  sss_m(ji,jj), 1.e-20   )
145                     IF( ln_sssr_bnd )   zerp = SIGN( 1., zerp ) * MIN( zerp_bnd, ABS(zerp) )
146                     emp(ji,jj) = emp (ji,jj) + zerp
147                     qns(ji,jj) = qns(ji,jj) - zerp * rcp * sst_m(ji,jj)
148                     erp(ji,jj) = zerp
149                  END DO
150               END DO
151               CALL iom_put( "erp", erp )                             ! freshwater flux damping
152               CALL wrk_dealloc( jpi,jpj,zsrp )
153            ENDIF
154            !
155         ENDIF
156         !
157      ENDIF
158      !
159      IF( nn_timing == 1 )  CALL timing_stop('sbc_ssr')
160      !
161   END SUBROUTINE sbc_ssr
162
163 
164   SUBROUTINE sbc_ssr_init
165      !!---------------------------------------------------------------------
166      !!                  ***  ROUTINE sbc_ssr_init  ***
167      !!
168      !! ** Purpose :   initialisation of surface damping term
169      !!
170      !! ** Method  : - Read namelist namsbc_ssr
171      !!              - Read observed SST and/or SSS if required
172      !!---------------------------------------------------------------------
173      INTEGER  ::   ji, jj   ! dummy loop indices
174      REAL(wp) ::   zerp     ! local scalar for evaporation damping
175      REAL(wp) ::   zqrp     ! local scalar for heat flux damping
176      REAL(wp) ::   zsrp     ! local scalar for unit conversion of rn_deds factor
177      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor
178      INTEGER  ::   ierror   ! return error code
179      !!
180      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
181      TYPE(FLD_N) ::   sn_sst, sn_sss        ! informations about the fields to be read
182      NAMELIST/namsbc_ssr/ cn_dir, nn_sstr, nn_sssr, rn_dqdt, rn_deds, sn_sst, sn_sss, ln_sssr_bnd, rn_sssr_bnd
183      INTEGER     ::  ios
184      !!----------------------------------------------------------------------
185      !
186 
187      REWIND( numnam_ref )              ! Namelist namsbc_ssr in reference namelist :
188      READ  ( numnam_ref, namsbc_ssr, IOSTAT = ios, ERR = 901)
189901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_ssr in reference namelist', lwp )
190
191      REWIND( numnam_cfg )              ! Namelist namsbc_ssr in configuration namelist :
192      READ  ( numnam_cfg, namsbc_ssr, IOSTAT = ios, ERR = 902 )
193902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_ssr in configuration namelist', lwp )
194      IF(lwm) WRITE ( numond, namsbc_ssr )
195
196      IF(lwp) THEN                 !* control print
197         WRITE(numout,*)
198         WRITE(numout,*) 'sbc_ssr : SST and/or SSS damping term '
199         WRITE(numout,*) '~~~~~~~ '
200         WRITE(numout,*) '   Namelist namsbc_ssr :'
201         WRITE(numout,*) '      SST restoring term (Yes=1)             nn_sstr     = ', nn_sstr
202         WRITE(numout,*) '      SSS damping term (Yes=1, salt flux)    nn_sssr     = ', nn_sssr
203         WRITE(numout,*) '                       (Yes=2, volume flux) '
204         WRITE(numout,*) '      dQ/dT (restoring magnitude on SST)     rn_dqdt     = ', rn_dqdt, ' W/m2/K'
205         WRITE(numout,*) '      dE/dS (restoring magnitude on SST)     rn_deds     = ', rn_deds, ' mm/day'
206         WRITE(numout,*) '      flag to bound erp term                 ln_sssr_bnd = ', ln_sssr_bnd
207         WRITE(numout,*) '      ABS(Max./Min.) erp threshold           rn_sssr_bnd = ', rn_sssr_bnd, ' mm/day'
208      ENDIF
209      !
210      !                            !* Allocate erp and qrp array
211      ALLOCATE( qrp(jpi,jpj), erp(jpi,jpj), STAT=ierror )
212      IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate erp and qrp array' )
213      !
214      IF( nn_sstr == 1 ) THEN      !* set sf_sst structure & allocate arrays
215         !
216         ALLOCATE( sf_sst(1), STAT=ierror )
217         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst structure' )
218         ALLOCATE( sf_sst(1)%fnow(jpi,jpj,1), STAT=ierror )
219         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst now array' )
220         !
221         ! fill sf_sst with sn_sst and control print
222         CALL fld_fill( sf_sst, (/ sn_sst /), cn_dir, 'sbc_ssr', 'SST restoring term toward SST data', 'namsbc_ssr' )
223         IF( sf_sst(1)%ln_tint )   ALLOCATE( sf_sst(1)%fdta(jpi,jpj,1,2), STAT=ierror )
224         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sst data array' )
225         !
226      ENDIF
227      !
228      IF( nn_sssr >= 1 ) THEN      !* set sf_sss structure & allocate arrays
229         !
230         ALLOCATE( sf_sss(1), STAT=ierror )
231         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss structure' )
232         ALLOCATE( sf_sss(1)%fnow(jpi,jpj,1), STAT=ierror )
233         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss now array' )
234         !
235         ! fill sf_sss with sn_sss and control print
236         CALL fld_fill( sf_sss, (/ sn_sss /), cn_dir, 'sbc_ssr', 'SSS restoring term toward SSS data', 'namsbc_ssr' )
237         IF( sf_sss(1)%ln_tint )   ALLOCATE( sf_sss(1)%fdta(jpi,jpj,1,2), STAT=ierror )
238         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_ssr: unable to allocate sf_sss data array' )
239         !
240      ENDIF
241      !
242      !                            !* Initialize qrp and erp if no restoring
243      IF( nn_sstr /= 1                   )   qrp(:,:) = 0._wp
244      IF( nn_sssr /= 1 .OR. nn_sssr /= 2 )   erp(:,:) = 0._wp
245      !
246   END SUBROUTINE sbc_ssr_init
247     
248   !!======================================================================
249END MODULE sbcssr
Note: See TracBrowser for help on using the repository browser.