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.
sbcblk_skin.F90 in NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC – NEMO

source: NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk_skin.F90 @ 11111

Last change on this file since 11111 was 11111, checked in by laurent, 5 years ago

LB: skin temperature now used by ECMWF bulk algorithm, and can be xios-ed ; new module "SBC/sbcblk_skin.F90" ; all physical constants defined in SBC now moved to "DOM/phycst.F90"; all air-properties and other ABL functions gathered in a new module: "SBC/sbcblk_phymbl.F90".
File size: 9.7 KB
Line 
1MODULE sbcblk_skin
2   !!======================================================================
3   !!                   ***  MODULE  sbcblk_skin  ***
4   !! Computes:
5   !!   * the surface skin temperature (aka SSST) based on the cool-skin/warm-layer
6   !!     scheme used at ECMWF
7   !!    Using formulation/param. of IFS of ECMWF (cycle 40r1)
8   !!   based on IFS doc (avaible online on the ECMWF's website)
9   !!
10   !!   From routine turb_ecmwf maintained and developed in AeroBulk
11   !!            (https://github.com/brodeau/aerobulk)
12   !!
13   !! ** Author: L. Brodeau, june 2019 / AeroBulk (https://github.com/brodeau/aerobulk)
14   !!----------------------------------------------------------------------
15   !! History :  4.0  !  2016-02  (L.Brodeau)   Original code
16   !!----------------------------------------------------------------------
17
18   !!----------------------------------------------------------------------
19   !!  cswl_ecmwf : computes the surface skin temperature (aka SSST)
20   !!               based on the cool-skin/warm-layer scheme used at ECMWF
21   !!----------------------------------------------------------------------
22   USE oce             ! ocean dynamics and tracers
23   USE dom_oce         ! ocean space and time domain
24   USE phycst          ! physical constants
25   USE sbc_oce         ! Surface boundary condition: ocean fields
26
27   USE sbcblk_phymbl  !LOLO: all thermodynamics functions, rho_air, q_sat, etc...
28
29   USE lib_mpp         ! distribued memory computing library
30   USE in_out_manager  ! I/O manager
31   USE lib_fortran     ! to use key_nosignedzero
32
33
34   IMPLICIT NONE
35   PRIVATE
36   
37   PUBLIC :: CSWL_ECMWF   ! called by sbcblk_algo_*.F90
38
39   !! Cool-Skin / Warm-Layer related parameters:
40   REAL(wp), PARAMETER :: rdt0 = 3600.*1.5 !: time step
41   REAL(wp), PARAMETER :: rd0  = 3.        !: Depth scale [m], "d" in Eq.11 (Zeng & Beljaars 2005)
42   REAL(wp), PARAMETER :: rNu0 = 0.5       !: Nu (exponent of temperature profile) Eq.11
43   !                                       !: (Zeng & Beljaars 2005) !: set to 0.5 instead of
44   !                                       !: 0.3 to respect a warming of +3 K in calm
45   !                                       !: condition for the insolation peak of +1000W/m^2
46   
47   INTEGER, PARAMETER :: nb_itt_wl = 10    !: number of sub-itterations for solving the differential equation in warm-layer part
48   !                                       !:  => use "nb_itt_wl = 1" for No itteration! => way cheaper !!!
49   !                                       !:    => assumes balance between the last 2 terms of Eq.11 (lhs of eq.11 = 0)
50   !                                       !:    => in that case no need for sub-itterations !
51   !                                       !:    => ACCEPTABLE IN MOST CONDITIONS ! (UNLESS: sunny + very calm/low-wind conditions)
52   !                                       !:  => Otherwize use "nb_itt_wl = 10"
53   !
54   !!----------------------------------------------------------------------
55CONTAINS
56
57
58   SUBROUTINE CSWL_ECMWF( pQsw, pQnsol, pustar, pSST, pTs )
59      !!---------------------------------------------------------------------
60      !!
61      !!  Cool-Skin Warm-Layer scheme according to Zeng & Beljaars, 2005 (GRL)
62      !!  " A prognostic scheme of sea surface skin temperature for modeling and data assimilation "
63      !!
64      !!    As included in IFS Cy40   /  E.C.M.W.F.
65      !!     ------------------------------------------------------------------
66      !!
67      !!  **   INPUT:
68      !!
69      !!     *pQsw*       net solar radiative flux to the ocean
70      !!     *pQnsol*     net non-solar heat flux to the ocean
71      !!     *pustar*     friction velocity u*
72      !!     *pSST*       SST
73      !!
74      !!   **  INPUT/OUTPUT:
75      !!     *pTs*  : as input  =>  previous estimate of skin temperature
76      !!             as output =>  new estimate of skin temperature
77      !!
78      !!------------------------------------------------------------------
79      REAL(wp), DIMENSION(jpi,jpj), INTENT(IN)    :: pQsw
80      REAL(wp), DIMENSION(jpi,jpj), INTENT(IN)    :: pQnsol
81      REAL(wp), DIMENSION(jpi,jpj), INTENT(IN)    :: pustar
82      REAL(wp), DIMENSION(jpi,jpj), INTENT(IN)    :: pSST
83      !
84      REAL(wp), DIMENSION(jpi,jpj), INTENT(INOUT) :: pTs
85      !
86      INTEGER :: ji,jj
87      !
88      REAL(wp) :: rmult, &
89         & zRhoCp_w, &
90         & ZCON2,ZCON3,ZCON4,ZCON5, zQnsol ,zQnet, zlamb, zdelta,&
91         & ZSRD,ZDSST,ZZ,ZEPDU2,&
92         & ZFI,zdL,zdL2, ztmp, &
93         & ZPHI,ZROADRW, &
94         & zus_a, &
95         & zdt, zfs, zsgn
96      !
97      REAL(wp), DIMENSION(jpi,jpj) :: &
98         &  zalpha_w, &       !: thermal expansion coefficient of seawater
99         & zus_w, zus_w2, &   !: u* and u*^2 in water
100         &       zdT_c,   &   !: cool skin temperature increment !lolo rm!
101         &       zdT_w        !: warm skin temperature increment
102      !
103      INTEGER :: nbi, jwl
104      !!------------------------------------------------------------------
105      !
106      !     1. Initialize constants for ocean warm layer and cool skin
107      !
108      !     1.1 General
109      !
110      ZEPDU2  = 0.01_wp   !    security constant for velocity**2   (m2/s2)
111      ZROADRW = rho0_a/rho0_w          ! Density ratio                      (-)
112      zRhoCp_w = rho0_w*rCp0_w
113      !
114      !     1.2C Warm layer parametrization constants
115      !
116      !    ZFI = Fraction of solar radiation absorbed in warm layer (-)
117      ZFI = 1._wp -0.28_wp*EXP(-71.5_wp*rd0) -0.27_wp*EXP(-2.8_wp*rd0) - 0.45_wp*EXP(-0.07_wp*rd0)  !: Eq. 8.135
118      !
119      ZCON3 = rd0*vkarmn*grav/(ZROADRW)**1.5_wp
120      ZCON4 = (rNu0 + 1._wp)*vkarmn/rd0
121      ZCON5 = (rNu0 + 1._wp)/(rNu0*rd0)
122      !
123      !     1.3 Cool skin parametrization constants
124      ZCON2 = 16._wp*grav*zRhoCp_w*rnu0_w**3/(rk0_w*rk0_w)
125      !
126      ! Friction velocities
127      ! "MAX( pustar(:,:), 1.E-4)" is u* in the air !
128      zus_w(:,:)  = MAX( pustar(:,:), 1.E-4_wp)*SQRT(ZROADRW)       ! u* in the water
129      zus_w2(:,:) = zus_w(:,:)*zus_w(:,:)
130      !
131      ! Ocean buoyancy
132      zalpha_w(:,:) = MAX( 1.E-5_wp , 1.E-5_wp*(pTs(:,:) - rt0) ) ! thermal expansion coefficient of water
133      !
134      zdT_c = 0._wp
135      zdT_w = 0._wp
136      !
137      !  3. Cool skin (Fairall et al. 1996)
138      !------------------------------------
139      DO jj = 1, jpj
140         DO ji = 1, jpi
141            !
142            ! Non-solar heat loss to the atmosphere:
143            zQnsol = MAX( 1._wp , - pQnsol(ji,jj) )
144
145            zlamb = 6._wp*(1._wp + (zQnsol*zalpha_w(ji,jj)*ZCON2/(zus_w2(ji,jj)*zus_w2(ji,jj)))**0.75)**(-1._wp/3._wp)
146
147            zdelta = zlamb*rnu0_w/zus_w(ji,jj)
148
149            !    Solar absorption
150            zfs   = 0.065_wp + 11._wp*zdelta - (6.6E-5_wp/zdelta)*(1._wp - EXP(-zdelta/8.E-4_wp))  ! Eq. 8.131 / IFS cy40r1, doc, Part IV,
151            zfs   = MAX(zfs , 0.01_wp)
152            zQnet = MAX( 1._wp , -zfs*pQsw(ji,jj) + zQnsol )
153            zdT_c(ji,jj) = -zdelta*zQnet/rk0_w
154
155         END DO
156      END DO
157
158      ! 2.2 Warm layer; formulation C (Xubin Zeng)
159      !!--------------------------------------------
160      nbi = MAX( 1 , nb_itt_wl)
161
162      IF( nbi > 1 ) THEN
163         !! Itterating for warm-layer solution
164         zdt   = rdt0/REAL(nbi)
165         rmult = 1._wp
166      ELSE
167         !! No itteration! Way cheaper !!!
168         !! Assuming balance between the last 2 terms of Eq. 11 (lhs of 11 = 0):
169         !!   => in that case no need for sub-itterations !
170         !! Acceptable in most conditions (UNLESS: sunny + no-wind conditions)
171         !!
172         zdt   = 1._wp
173         rmult = 0._wp
174      END IF
175
176      DO jwl = 1, nbi  ! itteration to solve implicitely equation for warm layer
177
178         DO jj = 1, jpj
179            DO ji = 1, jpi
180
181               ZDSST = pTs(ji,jj) - pSST(ji,jj) - zdT_c(ji,jj)
182
183               !! Buoyancy flux and stability parameter (zdl = -z/L) in water
184               !
185               !! Qt/(rho_w*Cpw):
186               ZSRD = ( pQsw(ji,jj)*ZFI + pQnsol(ji,jj) )/zRhoCp_w
187               !
188               zsgn = 0.5_wp + SIGN(0.5_wp, ZSRD)  ! ZSRD > 0. => 1.  / ZSRD < 0. => 0.
189               ztmp = MAX(ZDSST,0._wp)
190               zdl = (zsgn+1._wp) * ( zus_w2(ji,jj) * SQRT(ztmp/(5._wp*rd0*grav*zalpha_w(ji,jj)/rNu0)) ) & ! (ZDSST > 0.0 .AND. ZSRD < 0.0)
191                  &  +    zsgn    *  ZSRD                                                                  !   otherwize
192               !
193               zus_a = MAX( pustar(ji,jj), 1.E-4_wp )
194               zdL = ZCON3*zalpha_w(ji,jj)*zdL/(zus_a*zus_a*zus_a)
195
196               !! Stability function Phi_t(-z/L) (zdL is -z/L) :
197               zsgn = 0.5_wp + SIGN(0.5_wp, zdL)  ! zdl > 0. => 1.  / zdl < 0. => 0.
198               zdL2 = zdL*zdL
199               ZPHI =    zsgn      * ( 1._wp + (5._wp*zdL + 4._wp*zdL2)/(1._wp + 3._wp*zdL + 0.25_wp*zdL2) ) &  ! (zdL > 0) Takaya et al.
200                  & + (zsgn+1._wp) * ( 1._wp/SQRT(1._wp - 16._wp*(-ABS(zdL))) )        ! (zdl < 0) Eq. 8.136
201               !
202               !! FOR zdL > 0.0, old relations:
203               !         ZPHI = 1.+5._wp*zdL                                ! Eq. 8.136 (Large et al. 1994)
204               !         ZPHI = 1.+5.0*(zdL+zdL**2)/(1.0+3.0*zdL+zdL**2) ! SHEBA, Grachev et al. 2007
205
206               !! Solving 11 by itteration with time step of zdt...
207               ZZ = rmult*1._wp + ZCON4*zdt*zus_w(ji,jj)/ZPHI
208               ZZ = SIGN( MAX(ABS(ZZ) , 1e-4_wp), ZZ )
209               zdT_w(ji,jj) = MAX( 0._wp , (rmult*ZDSST + ZCON5*ZSRD*zdt)/ZZ )
210
211            END DO ! DO ji = 1, jpi
212         END DO ! DO jj = 1, jpj
213
214         ! 3. Apply warm layer and cool skin effects
215         !------------------------------------------
216         pTs(:,:) = pSST(:,:) + zdT_w(:,:) + zdT_c(:,:)
217
218      END DO  ! DO jwl = 1, nbi !: sub-itteration
219
220   END SUBROUTINE CSWL_ECMWF
221
222   !!======================================================================
223END MODULE sbcblk_skin
Note: See TracBrowser for help on using the repository browser.