source: LMDZ6/branches/Ocean_skin/libf/phylmd/surf_ocean_mod.F90 @ 3687

Last change on this file since 3687 was 3687, checked in by lguez, 4 years ago

Modify sensible heat due to rain sent to the ocean

Modify the sensible heat flux due to rain which is sent to the
ocean. Replace the computation of sens_prec_liq in procedure
calcul_fluxs by a call to sens_heat_rain. Set sens_prec_sol in
procedure calcul_fluxs to 0 because, for now, sens_heat_rain is
supposed to account for both rain and snow.

For the call to sens_heat_rain in procedure calcul_fluxs, we need
an additional dummy argument rhoa of calcul_fluxs. Add dummy
argument rhoa to ocean_cpl_noice, ocean_forced_noice,
ocean_forced_ice and ocean_cpl_ice because we need to pass it down
to calcul_fluxs.

Change the dimension of sens_prec_liq and sens_prec_sol in
procedures calcul_fluxs, ocean_cpl_noice, ocean_cpl_ice,
ocean_forced_noice, ocean_forced_ice, cpl_send_ocean_fields and
cpl_send_seaice_fields to knon.

In procedures ocean_forced_noice and ocean_cpl_noice, promote
sens_prec_liq from local variable to dummy argument because we need
it in surf_ocean. Remove useless initialization of sens_prec_liq
and sens_prec_sol in ocean_cpl_noice, ocean_cpl_ice,
ocean_forced_ice and ocean_forced_noice: they are intent out in
calcul_fluxs.

Remove variable rf of module phys_output_var_mod, we use
sens_prec_liq instead. Remove local variable yrf of procedure
pbl_surface. rf and yrf appeared in pbl_surface only to be output.
Remove variable o_rf of module phys_output_ctrlout_mod. Remove
dummy argument rf of procedure surf_ocean.

Do not call sens_heat_rain in surf_ocean since we now call it
from calcul_fluxs. Move the computation of rhoa in surf_ocean
before the calls to ocean_cpl_noice and ocean_forced_noice.

Add the computation of rhoa in surf_seaice, to pass it down to
ocean_cpl_ice and ocean_forced_ice.

If activate_ocean_skin == 1 then the results are changed because the
call to sens_heat_rain in calcul_fluxs now uses the SST from the
current time step of physics. On this point, the present revision
reverses revision [3463].
  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 13.9 KB
Line 
1!
2! $Id: surf_ocean_mod.F90 3687 2020-05-27 14:52:16Z lguez $
3!
4MODULE surf_ocean_mod
5
6  IMPLICIT NONE
7
8CONTAINS
9  !
10  !******************************************************************************
11  !
12  SUBROUTINE surf_ocean(rlon, rlat, swnet, lwnet, alb1, &
13       windsp, rmu0, fder, tsurf_in, &
14       itime, dtime, jour, knon, knindex, &
15       p1lay, z1lay, cdragh, cdragm, precip_rain, precip_snow, temp_air, spechum, &
16       AcoefH, AcoefQ, BcoefH, BcoefQ, &
17       AcoefU, AcoefV, BcoefU, BcoefV, &
18       ps, u1, v1, gustiness, rugoro, pctsrf, &
19       snow, qsurf, agesno, &
20       z0m, z0h, SFRWL, alb_dir_new, alb_dif_new, evap, fluxsens, fluxlat, &
21       tsurf_new, dflux_s, dflux_l, lmt_bils, &
22       flux_u1, flux_v1, t_int, s_int, ds_ns, dt_ns, dter, dser, tkt, tks, taur)
23
24    use albedo, only: alboc, alboc_cd
25    use bulk_flux_m, only: bulk_flux
26    USE dimphy, ONLY: klon, zmasq
27    USE surface_data, ONLY     : type_ocean
28    USE ocean_forced_mod, ONLY : ocean_forced_noice
29    USE ocean_slab_mod, ONLY   : ocean_slab_noice
30    USE ocean_cpl_mod, ONLY    : ocean_cpl_noice
31    USE indice_sol_mod, ONLY : nbsrf, is_oce
32    USE limit_read_mod
33    use config_ocean_skin_m, only: activate_ocean_skin
34    !
35    ! This subroutine will make a call to ocean_XXX_noice according to the ocean mode (force,
36    ! slab or couple). The calculations of albedo and rugosity for the ocean surface are
37    ! done in here because they are identical for the different modes of ocean.
38
39
40    INCLUDE "YOMCST.h"
41
42    include "clesphys.h"
43    ! for cycle_diurne and for iflag_z0_oce==-1 (prescribed z0)
44
45    ! Input variables
46    !******************************************************************************
47    INTEGER, INTENT(IN)                      :: itime, jour, knon
48    INTEGER, DIMENSION(klon), INTENT(IN)     :: knindex
49    REAL, INTENT(IN)                         :: dtime
50    REAL, DIMENSION(klon), INTENT(IN)        :: rlon, rlat
51    REAL, DIMENSION(klon), INTENT(IN)        :: swnet  ! net shortwave radiation at surface 
52    REAL, DIMENSION(klon), INTENT(IN)        :: lwnet  ! net longwave radiation at surface 
53    REAL, DIMENSION(klon), INTENT(IN)        :: alb1   ! albedo in visible SW interval
54    REAL, DIMENSION(klon), INTENT(IN)        :: windsp ! wind at 10 m, in m s-1
55    REAL, DIMENSION(klon), INTENT(IN)        :: rmu0 
56    REAL, DIMENSION(klon), INTENT(IN)        :: fder
57    REAL, INTENT(IN):: tsurf_in(klon) ! defined only for subscripts 1:knon
58    REAL, DIMENSION(klon), INTENT(IN)        :: p1lay,z1lay ! pression (Pa) et altitude (m) du premier niveau
59    REAL, DIMENSION(klon), INTENT(IN)        :: cdragh
60    REAL, DIMENSION(klon), INTENT(IN)        :: cdragm
61    REAL, DIMENSION(klon), INTENT(IN)        :: precip_rain, precip_snow
62    REAL, DIMENSION(klon), INTENT(IN)        :: temp_air, spechum
63    REAL, DIMENSION(klon), INTENT(IN)        :: AcoefH, AcoefQ, BcoefH, BcoefQ
64    REAL, DIMENSION(klon), INTENT(IN)        :: AcoefU, AcoefV, BcoefU, BcoefV
65    REAL, DIMENSION(klon), INTENT(IN)        :: ps
66    REAL, DIMENSION(klon), INTENT(IN)        :: u1, v1, gustiness
67    REAL, DIMENSION(klon), INTENT(IN)        :: rugoro
68    REAL, DIMENSION(klon,nbsrf), INTENT(IN)  :: pctsrf
69
70    ! In/Output variables
71    !******************************************************************************
72    REAL, DIMENSION(klon), INTENT(INOUT)     :: snow
73    REAL, DIMENSION(klon), INTENT(INOUT)     :: qsurf
74    REAL, DIMENSION(klon), INTENT(INOUT)     :: agesno
75    REAL, DIMENSION(klon), INTENT(inOUT):: z0h
76
77    real, intent(inout):: s_int(:) ! (knon) ocean-air interface salinity, in ppt
78
79    REAL, intent(inout):: ds_ns(:) ! (knon)
80    ! "delta salinity near surface". Salinity variation in the
81    ! near-surface turbulent layer. That is subskin salinity minus
82    ! foundation salinity. In ppt.
83
84    REAL, intent(inout):: dt_ns(:) ! (knon)
85    ! "delta temperature near surface". Temperature variation in the
86    ! near-surface turbulent layer. That is subskin temperature
87    ! minus foundation temperature. (Can be negative.) In K.
88
89    ! Output variables
90    !******************************************************************************
91    REAL, DIMENSION(klon), INTENT(OUT)       :: z0m
92    !albedo SB >>>
93    !    REAL, DIMENSION(klon), INTENT(OUT)       :: alb1_new  ! new albedo in visible SW interval
94    !    REAL, DIMENSION(klon), INTENT(OUT)       :: alb2_new  ! new albedo in near IR interval
95    REAL, DIMENSION(6), INTENT(IN)          :: SFRWL
96    REAL, DIMENSION(klon,nsw), INTENT(OUT)       :: alb_dir_new,alb_dif_new
97    !albedo SB <<<     
98    REAL, DIMENSION(klon), INTENT(OUT)       :: evap, fluxsens, fluxlat
99    REAL, INTENT(OUT):: tsurf_new(klon) ! sea surface temperature, in K
100    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
101    REAL, DIMENSION(klon), INTENT(OUT)       :: lmt_bils
102    REAL, DIMENSION(klon), INTENT(OUT)       :: flux_u1, flux_v1
103
104    REAL, intent(out):: t_int(:) ! (knon) ocean-air interface temperature, in K
105
106    REAL, intent(out):: dter(:) ! (knon)
107    ! Temperature variation in the diffusive microlayer, that is
108    ! ocean-air interface temperature minus subskin temperature. In
109    ! K.
110
111    REAL, intent(out):: dser(:) ! (knon)
112    ! Salinity variation in the diffusive microlayer, that is
113    ! ocean-air interface salinity minus subskin salinity. In ppt.
114
115    REAL, intent(out):: tkt(:) ! (knon)
116    ! épaisseur (m) de la couche de diffusion thermique (microlayer)
117    ! cool skin thickness
118
119    REAL, intent(out):: tks(:) ! (knon)
120    ! épaisseur (m) de la couche de diffusion de masse (microlayer)
121
122    REAL, intent(out):: taur(:) ! (knon)
123    ! momentum flux due to rain, in Pa
124
125    ! Local variables
126    !******************************************************************************
127    INTEGER               :: i, k
128    REAL                  :: tmp
129    REAL, PARAMETER       :: cepdu2=(0.1)**2
130    REAL, DIMENSION(klon) :: alb_eau, z0_lim
131    REAL, DIMENSION(klon) :: radsol
132    REAL, DIMENSION(klon) :: cdragq ! Cdrag pour l'evaporation
133    CHARACTER(len=20),PARAMETER :: modname="surf_ocean"
134    real rhoa(knon) ! density of moist air  (kg / m3)
135    REAL sens_prec_liq(knon)
136
137    real sss(klon)
138    ! Bulk salinity of the surface layer of the ocean, in ppt. (Only
139    ! defined for subscripts 1:knon, but we have to declare it with
140    ! size klon because of the coupling machinery.)
141
142    ! End definition
143    !******************************************************************************
144
145
146    !******************************************************************************
147    ! Calculate total net radiance at surface
148    !
149    !******************************************************************************
150    radsol(1:klon) = 0.0 ! initialisation a priori inutile
151    radsol(1:knon) = swnet(1:knon) + lwnet(1:knon)
152
153    !******************************************************************************
154    ! Cdragq computed from cdrag
155    ! The difference comes only from a factor (f_z0qh_oce) on z0, so that
156    ! it can be computed inside surf_ocean
157    ! More complicated appraches may require the propagation through
158    ! pbl_surface of an independant cdragq variable.
159    !******************************************************************************
160
161    IF ( f_z0qh_oce .ne. 1.) THEN
162       ! Si on suit les formulations par exemple de Tessel, on
163       ! a z0h=0.4*nu/u*, z0q=0.62*nu/u*, d'ou f_z0qh_oce=0.62/0.4=1.55
164       cdragq(1:knon)=cdragh(1:knon)*                                      &
165            log(z1lay(1:knon)/z0h(1:knon))/log(z1lay(1:knon)/(f_z0qh_oce*z0h(1:knon)))
166    ELSE
167       cdragq(1:knon)=cdragh(1:knon)
168    ENDIF
169
170    rhoa = PS(:KNON) / (Rd * temp_air(:knon) * (1. + retv * spechum(:knon)))
171    !******************************************************************************
172    ! Switch according to type of ocean (couple, slab or forced)
173    !******************************************************************************
174    SELECT CASE(type_ocean)
175    CASE('couple')
176       CALL ocean_cpl_noice( &
177            swnet, lwnet, alb1, &
178            windsp, fder, &
179            itime, dtime, knon, knindex, &
180            p1lay, cdragh, cdragq, cdragm, precip_rain, precip_snow,temp_air,spechum,&
181            AcoefH, AcoefQ, BcoefH, BcoefQ, &
182            AcoefU, AcoefV, BcoefU, BcoefV, &
183            ps, u1, v1, gustiness, tsurf_in, &
184            radsol, snow, agesno, &
185            qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
186            tsurf_new, dflux_s, dflux_l, sens_prec_liq, sss, s_int, rhoa)
187
188    CASE('slab')
189       CALL ocean_slab_noice( &
190            itime, dtime, jour, knon, knindex, &
191            p1lay, cdragh, cdragq, cdragm, precip_rain, precip_snow, temp_air, spechum,&
192            AcoefH, AcoefQ, BcoefH, BcoefQ, &
193            AcoefU, AcoefV, BcoefU, BcoefV, &
194            ps, u1, v1, gustiness, tsurf_in, &
195            radsol, snow, &
196            qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
197            tsurf_new, dflux_s, dflux_l, lmt_bils)
198
199    CASE('force')
200       CALL ocean_forced_noice( &
201            itime, dtime, jour, knon, knindex, &
202            p1lay, cdragh, cdragq, cdragm, precip_rain, precip_snow, &
203            temp_air, spechum, &
204            AcoefH, AcoefQ, BcoefH, BcoefQ, &
205            AcoefU, AcoefV, BcoefU, BcoefV, &
206            ps, u1, v1, gustiness, tsurf_in, &
207            radsol, snow, agesno, &
208            qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
209            tsurf_new, dflux_s, dflux_l, sens_prec_liq, rhoa)
210    END SELECT
211
212    !******************************************************************************
213    ! fcodron: compute lmt_bils  forced case (same as wfbils_oce / 1.-contfracatm)
214    !******************************************************************************
215    IF (type_ocean.NE.'slab') THEN
216       lmt_bils(1:klon)=0.
217       DO i=1,knon
218          lmt_bils(knindex(i))=(swnet(i)+lwnet(i)+fluxsens(i)+fluxlat(i)) &
219               *pctsrf(knindex(i),is_oce)/(1.-zmasq(knindex(i)))
220       END DO
221    END IF
222
223    !******************************************************************************
224    ! Calculate ocean surface albedo
225    !******************************************************************************
226    !albedo SB >>>
227    IF (iflag_albedo==0) THEN
228       !--old parametrizations of ocean surface albedo
229       !
230       IF (iflag_cycle_diurne.GE.1) THEN
231          !
232          CALL alboc_cd(rmu0,alb_eau)
233          !
234          !--ad-hoc correction for model radiative balance tuning
235          !--now outside alboc_cd routine
236          alb_eau(1:klon) = fmagic*alb_eau(1:klon) + pmagic
237          alb_eau(1:klon)=MIN(MAX(alb_eau(1:klon),0.0),1.0)
238          !
239       ELSE
240          !
241          CALL alboc(REAL(jour),rlat,alb_eau)
242          !--ad-hoc correction for model radiative balance tuning
243          !--now outside alboc routine
244          alb_eau(1:klon) = fmagic*alb_eau(1:klon) + pmagic
245          alb_eau(1:klon)=MIN(MAX(alb_eau(1:klon),0.04),0.60)
246          !
247       ENDIF
248       !
249       DO i =1, knon
250          DO  k=1,nsw
251             alb_dir_new(i,k) = alb_eau(knindex(i))
252          ENDDO
253       ENDDO
254       !IM 09122015 next line corresponds to the old way of doing in LMDZ5A/IPSLCM5A versions
255       !albedo for diffuse radiation is taken the same as for direct radiation
256       alb_dif_new(1:knon,:)=alb_dir_new(1:knon,:)
257       !IM 09122015 end
258       !
259    ELSE IF (iflag_albedo==1) THEN
260       !--new parametrization of ocean surface albedo by Sunghye Baek
261       !--albedo for direct and diffuse radiation are different
262       !
263       CALL ocean_albedo(knon,rmu0,knindex,windsp,SFRWL,alb_dir_new,alb_dif_new)
264       !
265       !--ad-hoc correction for model radiative balance tuning
266       alb_dir_new(1:knon,:) = fmagic*alb_dir_new(1:knon,:) + pmagic
267       alb_dif_new(1:knon,:) = fmagic*alb_dif_new(1:knon,:) + pmagic
268       alb_dir_new(1:knon,:)=MIN(MAX(alb_dir_new(1:knon,:),0.0),1.0)
269       alb_dif_new(1:knon,:)=MIN(MAX(alb_dif_new(1:knon,:),0.0),1.0)
270       !
271    ELSE IF (iflag_albedo==2) THEN
272       ! F. Codron albedo read from limit.nc
273       CALL limit_read_rug_alb(itime, dtime, jour,&
274            knon, knindex, z0_lim, alb_eau)
275       DO i =1, knon
276          DO  k=1,nsw
277             alb_dir_new(i,k) = alb_eau(i)
278          ENDDO
279       ENDDO
280       alb_dif_new=alb_dir_new
281    ENDIF
282    !albedo SB <<<
283
284    !******************************************************************************
285    ! Calculate the rugosity
286    !******************************************************************************
287    IF (iflag_z0_oce==0) THEN
288       DO i = 1, knon
289          tmp = MAX(cepdu2,gustiness(i)+u1(i)**2+v1(i)**2)
290          z0m(i) = 0.018*cdragm(i) * (gustiness(i)+u1(i)**2+v1(i)**2)/RG  &
291               +  0.11*14e-6 / SQRT(cdragm(i) * tmp)
292          z0m(i) = MAX(1.5e-05,z0m(i))
293       ENDDO
294       z0h(1:knon)=z0m(1:knon) ! En attendant mieux
295
296    ELSE IF (iflag_z0_oce==1) THEN
297       DO i = 1, knon
298          tmp = MAX(cepdu2,gustiness(i)+u1(i)**2+v1(i)**2)
299          z0m(i) = 0.018*cdragm(i) * (gustiness(i)+u1(i)**2+v1(i)**2)/RG  &
300               + 0.11*14e-6 / SQRT(cdragm(i) * tmp)
301          z0m(i) = MAX(1.5e-05,z0m(i))
302          z0h(i)=0.4*14e-6 / SQRT(cdragm(i) * tmp)
303       ENDDO
304    ELSE IF (iflag_z0_oce==-1) THEN
305       DO i = 1, knon
306          z0m(i) = z0min
307          z0h(i) = z0min
308       ENDDO
309    ELSE
310       CALL abort_physic(modname,'version non prevue',1)
311    ENDIF
312
313    if (activate_ocean_skin >= 1) then
314       if (type_ocean /= 'couple') sss(:knon) = 35.
315       call bulk_flux(tkt, tks, taur, dter, dser, t_int, s_int, ds_ns, dt_ns, &
316            u = windsp(:knon), t_ocean_1 = tsurf_new(:knon), s1 = sss(:knon), &
317            rain = precip_rain(:knon) + precip_snow(:knon), &
318            hf = - fluxsens(:knon), hlb = - fluxlat(:knon), &
319            rnl = - lwnet(:knon), &
320            tau = sqrt(flux_u1(:knon)**2 + flux_v1(:knon)**2), rhoa = rhoa, &
321            xlv = [(rlvtt, i = 1, knon)], rf = - sens_prec_liq, dtime = dtime, &
322            rns = swnet(:knon))
323       if (activate_ocean_skin == 2) tsurf_new(:knon) = t_int
324    end if
325   
326  END SUBROUTINE surf_ocean
327  !****************************************************************************
328  !
329END MODULE surf_ocean_mod
Note: See TracBrowser for help on using the repository browser.