Index: LMDZ6/trunk/libf/phylmdiso/phys_state_var_mod.F90
===================================================================
--- LMDZ6/trunk/libf/phylmdiso/phys_state_var_mod.F90	(revision 3927)
+++ LMDZ6/trunk/libf/phylmdiso/phys_state_var_mod.F90	(revision 3940)
@@ -1,4 +1,4 @@
 !
-! $Id: phys_state_var_mod.F90 3496 2019-05-10 10:17:35Z jyg $
+! $Id: phys_state_var_mod.F90 3888 2021-05-05 10:50:37Z jyg $
 !
       MODULE phys_state_var_mod
@@ -32,4 +32,6 @@
       REAL, ALLOCATABLE, SAVE :: ftsol(:,:)
 !$OMP THREADPRIVATE(ftsol)
+      REAL, ALLOCATABLE, SAVE :: beta_aridity(:,:)
+!$OMP THREADPRIVATE(beta_aridity)
       REAL,ALLOCATABLE,SAVE :: qsol(:),fevap(:,:),z0m(:,:),z0h(:,:),agesno(:,:)
 !$OMP THREADPRIVATE(qsol,fevap,z0m,z0h,agesno)
@@ -47,5 +49,5 @@
 !albedo SB >>>
       REAL, ALLOCATABLE, SAVE :: falb_dif(:,:,:), falb_dir(:,:,:)
-      real, allocatable, save :: chl_con(:)
+      REAL, ALLOCATABLE, SAVE :: chl_con(:)
 !$OMP THREADPRIVATE(falb_dir,falb_dif,chl_con)
 !albedo SB <<<
@@ -106,8 +108,4 @@
       REAL, ALLOCATABLE, SAVE :: coefm(:,:,:) ! Kz momentum
 !$OMP THREADPRIVATE(pbl_tke, coefh,coefm)
-!nrlmd<
-      REAL, ALLOCATABLE, SAVE :: delta_tsurf(:,:) ! Surface temperature difference inside-outside cold pool
-!$OMP THREADPRIVATE(delta_tsurf)
-!>nrlmd
       REAL, ALLOCATABLE, SAVE :: zmax0(:), f0(:) ! 
 !$OMP THREADPRIVATE(zmax0,f0)
@@ -294,4 +292,8 @@
       REAL,ALLOCATABLE,SAVE :: wake_delta_pbl_TKE(:,:,:)
 !$OMP THREADPRIVATE(wake_delta_pbl_TKE)
+!nrlmd<
+      REAL, ALLOCATABLE, SAVE :: delta_tsurf(:,:) ! Surface temperature difference inside-outside cold pool
+!$OMP THREADPRIVATE(delta_tsurf)
+!>nrlmd
 !>jyg
 !
@@ -442,8 +444,9 @@
 !$OMP THREADPRIVATE(ccm)
 
-!!! nrlmd le 10/04/2012
       REAL,SAVE,ALLOCATABLE :: ale_bl_trig(:)
 !$OMP THREADPRIVATE(ale_bl_trig)
-!!! fin nrlmd le 10/04/2012
+
+      REAL,SAVE,ALLOCATABLE :: ratqs_inter(:,:)
+!$OMP THREADPRIVATE(ratqs_inter)
 
 #ifdef ISO
@@ -462,5 +465,26 @@
 !$OMP THREADPRIVATE(is_initialized)    
 
-CONTAINS
+      ! Ocean-atmosphere interface:
+
+      REAL, ALLOCATABLE, SAVE:: ds_ns(:) ! (klon)
+      ! "delta salinity near surface". Salinity variation in the
+      ! near-surface turbulent layer. That is subskin salinity minus
+      ! foundation salinity. In ppt.
+
+      REAL, ALLOCATABLE, SAVE:: dt_ns(:) ! (klon)
+      ! "delta temperature near surface". Temperature variation in the
+      ! near-surface turbulent layer. That is subskin temperature
+      ! minus foundation temperature. (Can be negative.) In K.
+      
+      REAL, ALLOCATABLE, SAVE:: delta_sst(:) ! (klon)
+      ! Ocean-air interface temperature minus bulk SST, in
+      ! K. Allocated and defined only if activate_ocean_skin >= 1.
+
+      REAL, ALLOCATABLE, SAVE:: delta_sal(:) ! (klon)
+      ! Ocean-air interface salinity minus bulk salinity, in ppt
+      
+      !$OMP THREADPRIVATE(delta_sal, ds_ns, dt_ns, delta_sst)
+
+    CONTAINS
 
 !======================================================================
@@ -473,4 +497,5 @@
 #endif
 USE indice_sol_mod
+use config_ocean_skin_m, only: activate_ocean_skin
 IMPLICIT NONE
 
@@ -485,8 +510,10 @@
 include "clesphys.h"
 
+      print*, 'is_initialized', is_initialized
       IF (is_initialized) RETURN
       is_initialized=.TRUE.
       ALLOCATE(pctsrf(klon,nbsrf))
       ALLOCATE(ftsol(klon,nbsrf))
+      ALLOCATE(beta_aridity(klon,nbsrf))
       ALLOCATE(qsol(klon),fevap(klon,nbsrf))
       ALLOCATE(z0m(klon,nbsrf+1),z0h(klon,nbsrf+1),agesno(klon,nbsrf))
@@ -496,5 +523,7 @@
       ALLOCATE(falb2(klon,nbsrf))
 !albedo SB >>> 
+      print*, 'allocate falb'
       ALLOCATE(falb_dir(klon,nsw,nbsrf),falb_dif(klon,nsw,nbsrf))
+      print*, 'allocate falb good', falb_dir(1,1,1)
       ALLOCATE(chl_con(klon))
 !albedo SB <<<
@@ -675,9 +704,9 @@
 #endif      
 #endif
-!!! nrlmd le 10/04/2012
+
       ALLOCATE(ale_bl_trig(klon))
-!!! fin nrlmd le 10/04/2012
+      ALLOCATE(ratqs_inter(klon,klev))
       IF (ok_gwd_rando) THEN
-        allocate(du_gwd_rando(klon, klev))
+        ALLOCATE(du_gwd_rando(klon, klev))
         du_gwd_rando(:,:)=0.
       ENDIF
@@ -686,147 +715,155 @@
         du_gwd_front(:,:) = 0 !ym missing init   
       ENDIF
-END SUBROUTINE phys_state_var_init
+      if (activate_ocean_skin >= 1) ALLOCATE(delta_sal(klon), ds_ns(klon), &
+           dt_ns(klon), delta_sst(klon))
+
+    END SUBROUTINE phys_state_var_init
 
 !======================================================================
-SUBROUTINE phys_state_var_end
+    SUBROUTINE phys_state_var_end
+      ! Useful only for lmdz1d.
 !USE dimphy
 USE indice_sol_mod
+use config_ocean_skin_m, only: activate_ocean_skin
 IMPLICIT NONE
 include "clesphys.h"
 
-      deallocate(pctsrf, ftsol, falb1, falb2)
-      deallocate(qsol,fevap,z0m,z0h,agesno)
+      DEALLOCATE(pctsrf, ftsol, falb1, falb2)
+      DEALLOCATE(beta_aridity)
+      DEALLOCATE(qsol,fevap,z0m,z0h,agesno)
 !FC
-      deallocate(treedrg)
-      deallocate(rain_fall, snow_fall, solsw, solswfdiff, sollw, radsol, swradcorr)
-      deallocate(zmea, zstd, zsig, zgam)
-      deallocate(zthe, zpic, zval)
-      deallocate(rugoro, t_ancien, q_ancien, clwcon, rnebcon)
-      deallocate(qs_ancien, ql_ancien)
-      deallocate(prw_ancien, prlw_ancien, prsw_ancien)
-      deallocate(qtc_cv,sigt_cv)
-      deallocate(u_ancien, v_ancien)
-      deallocate(tr_ancien)                           !RomP
-      deallocate(ratqs, pbl_tke,coefh,coefm)
+      DEALLOCATE(treedrg)
+      DEALLOCATE(rain_fall, snow_fall, solsw, solswfdiff, sollw, radsol, swradcorr)
+      DEALLOCATE(zmea, zstd, zsig, zgam)
+      DEALLOCATE(zthe, zpic, zval)
+      DEALLOCATE(rugoro, t_ancien, q_ancien, clwcon, rnebcon)
+      DEALLOCATE(qs_ancien, ql_ancien)
+      DEALLOCATE(prw_ancien, prlw_ancien, prsw_ancien)
+      DEALLOCATE(qtc_cv,sigt_cv)
+      DEALLOCATE(u_ancien, v_ancien)
+      DEALLOCATE(tr_ancien)                           !RomP
+      DEALLOCATE(ratqs, pbl_tke,coefh,coefm)
+      DEALLOCATE(zmax0, f0)
+      DEALLOCATE(sig1, w01)
+      DEALLOCATE(entr_therm, fm_therm)
+      DEALLOCATE(detr_therm)
+      DEALLOCATE(clwcon0th, rnebcon0th)
+! radiation outputs
+      DEALLOCATE(swdnc0, swdn0, swdn)
+      DEALLOCATE(swupc0, swup0, swup)
+      DEALLOCATE(lwdnc0, lwdn0, lwdn)
+      DEALLOCATE(lwupc0, lwup0, lwup)
+      DEALLOCATE(SWdn200clr, SWdn200)
+      DEALLOCATE(SWup200clr, SWup200)
+      DEALLOCATE(LWdn200clr, LWdn200)
+      DEALLOCATE(LWup200clr, LWup200)
+      DEALLOCATE(LWdnTOA, LWdnTOAclr)
+! pressure level
+      DEALLOCATE(tsumSTD)
+      DEALLOCATE(usumSTD, vsumSTD)
+      DEALLOCATE(wsumSTD, phisumSTD)
+      DEALLOCATE(tnondef)
+      DEALLOCATE(qsumSTD, rhsumSTD)
+      DEALLOCATE(uvsumSTD)
+      DEALLOCATE(vqsumSTD)
+      DEALLOCATE(vTsumSTD)
+      DEALLOCATE(wqsumSTD)
+      DEALLOCATE(vphisumSTD)
+      DEALLOCATE(wTsumSTD)
+      DEALLOCATE(u2sumSTD)
+      DEALLOCATE(v2sumSTD)
+      DEALLOCATE(T2sumSTD)
+      DEALLOCATE(O3sumSTD)
+      DEALLOCATE(O3daysumSTD)
+!IM beg
+      DEALLOCATE(wlevSTD,ulevSTD,vlevSTD,tlevSTD,qlevSTD,rhlevSTD,philevSTD)
+      DEALLOCATE(uvSTD,vqSTD,vTSTD,wqSTD,vphiSTD,wTSTD,u2STD,v2STD,T2STD,O3STD,O3daySTD)
+!IM end
+      DEALLOCATE(seed_old)
+      DEALLOCATE(zuthe, zvthe)
+      DEALLOCATE(alb_neig)
+      DEALLOCATE(ema_cbmf)
+      DEALLOCATE(ema_pcb, ema_pct)
+      DEALLOCATE(Mipsh, Ma, qcondc)
+      DEALLOCATE(wd, sigd)
+      DEALLOCATE(cin, ALE, ALP)
+      DEALLOCATE(ftd, fqd)
+      DEALLOCATE(Ale_bl, Alp_bl)
+      DEALLOCATE(ale_wake)
+      DEALLOCATE(ale_bl_stat)
+      DEALLOCATE(lalim_conv, wght_th)
+      DEALLOCATE(wake_deltat, wake_deltaq)
+      DEALLOCATE(wake_s, awake_dens, wake_dens)
+      DEALLOCATE(wake_Cstar, wake_pe, wake_fip)
+!jyg<
+      DEALLOCATE(wake_delta_pbl_TKE)
 !nrlmd<
-      deallocate(delta_tsurf)
+      DEALLOCATE(delta_tsurf)
 !>nrlmd
-      deallocate(zmax0, f0)
-      deallocate(sig1, w01)
-      deallocate(entr_therm, fm_therm)
-      deallocate(detr_therm)
-      deallocate(clwcon0th, rnebcon0th)
-! radiation outputs
-      deallocate(swdnc0, swdn0, swdn)
-      deallocate(swupc0, swup0, swup)
-      deallocate(lwdnc0, lwdn0, lwdn)
-      deallocate(lwupc0, lwup0, lwup)
-      deallocate(SWdn200clr, SWdn200)
-      deallocate(SWup200clr, SWup200)
-      deallocate(LWdn200clr, LWdn200)
-      deallocate(LWup200clr, LWup200)
-      deallocate(LWdnTOA, LWdnTOAclr)
-! pressure level
-      deallocate(tsumSTD)
-      deallocate(usumSTD, vsumSTD)
-      deallocate(wsumSTD, phisumSTD)
-      deallocate(tnondef)
-      deallocate(qsumSTD, rhsumSTD)
-      deallocate(uvsumSTD)
-      deallocate(vqsumSTD)
-      deallocate(vTsumSTD)
-      deallocate(wqsumSTD)
-      deallocate(vphisumSTD)
-      deallocate(wTsumSTD)
-      deallocate(u2sumSTD)
-      deallocate(v2sumSTD)
-      deallocate(T2sumSTD)
-      deallocate(O3sumSTD)
-      deallocate(O3daysumSTD)
-!IM beg
-      deallocate(wlevSTD,ulevSTD,vlevSTD,tlevSTD,qlevSTD,rhlevSTD,philevSTD)
-      deallocate(uvSTD,vqSTD,vTSTD,wqSTD,vphiSTD,wTSTD,u2STD,v2STD,T2STD,O3STD,O3daySTD)
-!IM end
-      deallocate(seed_old)
-      deallocate(zuthe, zvthe)
-      deallocate(alb_neig)
-      deallocate(ema_cbmf)
-      deallocate(ema_pcb, ema_pct)
-      deallocate(Mipsh, Ma, qcondc)
-      deallocate(wd, sigd)
-      deallocate(cin, ALE, ALP)
-      deallocate(ftd, fqd)
-      deallocate(Ale_bl, Alp_bl)
-      deallocate(ale_wake)
-      deallocate(ale_bl_stat)
-      deallocate(lalim_conv, wght_th)
-      deallocate(wake_deltat, wake_deltaq)
-      deallocate(wake_s, awake_dens, wake_dens)
-      deallocate(wake_Cstar, wake_pe, wake_fip)
-!jyg<
-      deallocate(wake_delta_pbl_TKE)
 !>jyg
-      deallocate(pfrac_impa, pfrac_nucl)
-      deallocate(pfrac_1nucl)
-      deallocate(total_rain, nday_rain)
-      deallocate(paire_ter)
-      deallocate(albsol1, albsol2)
+      DEALLOCATE(pfrac_impa, pfrac_nucl)
+      DEALLOCATE(pfrac_1nucl)
+      DEALLOCATE(total_rain, nday_rain)
+      DEALLOCATE(paire_ter)
+      DEALLOCATE(albsol1, albsol2)
 !albedo SB >>>
-      deallocate(albsol_dir,albsol_dif,falb_dir,falb_dif,chl_con)
+      DEALLOCATE(albsol_dir,albsol_dif,falb_dir,falb_dif,chl_con)
 !albedo SB <<<
-      deallocate(wo)
-      deallocate(clwcon0,rnebcon0)
-      deallocate(heat, heat0) 
-      deallocate(cool, cool0)
-      deallocate(heat_volc, cool_volc) 
-      deallocate(topsw, toplw)
-      deallocate(sollwdown, sollwdownclr)
-      deallocate(gustiness)
-      deallocate(toplwdown, toplwdownclr)
-      deallocate(topsw0,toplw0,solsw0,sollw0)
-      deallocate(albpla)
+      DEALLOCATE(wo)
+      DEALLOCATE(clwcon0,rnebcon0)
+      DEALLOCATE(heat, heat0) 
+      DEALLOCATE(cool, cool0)
+      DEALLOCATE(heat_volc, cool_volc) 
+      DEALLOCATE(topsw, toplw)
+      DEALLOCATE(sollwdown, sollwdownclr)
+      DEALLOCATE(gustiness)
+      DEALLOCATE(toplwdown, toplwdownclr)
+      DEALLOCATE(topsw0,toplw0,solsw0,sollw0)
+      DEALLOCATE(albpla)
 !IM ajout variables CFMIP2/CMIP5
-      deallocate(heatp, coolp)
-      deallocate(heat0p, cool0p)
-      deallocate(radsolp, topswp, toplwp)
-      deallocate(albplap)
-      deallocate(solswp, solswfdiffp, sollwp)
-      deallocate(sollwdownp)
-      deallocate(topsw0p,toplw0p)
-      deallocate(solsw0p,sollw0p)
-      deallocate(lwdnc0p, lwdn0p, lwdnp)
-      deallocate(lwupc0p, lwup0p, lwupp)
-      deallocate(swdnc0p, swdn0p, swdnp)
-      deallocate(swupc0p, swup0p, swupp)
-      deallocate(cape)
-      deallocate(pbase,bbase)
-      deallocate(zqasc)
-      deallocate(ibas_con, itop_con)
-      deallocate(rain_con, snow_con)
-      deallocate(rlonPOS)
-      deallocate(newsst)
-      deallocate(ustar,u10m, v10m,wstar)
-      deallocate(topswad, solswad)
-      deallocate(topswai, solswai)
-      deallocate(tau_aero,piz_aero,cg_aero)
-      deallocate(tau_aero_sw_rrtm,piz_aero_sw_rrtm,cg_aero_sw_rrtm)
-      deallocate(tau_aero_lw_rrtm,piz_aero_lw_rrtm,cg_aero_lw_rrtm)
-      deallocate(ccm)
-      if (ok_gwd_rando) deallocate(du_gwd_rando)
-      if (.not. ok_hines .and. ok_gwd_rando) deallocate(du_gwd_front)
-       
-!!! nrlmd le 10/04/2012
-      deallocate(ale_bl_trig)
-!!! fin nrlmd le 10/04/2012
+      DEALLOCATE(heatp, coolp)
+      DEALLOCATE(heat0p, cool0p)
+      DEALLOCATE(radsolp, topswp, toplwp)
+      DEALLOCATE(albplap)
+      DEALLOCATE(solswp, solswfdiffp, sollwp)
+      DEALLOCATE(sollwdownp)
+      DEALLOCATE(topsw0p,toplw0p)
+      DEALLOCATE(solsw0p,sollw0p)
+      DEALLOCATE(lwdnc0p, lwdn0p, lwdnp)
+      DEALLOCATE(lwupc0p, lwup0p, lwupp)
+      DEALLOCATE(swdnc0p, swdn0p, swdnp)
+      DEALLOCATE(swupc0p, swup0p, swupp)
+      DEALLOCATE(cape)
+      DEALLOCATE(pbase,bbase)
+      DEALLOCATE(zqasc)
+      DEALLOCATE(ibas_con, itop_con)
+      DEALLOCATE(rain_con, snow_con)
+      DEALLOCATE(rlonPOS)
+      DEALLOCATE(newsst)
+      DEALLOCATE(ustar,u10m, v10m,wstar)
+      DEALLOCATE(topswad, solswad)
+      DEALLOCATE(topswai, solswai)
+      DEALLOCATE(tau_aero,piz_aero,cg_aero)
+      DEALLOCATE(tau_aero_sw_rrtm,piz_aero_sw_rrtm,cg_aero_sw_rrtm)
+      DEALLOCATE(tau_aero_lw_rrtm,piz_aero_lw_rrtm,cg_aero_lw_rrtm)
+      DEALLOCATE(ccm)
+      if (ok_gwd_rando) DEALLOCATE(du_gwd_rando)
+      if (.not. ok_hines .and. ok_gwd_rando) DEALLOCATE(du_gwd_front)
+      DEALLOCATE(ale_bl_trig)
+      DEALLOCATE(ratqs_inter)
+
+      if (activate_ocean_skin >= 1) deALLOCATE(delta_sal, ds_ns, dt_ns, &
+           delta_sst)
 
 #ifdef ISO    
-      deallocate(xtsol,fxtevap)  
-      deallocate(xt_ancien,xtl_ancien,xts_ancien, fxtd, wake_deltaxt)
-      deallocate(xtrain_fall, xtsnow_fall, xtrain_con, xtsnow_con)
+      DEALLOCATE(xtsol,fxtevap)  
+      DEALLOCATE(xt_ancien,xtl_ancien,xts_ancien, fxtd, wake_deltaxt)
+      DEALLOCATE(xtrain_fall, xtsnow_fall, xtrain_con, xtsnow_con)
 #ifdef ISOTRAC 
-      deallocate(bassin_map,boite_map) 
+      DEALLOCATE(bassin_map,boite_map) 
 #endif        
 #endif
       is_initialized=.FALSE.
+      
 END SUBROUTINE phys_state_var_end
 
