Changeset 3940 for LMDZ6


Ignore:
Timestamp:
Jun 15, 2021, 1:18:14 PM (3 years ago)
Author:
crisi
Message:

replace files by symbloic liks from phylmdiso towards phylmd.
Many files at once

Location:
LMDZ6/trunk/libf/phylmdiso
Files:
271 added
17 edited

Legend:

Unmodified
Added
Removed
  • LMDZ6/trunk/libf/phylmdiso/fonte_neige_mod.F90

    r3927 r3940  
    2828  REAL, PRIVATE                               :: tau_calv 
    2929  !$OMP THREADPRIVATE(tau_calv)
    30   REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE  :: ffonte_global
     30  REAL, ALLOCATABLE, DIMENSION(:,:)           :: ffonte_global
    3131  !$OMP THREADPRIVATE(ffonte_global)
    32   REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE  :: fqfonte_global
     32  REAL, ALLOCATABLE, DIMENSION(:,:)           :: fqfonte_global
    3333  !$OMP THREADPRIVATE(fqfonte_global)
    34   REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE  :: fqcalving_global
     34  REAL, ALLOCATABLE, DIMENSION(:,:)           :: fqcalving_global
    3535  !$OMP THREADPRIVATE(fqcalving_global)
    36   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE  :: runofflic_global
     36  REAL, ALLOCATABLE, DIMENSION(:)             :: runofflic_global
    3737  !$OMP THREADPRIVATE(runofflic_global)
    3838#ifdef ISO
     
    6060! The variable run_off_lic_0 is initialized to the field read from
    6161! restart file. The other variables are initialized to zero.
    62 
    6362!
    6463!****************************************************************************************
     
    8483    run_off_lic_0(:) = restart_runoff(:)
    8584
    86 
    8785!****************************************************************************************
    8886! Allocate other variables and initilize to zero
     
    130128    ENDIF
    131129    runofflic_global(:) = 0.0
    132 
    133 
    134130
    135131!****************************************************************************************
     
    348344
    349345    snow_evap = 0.
    350    
     346 
    351347#ifdef ISOVERIF
    352348        write(*,*) 'klon,snow_evap(413)=',klon,snow_evap(413)
     
    406402       ! Y'a-t-il fonte de neige?
    407403       neige_fond = (snow(i)>epsfra .OR. nisurf==is_sic .OR. nisurf==is_lic) .AND. tsurf_new(i)>=RTT
    408        IF (neige_fond) THEN   
     404       IF (neige_fond) THEN
    409405          fq_fonte     = MIN( MAX((tsurf_new(i)-RTT )/chasno,0.0),snow(i))
    410406          ffonte(i)    = fq_fonte * RLMLT/dtime
     
    412408          snow(i)      = MAX(0., snow(i) - fq_fonte)
    413409          bil_eau_s(i) = bil_eau_s(i) + fq_fonte
    414           tsurf_new(i) = tsurf_new(i) - fq_fonte * chasno   
     410          tsurf_new(i) = tsurf_new(i) - fq_fonte * chasno 
    415411#ifdef ISO
    416412        fq_fonte_diag(i)=fq_fonte
  • LMDZ6/trunk/libf/phylmdiso/infotrac_phy.F90

    r3927 r3940  
    2323  INTEGER, SAVE :: nqtottr
    2424!$OMP THREADPRIVATE(nqtottr)
     25
     26! ThL : number of CO2 tracers                   ModThL
     27  INTEGER, SAVE :: nqCO2
     28!$OMP THREADPRIVATE(nqCO2)
    2529
    2630#ifdef CPP_StratAer
     
    194198    ALLOCATE(solsym(nbtr))
    195199    solsym(:)=solsym_(:)
    196  
     200     
    197201    IF(prt_level.ge.1) THEN
    198       write(lunout,*) TRIM(modname)//": nqtot,nqo,nbtr",nqtot,nqo,nbtr
     202      write(lunout,*) TRIM(modname)//": nqtot,nqo,nbtr,nqCO2",nqtot,nqo,nbtr,nqCO2
    199203    ENDIF
    200204   
  • LMDZ6/trunk/libf/phylmdiso/ocean_forced_mod.F90

    r3927 r3940  
    1919       AcoefH, AcoefQ, BcoefH, BcoefQ, &
    2020       AcoefU, AcoefV, BcoefU, BcoefV, &
    21        ps, u1, v1, gustiness, &
     21       ps, u1, v1, gustiness, tsurf_in, &
    2222       radsol, snow, agesno, &
    2323       qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    24        tsurf_new, dflux_s, dflux_l &
     24       tsurf_new, dflux_s, dflux_l, sens_prec_liq, rhoa &
    2525#ifdef ISO
    2626            ,xtprecip_rain, xtprecip_snow, xtspechum,Roce,rlat, &
     
    4040    USE indice_sol_mod
    4141    USE phys_output_var_mod, ONLY : sens_prec_liq_o, sens_prec_sol_o, lat_prec_liq_o, lat_prec_sol_o
     42    use config_ocean_skin_m, only: activate_ocean_skin
    4243#ifdef ISO
    4344  USE infotrac_phy, ONLY: ntraciso,niso
     
    5354    INCLUDE "YOMCST.h"
    5455    INCLUDE "clesphys.h"
    55 
     56    INCLUDE "flux_arp.h"
    5657
    5758! Input arguments
     
    6869    REAL, DIMENSION(klon), INTENT(IN)        :: ps
    6970    REAL, DIMENSION(klon), INTENT(IN)        :: u1, v1, gustiness
     71    REAL, DIMENSION(klon), INTENT(IN)        :: tsurf_in
     72    real, intent(in):: rhoa(:) ! (knon) density of moist air  (kg / m3)
     73
    7074#ifdef ISO
    7175    REAL, DIMENSION(ntraciso,klon), INTENT(IN)    :: xtprecip_rain, xtprecip_snow
     
    9195    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    9296    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l
     97    REAL, intent(out):: sens_prec_liq(:) ! (knon)
     98
    9399#ifdef ISO     
    94100    REAL, DIMENSION(ntraciso,klon), INTENT(OUT)    :: xtevap ! isotopes in evaporation flux
     
    104110    REAL, DIMENSION(klon)       :: u1_lay, v1_lay
    105111    LOGICAL                     :: check=.FALSE.
    106     REAL, DIMENSION(klon) :: sens_prec_liq, sens_prec_sol   
     112    REAL sens_prec_sol(knon)
    107113    REAL, DIMENSION(klon) :: lat_prec_liq, lat_prec_sol   
    108114
     
    139145!!jyg    if (knon.eq.1) then ! single-column model
    140146    if (klon_glo.eq.1) then ! single-column model
    141       CALL read_tsurf1d(knon,tsurf_lim) ! new
    142 #ifdef ISO
    143       write(*,*) 'ocean_forced_mod 143: isotopes pas prévus ici'
    144       stop
    145 #endif
     147      ! EV: now surface Tin flux_arp.h
     148      !CALL read_tsurf1d(knon,tsurf_lim) ! new
     149       DO i = 1, knon
     150        tsurf_lim(i) = tg
     151       ENDDO
     152
    146153    else ! GCM
    147154      CALL limit_read_sst(knon,knindex,tsurf_lim &
     
    159166!****************************************************************************************
    160167! Set some variables for calcul_fluxs
    161     cal = 0.
    162     beta = 1.
    163     dif_grnd = 0.
     168    !cal = 0.
     169    !beta = 1.
     170    !dif_grnd = 0.
     171   
     172   
     173    ! EV: use calbeta to calculate beta
     174    ! Need to initialize qsurf for calbeta but it is not modified by this routine
     175    qsurf(:)=0.
     176    CALL calbeta(dtime, is_oce, knon, snow, qsurf, beta, cal, dif_grnd)
     177
     178
    164179    alb_neig(:) = 0.
    165180    agesno(:) = 0.
    166     sens_prec_liq = 0.; sens_prec_sol = 0.; lat_prec_liq = 0.; lat_prec_sol = 0.
     181    lat_prec_liq = 0.; lat_prec_sol = 0.
    167182
    168183! Suppose zero surface speed
     
    174189! Calcul de tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l and qsurf
    175190    CALL calcul_fluxs(knon, is_oce, dtime, &
    176          tsurf_lim, p1lay, cal, beta, cdragh, cdragq, ps, &
     191         merge(tsurf_in, tsurf_lim, activate_ocean_skin == 2), p1lay, cal, &
     192         beta, cdragh, cdragq, ps, &
    177193         precip_rain, precip_snow, snow, qsurf,  &
    178194         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, gustiness, &
    179195         f_qsat_oce,AcoefH, AcoefQ, BcoefH, BcoefQ, &
    180196         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
    181          sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol)
     197         sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol, rhoa)
     198    if (activate_ocean_skin == 2) tsurf_new = tsurf_lim
    182199
    183200    do j = 1, knon
     
    233250       radsol, snow, qsol, agesno, tsoil, &
    234251       qsurf, alb1_new, alb2_new, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    235        tsurf_new, dflux_s, dflux_l &
     252       tsurf_new, dflux_s, dflux_l, rhoa &
    236253#ifdef ISO
    237254            ,xtprecip_rain, xtprecip_snow, xtspechum,Roce, &
     
    262279#endif
    263280
    264 !    INCLUDE "indicesol.h"
     281!   INCLUDE "indicesol.h"
    265282    INCLUDE "dimsoil.h"
    266283    INCLUDE "YOMCST.h"
    267284    INCLUDE "clesphys.h"
     285    INCLUDE "flux_arp.h"
    268286
    269287! Input arguments
     
    281299    REAL, DIMENSION(klon), INTENT(IN)    :: ps
    282300    REAL, DIMENSION(klon), INTENT(IN)    :: u1, v1, gustiness
     301    real, intent(in):: rhoa(:) ! (knon) density of moist air  (kg / m3)
    283302#ifdef ISO
    284303    REAL, DIMENSION(ntraciso,klon), INTENT(IN)    :: xtprecip_rain, xtprecip_snow
     
    323342    REAL, DIMENSION(klon)       :: u0, v0
    324343    REAL, DIMENSION(klon)       :: u1_lay, v1_lay
    325     REAL, DIMENSION(klon)       :: sens_prec_liq, sens_prec_sol   
     344    REAL sens_prec_liq(knon), sens_prec_sol (knon)
    326345    REAL, DIMENSION(klon)       :: lat_prec_liq, lat_prec_sol   
    327346
     
    354373    tsurf_tmp(:) = tsurf_in(:)
    355374
    356 ! calculate the parameters cal, beta, capsol and dif_grnd
     375! calculate the parameters cal, beta, capsol and dif_grnd and then recalculate cal
    357376    CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, capsol, dif_grnd)
    358377
     
    370389    ENDIF
    371390
    372     beta = 1.0
    373     sens_prec_liq = 0.; sens_prec_sol = 0.; lat_prec_liq = 0.; lat_prec_sol = 0.
     391!    beta = 1.0
     392    lat_prec_liq = 0.; lat_prec_sol = 0.
    374393
    375394! Suppose zero surface speed
     
    384403         f_qsat_oce,AcoefH, AcoefQ, BcoefH, BcoefQ, &
    385404         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
    386          sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol)
     405         sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol, rhoa)
    387406    do j = 1, knon
    388407      i = knindex(j)
     
    489508! 1D case
    490509!************************************************************************
    491   SUBROUTINE read_tsurf1d(knon,sst_out)
    492 
     510!  SUBROUTINE read_tsurf1d(knon,sst_out)
     511!
    493512! This subroutine specifies the surface temperature to be used in 1D simulations
    494 
    495       USE dimphy, ONLY : klon
    496 
    497       INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
    498       REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out  ! tsurf used to force the single-column model
    499 
    500        INTEGER :: i
     513!
     514!      USE dimphy, ONLY : klon
     515!
     516!      INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
     517!      REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out  ! tsurf used to force the single-column model
     518!
     519!       INTEGER :: i
    501520! COMMON defined in lmdz1d.F:
    502        real ts_cur
    503        common /sst_forcing/ts_cur
    504 
    505        DO i = 1, knon
    506         sst_out(i) = ts_cur
    507        ENDDO
    508 
    509       END SUBROUTINE read_tsurf1d
    510 
     521!       real ts_cur
     522!       common /sst_forcing/ts_cur
     523!
     524!       DO i = 1, knon
     525!        sst_out(i) = ts_cur
     526!       ENDDO
     527!
     528!      END SUBROUTINE read_tsurf1d
     529!
    511530!
    512531!************************************************************************
    513 !
    514532END MODULE ocean_forced_mod
    515533
  • LMDZ6/trunk/libf/phylmdiso/pbl_surface_mod.F90

    r3927 r3940  
    11!
    2 ! $Id: pbl_surface_mod.F90 3435 2019-01-22 15:21:59Z fairhead $
     2! $Id: pbl_surface_mod.F90 3906 2021-05-19 10:35:18Z jyg $
    33!
    44MODULE pbl_surface_mod
     
    2323  USE climb_wind_mod,      ONLY : climb_wind_down, climb_wind_up
    2424  USE coef_diff_turb_mod,  ONLY : coef_diff_turb
    25   USE wx_pbl_mod,          ONLY : wx_pbl_init, wx_pbl_final, &
    26 !!                                  wx_pbl_fuse_no_dts, wx_pbl_split_no_dts, &
    27 !!                                  wx_pbl_fuse, wx_pbl_split
    28                                   wx_pbl0_fuse, wx_pbl0_split
     25  USE ioipsl_getin_p_mod,  ONLY : getin_p
     26  USE cdrag_mod
     27  USE stdlevvar_mod
     28  USE wx_pbl_var_mod,      ONLY : wx_pbl_init, wx_pbl_final, &
     29                                  wx_pbl_prelim_0, wx_pbl_prelim_beta
     30  USE wx_pbl_mod,          ONLY : wx_pbl0_merge, wx_pbl_split, wx_pbl_dts_merge, &
     31                                  wx_pbl_check, wx_pbl_dts_check, wx_evappot
     32  use config_ocean_skin_m, only: activate_ocean_skin
    2933
    3034  IMPLICIT NONE
     
    3337  REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE     :: fder   ! flux drift
    3438  !$OMP THREADPRIVATE(fder)
    35   REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC, SAVE   :: snow   ! snow at surface
     39  REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC, SAVE    :: snow   ! snow at surface
    3640  !$OMP THREADPRIVATE(snow)
    3741  REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE   :: qsurf  ! humidity at surface
    3842  !$OMP THREADPRIVATE(qsurf)
    39   REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: ftsoil ! soil temperature
     43  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE          :: ftsoil ! soil temperature
    4044  !$OMP THREADPRIVATE(ftsoil)
     45  REAL, ALLOCATABLE, DIMENSION(:), SAVE              :: ydTs0, ydqs0 
     46                                                     ! nul forced temperature and humidity differences
     47  !$OMP THREADPRIVATE(ydTs0, ydqs0)
    4148
    4249#ifdef ISO
     
    5158  INTEGER, SAVE :: iflag_pbl_surface_t2m_bug
    5259  !$OMP THREADPRIVATE(iflag_pbl_surface_t2m_bug)
     60  INTEGER, SAVE :: iflag_new_t2mq2m
     61  !$OMP THREADPRIVATE(iflag_new_t2mq2m)
     62
    5363!FC
    5464!  integer, save :: iflag_frein
     
    7888    REAL, DIMENSION(klon, nbsrf), INTENT(IN)          :: qsurf_rst
    7989    REAL, DIMENSION(klon, nsoilmx, nbsrf), INTENT(IN) :: ftsoil_rst
    80 
    8190 
    8291! Local variables
     
    102111    IF (ierr /= 0) CALL abort_physic('pbl_surface_init', 'pb in allocation',1)
    103112
     113    ALLOCATE(ydTs0(klon), stat=ierr)
     114    IF (ierr /= 0) CALL abort_physic('pbl_surface_init', 'pb in allocation',1)
     115
     116    ALLOCATE(ydqs0(klon), stat=ierr)
     117    IF (ierr /= 0) CALL abort_physic('pbl_surface_init', 'pb in allocation',1)
     118
    104119    fder(:)       = fder_rst(:)
    105120    snow(:,:)     = snow_rst(:,:)
    106121    qsurf(:,:)    = qsurf_rst(:,:)
    107122    ftsoil(:,:,:) = ftsoil_rst(:,:,:)
    108 
     123    ydTs0(:) = 0.
     124    ydqs0(:) = 0.
    109125
    110126!****************************************************************************************
     
    152168    iflag_pbl_surface_t2m_bug=0
    153169    CALL getin_p('iflag_pbl_surface_t2m_bug',iflag_pbl_surface_t2m_bug)
     170    WRITE(lunout,*) 'iflag_pbl_surface_t2m_bug=',iflag_pbl_surface_t2m_bug
    154171!FC
    155172!    iflag_frein = 0
     
    240257       debut,     lafin,                              &
    241258       rlon,      rlat,      rugoro,   rmu0,          &
    242        zsig,      lwdown_m,  pphi,     cldt,          &
    243        rain_f,    snow_f,    solsw_m,  sollw_m,       &
     259       lwdown_m,  cldt,          &
     260       rain_f,    snow_f,    solsw_m,  solswfdiff_m, sollw_m,       &
    244261       gustiness,                                     &
    245262       t,         q,         u,        v,             &
     
    252269       ts,SFRWL,   alb_dir, alb_dif,ustar, u10m, v10m,wstar, &
    253270       cdragh,    cdragm,   zu1,    zv1,              &
     271!jyg<   (26/09/2019)
     272       beta, &
     273!>jyg
    254274       alb_dir_m,    alb_dif_m,  zxsens,   zxevap,    &
    255275       alb3_lic,  runoff,    snowhgt,   qsnow,     to_ice,    sissnow,  &
    256        zxtsol,    zxfluxlat, zt2m,     qsat2m,       &
     276       zxtsol,    zxfluxlat, zt2m,     qsat2m, zn2mout, &
    257277       d_t,       d_q,       d_u,      d_v, d_t_diss, &
    258278!!! nrlmd+jyg le 02/05/2011 et le 20/02/2012
     
    275295       s_therm,   s_trmb1,   s_trmb2,  s_trmb3,       &
    276296       zustar,zu10m,  zv10m,    fder_print,    &
    277        zxqsurf,   rh2m,      zxfluxu,  zxfluxv,       &
     297       zxqsurf, delta_qsurf,                       &
     298       rh2m,      zxfluxu,  zxfluxv,               &
    278299       z0m, z0h,   agesno,  sollw,    solsw,         &
    279300       d_ts,      evap,    fluxlat,  t2m,           &
     
    283304!jyg<
    284305!!       zxfluxt,   zxfluxq,   q2m,      flux_q, tke,   &
    285        zxfluxt,   zxfluxq,   q2m,      flux_q, tke_x,   &
     306       zxfluxt,   zxfluxq,   q2m,      flux_q, tke_x,  &
    286307!>jyg
    287308!!! nrlmd+jyg le 02/05/2011 et le 20/02/2012
     
    338359! z0m, z0h ----input-R- longeur de rugosite (en m)
    339360! Martin
    340 ! zsig-----input-R- slope
    341361! cldt-----input-R- total cloud fraction
    342 ! pphi-----input-R- geopotentiel de chaque couche (g z) (reference sol)
    343362! Martin
    344363!
     
    370389    USE carbon_cycle_mod,   ONLY : carbon_cycle_cpl, carbon_cycle_tr, level_coupling_esm
    371390    USE carbon_cycle_mod,   ONLY : co2_send, nbcf_out, fields_out, yfields_out, cfname_out
     391    use hbtm_mod, only: hbtm
    372392    USE indice_sol_mod
    373393    USE time_phylmdz_mod,   ONLY : day_ini,annee_ref,itau_phy
     
    385405#endif
    386406    USE ioipsl_getin_p_mod, ONLY : getin_p
     407    use phys_state_var_mod, only: ds_ns, dt_ns, delta_sst, delta_sal, zsig, zmea
     408    use phys_output_var_mod, only: dter, dser, tkt, tks, taur, sss
     409#ifdef CPP_XIOS
     410    USE wxios, ONLY: missing_val
     411#else
     412    use netcdf, only: missing_val => nf90_fill_real
     413#endif
     414
     415     
     416
    387417
    388418    IMPLICIT NONE
     
    412442    REAL, DIMENSION(klon),        INTENT(IN)        :: snow_f  ! snow fall
    413443    REAL, DIMENSION(klon),        INTENT(IN)        :: solsw_m ! net shortwave radiation at mean surface
     444    REAL, DIMENSION(klon),        INTENT(IN)        :: solswfdiff_m ! diffuse fraction fordownward shortwave radiation at mean surface
    414445    REAL, DIMENSION(klon),        INTENT(IN)        :: sollw_m ! net longwave radiation at mean surface
    415446    REAL, DIMENSION(klon,klev),   INTENT(IN)        :: t       ! temperature (K)
     
    421452    REAL, DIMENSION(klon, nbsrf), INTENT(IN)        :: pctsrf  ! sub-surface fraction
    422453! Martin
    423     REAL, DIMENSION(klon),        INTENT(IN)        :: zsig    ! slope
    424454    REAL, DIMENSION(klon),        INTENT(IN)        :: lwdown_m ! downward longwave radiation at mean s   
    425455    REAL, DIMENSION(klon),        INTENT(IN)        :: gustiness ! gustiness
    426456
    427457    REAL, DIMENSION(klon),        INTENT(IN)        :: cldt    ! total cloud fraction
    428     REAL, DIMENSION(klon,klev),   INTENT(IN)        :: pphi    ! geopotential (m2/s2)
    429 ! Martin
    430458
    431459#ifdef ISO
     
    451479! Input/Output variables
    452480!****************************************************************************************
     481!jyg<
     482    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: beta    ! Aridity factor
     483!>jyg
    453484    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: ts      ! temperature at surface (K)
    454485    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: delta_tsurf !surface temperature difference between
     
    496527    REAL, DIMENSION(klon),        INTENT(OUT)       :: zxfluxlat  ! latent flux, mean for each grid point
    497528    REAL, DIMENSION(klon),        INTENT(OUT)       :: zt2m       ! temperature at 2m, mean for each grid point
     529    INTEGER, DIMENSION(klon, 6),  INTENT(OUT)       :: zn2mout    ! number of times the 2m temperature is out of the [tsol,temp]
    498530    REAL, DIMENSION(klon),        INTENT(OUT)       :: qsat2m
    499531    REAL, DIMENSION(klon, klev),  INTENT(OUT)       :: d_t        ! change in temperature
     
    560592    REAL, DIMENSION(klon),        INTENT(OUT)       :: fder_print ! fder for printing (=fder(i) + dflux_t(i) + dflux_q(i))
    561593    REAL, DIMENSION(klon),        INTENT(OUT)       :: zxqsurf    ! humidity at surface, mean for each grid point
     594    REAL, DIMENSION(klon),        INTENT(OUT)       :: delta_qsurf! humidity difference at surface, mean for each grid point
    562595    REAL, DIMENSION(klon),        INTENT(OUT)       :: rh2m       ! relative humidity at 2m
    563596    REAL, DIMENSION(klon, klev),  INTENT(OUT)       :: zxfluxu    ! u wind tension, mean for each grid point
    564597    REAL, DIMENSION(klon, klev),  INTENT(OUT)       :: zxfluxv    ! v wind tension, mean for each grid point
    565     REAL, DIMENSION(klon, nbsrf+1), INTENT(INOUT)     :: z0m,z0h      ! rugosity length (m)
    566     REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)       :: agesno   ! age of snow at surface
     598    REAL, DIMENSION(klon, nbsrf+1), INTENT(INOUT)   :: z0m,z0h      ! rugosity length (m)
     599    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: agesno   ! age of snow at surface
    567600    REAL, DIMENSION(klon, nbsrf), INTENT(OUT)       :: solsw      ! net shortwave radiation at surface
    568601    REAL, DIMENSION(klon, nbsrf), INTENT(OUT)       :: sollw      ! net longwave radiation at surface
    569602    REAL, DIMENSION(klon, nbsrf), INTENT(OUT)       :: d_ts       ! change in temperature at surface
    570     REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)       :: evap     ! evaporation at surface
     603    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: evap     ! evaporation at surface
    571604    REAL, DIMENSION(klon, nbsrf), INTENT(OUT)       :: fluxlat    ! latent flux
    572605    REAL, DIMENSION(klon, nbsrf), INTENT(OUT)       :: t2m        ! temperature at 2 meter height
     
    605638
    606639! Martin
    607 ! sisvat
     640! inlandsis
    608641    REAL, DIMENSION(klon),       INTENT(OUT)        :: qsnow      ! snow water content
    609642    REAL, DIMENSION(klon),       INTENT(OUT)        :: snowhgt    ! snow height
     
    632665    INTEGER                            :: n
    633666! << PC
    634     INTEGER                            :: iflag_split
     667    INTEGER                            :: iflag_split, iflag_split_ref
    635668    INTEGER                            :: i, k, nsrf
    636669    INTEGER                            :: knon, j
     
    643676    REAL, DIMENSION(klon)              :: r_co2_ppm     ! taux CO2 atmosphere
    644677    REAL, DIMENSION(klon)              :: yts, yz0m, yz0h, ypct
     678    REAL, DIMENSION(klon)              :: yz0h_old
    645679!albedo SB >>>
    646680    REAL, DIMENSION(klon)              :: yalb,yalb_vis
    647681!albedo SB <<<
    648682    REAL, DIMENSION(klon)              :: yt1, yq1, yu1, yv1
     683    REAL, DIMENSION(klon)              :: yqa
    649684    REAL, DIMENSION(klon)              :: ysnow, yqsurf, yagesno, yqsol
    650685    REAL, DIMENSION(klon)              :: yrain_f, ysnow_f
     
    670705    REAL, DIMENSION(klon)              :: y_flux_u1, y_flux_v1
    671706    REAL, DIMENSION(klon)              :: yt2m, yq2m, yu10m
     707    INTEGER, DIMENSION(klon, nbsrf, 6) :: yn2mout, yn2mout_x, yn2mout_w
     708    INTEGER, DIMENSION(klon, nbsrf, 6) :: n2mout, n2mout_x, n2mout_w
    672709    REAL, DIMENSION(klon)              :: yustar
    673710    REAL, DIMENSION(klon)              :: ywstar
     
    690727    REAL, DIMENSION(klon)              :: yz0h_oupas
    691728    REAL, DIMENSION(klon)              :: yfluxsens
     729    REAL, DIMENSION(klon)              :: AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0
    692730    REAL, DIMENSION(klon)              :: AcoefH, AcoefQ, BcoefH, BcoefQ
    693731#ifdef ISO
     
    696734    REAL, DIMENSION(klon)              :: AcoefU, AcoefV, BcoefU, BcoefV
    697735    REAL, DIMENSION(klon)              :: ypsref
    698     REAL, DIMENSION(klon)              :: yevap, ytsurf_new, yalb3_new
     736    REAL, DIMENSION(klon)              :: yevap, yevap_pot, ytsurf_new, yalb3_new
    699737!albedo SB >>>
    700738    REAL, DIMENSION(klon,nsw)          :: yalb_dir_new, yalb_dif_new
     
    708746    REAL, DIMENSION(klon,klev)         :: y_flux_u, y_flux_v
    709747    REAL, DIMENSION(klon,klev)         :: ycoefh, ycoefm,ycoefq
    710     REAL, DIMENSION(klon)              :: ycdragh, ycdragm
     748    REAL, DIMENSION(klon)              :: ycdragh, ycdragq, ycdragm
    711749    REAL, DIMENSION(klon,klev)         :: yu, yv
    712750    REAL, DIMENSION(klon,klev)         :: yt, yq
     
    746784    REAL, DIMENSION(klon,klev)         :: ycoefh_x, ycoefm_x, ycoefh_w, ycoefm_w
    747785    REAL, DIMENSION(klon,klev)         :: ycoefq_x, ycoefq_w
    748     REAL, DIMENSION(klon)              :: ycdragh_x, ycdragm_x, ycdragh_w, ycdragm_w
     786    REAL, DIMENSION(klon)              :: ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w
     787    REAL, DIMENSION(klon)              :: ycdragm_x, ycdragm_w
    749788    REAL, DIMENSION(klon)              :: AcoefH_x, AcoefQ_x, BcoefH_x, BcoefQ_x
    750789    REAL, DIMENSION(klon)              :: AcoefH_w, AcoefQ_w, BcoefH_w, BcoefQ_w
     
    767806    REAL                               :: zx_qs_surf, zcor_surf, zdelta_surf
    768807    REAL, DIMENSION(klon)              :: ytsurf_th, yqsatsurf
     808!jyg<
    769809    REAL, DIMENSION(klon)              :: ybeta
     810    REAL, DIMENSION(klon)              :: ybeta_prev
     811!>jyg
    770812    REAL, DIMENSION(klon, klev)        :: d_u_x
    771813    REAL, DIMENSION(klon, klev)        :: d_u_w
     
    915957!!! nrlmd le 13/06/2011
    916958    REAL, DIMENSION(klon)              :: y_delta_flux_t1, y_delta_flux_q1, y_delta_flux_u1, y_delta_flux_v1
    917     REAL, DIMENSION(klon)              :: y_delta_tsurf,delta_coef,tau_eq
     959    REAL, DIMENSION(klon)              :: y_delta_tsurf, y_delta_tsurf_new
     960    REAL, DIMENSION(klon)              :: delta_coef, tau_eq
     961    REAL, DIMENSION(klon)              :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn
     962    REAL, DIMENSION(klon)              :: phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0
     963    REAL, DIMENSION(klon)              :: y_delta_qsurf
     964    REAL, DIMENSION(klon)              :: y_delta_qsats
     965    REAL, DIMENSION(klon)              :: yg_T, yg_Q
     966    REAL, DIMENSION(klon)              :: yGamma_dTs_phiT, yGamma_dQs_phiQ
     967    REAL, DIMENSION(klon)              :: ydTs_ins, ydqs_ins
     968!
    918969    REAL, PARAMETER                    :: facteur=2./sqrt(3.14)
    919970    REAL, PARAMETER                    :: inertia=2000.
     
    928979    REAL, DIMENSION(klon)              :: Kech_m
    929980    REAL, DIMENSION(klon)              :: Kech_m_x, Kech_m_w
    930     REAL, DIMENSION(klon)              :: yts_x,yts_w
     981    REAL, DIMENSION(klon)              :: yts_x, yts_w
     982    REAL, DIMENSION(klon)              :: yqsatsrf0_x, yqsatsrf0_w
     983    REAL, DIMENSION(klon)              :: yqsurf_x, yqsurf_w
    931984!jyg<
    932985!!    REAL, DIMENSION(klon)              :: Kech_Hp, Kech_H_xp, Kech_H_wp
     
    935988!!    REAL, DIMENSION(klon)              :: Kech_Vp, Kech_V_xp, Kech_V_wp
    936989!>jyg
    937 !jyg<
    938     REAL, DIMENSION(klon)              :: ah, bh     ! coefficients of the delta_Tsurf equation
    939 !>jyg
     990
     991    REAL                               :: fact_cdrag
     992    REAL                               :: z1lay
    940993
    941994    REAL                               :: vent
     
    9711024    REAL, DIMENSION(klon)              :: ytoice
    9721025    REAL, DIMENSION(klon)              :: ysnowhgt, yqsnow, ysissnow, yrunoff
     1026    REAL, DIMENSION(klon)              :: yzmea
    9731027    REAL, DIMENSION(klon)              :: yzsig
    974     REAL, DIMENSION(klon,klev)         :: ypphi
    9751028    REAL, DIMENSION(klon)              :: ycldt
    9761029    REAL, DIMENSION(klon)              :: yrmu0
    9771030    ! Martin
     1031
     1032    REAL, DIMENSION(klon):: ydelta_sst, ydelta_sal, yds_ns, ydt_ns, ydter, ydser, &
     1033         ytkt, ytks, ytaur, ysss
     1034    ! compression of delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, tkt, tks,
     1035    ! taur, sss on ocean points
     1036
    9781037#ifdef ISO
    9791038    REAL, DIMENSION(klon)       :: h1
     
    9911050!
    9921051!!jyg      iflag_split = mod(iflag_pbl_split,2)
    993       iflag_split = mod(iflag_pbl_split,10)
     1052!!jyg      iflag_split = mod(iflag_pbl_split,10)
     1053      iflag_split_ref = mod(iflag_pbl_split,10)
    9941054
    9951055#ifdef ISO     
     
    10381098
    10391099    IF (first_call) THEN
     1100
     1101       iflag_new_t2mq2m=1
     1102       CALL getin_p('iflag_new_t2mq2m',iflag_new_t2mq2m)
     1103       WRITE(lunout,*) 'pbl_iflag_new_t2mq2m=',iflag_new_t2mq2m
     1104
    10401105       print*,'PBL SURFACE AVEC GUSTINESS'
    10411106       first_call=.FALSE.
    10421107     
    10431108       ! Initialize ok_flux_surf (for 1D model)
    1044        if (klon_glo>1) ok_flux_surf=.FALSE.
     1109       IF (klon_glo>1) ok_flux_surf=.FALSE.
     1110       IF (klon_glo>1) ok_forc_tsurf=.FALSE.
    10451111
    10461112       ! intialize beta_land
     
    11131179 zxfluxlat(:)=0.
    11141180 zt2m(:)=0. ; zq2m(:)=0. ; qsat2m(:)=0. ; rh2m(:)=0.
     1181 zn2mout(:,:)=0 ;
    11151182 d_t(:,:)=0. ; d_t_diss(:,:)=0. ; d_q(:,:)=0. ; d_u(:,:)=0. ; d_v(:,:)=0.
    11161183 zcoefh(:,:,:)=0. ; zcoefm(:,:,:)=0.
     
    11281195 fder_print(:)=0.
    11291196 zxqsurf(:)=0.
     1197 delta_qsurf(:) = 0.
    11301198 zxfluxu(:,:)=0. ; zxfluxv(:,:)=0.
    11311199 solsw(:,:)=0. ; sollw(:,:)=0.
     
    11641232!!    tke(:,:,is_ave)=0.
    11651233    tke_x(:,:,is_ave)=0.
     1234
    11661235    wake_dltke(:,:,is_ave)=0.
    11671236!>jyg
     
    11841253    yqsurf = 0.0  ; yalb = 0.0 ; yalb_vis = 0.0
    11851254!albedo SB <<<
    1186     yrain_f = 0.0 ; ysnow_f = 0.0    ; yfder = 0.0     ; ysolsw = 0.0   
     1255    yrain_f = 0.0 ; ysnow_f = 0.0    ; yfder = 0.0     ; ysolsw = 0.0
    11871256    ysollw = 0.0  ; yz0m = 0.0 ; yz0h = 0.0    ; yu1 = 0.0   
    11881257    yv1 = 0.0     ; ypaprs = 0.0     ; ypplay = 0.0
     
    11951264!!    d_t_diss= 0.0 ;d_u = 0.0     ; d_v = 0.0
    11961265    yqsol = 0.0   
    1197     ytherm = 0.0  ; ytke=0.
     1266
     1267    ytke=0.
    11981268!FC
    11991269    y_treedrg=0.
     
    12021272    ysnowhgt = 0.0; yqsnow = 0.0     ; yrunoff = 0.0   ; ytoice =0.0
    12031273    yalb3_new = 0.0  ; ysissnow = 0.0
    1204     ypphi = 0.0   ; ycldt = 0.0      ; yrmu0 = 0.0
     1274    ycldt = 0.0      ; yrmu0 = 0.0
    12051275    ! Martin
    12061276
     
    12181288    y_delta_flux_t1=0.
    12191289    ydtsurf_th=0.
    1220     yts_x=0.      ; yts_w=0.
    1221     y_delta_tsurf=0.
     1290    yts_x(:)=0.      ; yts_w(:)=0.
     1291    y_delta_tsurf(:)=0. ; y_delta_qsurf(:)=0.
     1292    yqsurf_x(:)=0.      ; yqsurf_w(:)=0.
     1293    yg_T(:) = 0. ;        yg_Q(:) = 0.
     1294    yGamma_dTs_phiT(:) = 0. ; yGamma_dQs_phiQ(:) = 0.
     1295    ydTs_ins(:) = 0. ; ydqs_ins(:) = 0.
     1296
    12221297!!!
    12231298    ytsoil = 999999.
     
    14101485       DO i = 1, klon
    14111486          sollw(i,nsrf) = sollw_m(i) + 4.0*RSIGMA*ztsol(i)**3 * (ztsol(i)-ts(i,nsrf))
    1412 
    1413 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    1414 !         ! Martin
    1415 ! Apparently introduced for sisvat but not used
    1416 !         sollwd(i,nsrf)= sollwd_m(i)
    1417 !         ! Martin
    1418 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    1419 
     1487!--OB this line is not satisfactory because alb is the direct albedo not total albedo
    14201488          solsw(i,nsrf) = solsw_m(i) * (1.-alb(i,nsrf)) / (1.-alb_m(i))
    14211489       ENDDO
     
    14401508!>al1
    14411509
     1510!--OB add diffuse fraction of SW down
     1511   DO n=1,nbcf_out
     1512       IF (cfname_out(n) == "swdownfdiff" ) fields_out(:,n) = solswfdiff_m(:)
     1513   ENDDO
    14421514! >> PC
    14431515   IF (carbon_cycle_cpl .AND. carbon_cycle_tr .AND. nbcf_out.GT.0 ) THEN
     
    14611533!
    14621534!****************************************************************************************
    1463    
    1464     loop_nbsrf: DO nsrf = 1, nbsrf
     1535                                                                          !<<<<<<<<<<<<<
     1536    loop_nbsrf: DO nsrf = 1, nbsrf                                        !<<<<<<<<<<<<<
     1537                                                                          !<<<<<<<<<<<<<
    14651538       IF (prt_level >=10) print *,' Loop nsrf ',nsrf
     1539!
     1540       IF (iflag_split_ref == 3) THEN
     1541         IF (nsrf == is_oce) THEN
     1542            iflag_split = 1
     1543         ELSE
     1544            iflag_split=0
     1545         ENDIF   !! (nsrf == is_oce)
     1546       ELSE                     
     1547         iflag_split = iflag_split_ref
     1548       ENDIF   !! (iflag_split_ref == 3)
    14661549
    14671550! Search for index(ni) and size(knon) of domaine to treat
     
    14991582!****************************************************************************************
    15001583
     1584!
     1585!jyg<    (20190926)
     1586!   Provisional : set ybeta to standard values
     1587       IF (nsrf .NE. is_ter) THEN
     1588           ybeta(:) = 1.
     1589       ELSE
     1590           IF (iflag_split .EQ. 0) THEN
     1591              ybeta(:) = 1.
     1592           ELSE
     1593             DO j = 1, knon
     1594                i = ni(j)
     1595                ybeta(j)   = beta(i,nsrf)
     1596             ENDDO
     1597           ENDIF  ! (iflag_split .LE.1)
     1598       ENDIF !  (nsrf .NE. is_ter)
     1599!>jyg
     1600!
    15011601       DO j = 1, knon
    15021602          i = ni(j)
     
    15311631          ywindsp(j) = windsp(i,nsrf)
    15321632!>jyg
    1533           ! Martin
     1633          ! Martin and Etienne
     1634          yzmea(j)   = zmea(i)
    15341635          yzsig(j)   = zsig(i)
    15351636          ycldt(j)   = cldt(i)
     
    16561757             ytke_w(j,k)      = tke_x(i,k,nsrf)+wake_dltke(i,k,nsrf)
    16571758             ywake_dltke(j,k) = wake_dltke(i,k,nsrf)
     1759           
    16581760!>jyg
    16591761          ENDDO
     
    16951797          ENDDO
    16961798       ENDIF
     1799
     1800       if (nsrf == is_oce .and. activate_ocean_skin >= 1) then
     1801          if (activate_ocean_skin == 2 .and. type_ocean == "couple") then
     1802             ydelta_sal(:knon) = delta_sal(ni(:knon))
     1803             ydelta_sst(:knon) = delta_sst(ni(:knon))
     1804          end if
     1805         
     1806          yds_ns(:knon) = ds_ns(ni(:knon))
     1807          ydt_ns(:knon) = dt_ns(ni(:knon))
     1808       end if
    16971809       
    16981810!****************************************************************************************
     
    17351847      ENDDO
    17361848     ENDIF
     1849
    17371850        IF (prt_level >=10) print *,'clcdrag -> ycdragh ', ycdragh
    17381851       ELSE  !(iflag_split .eq.0)
     
    17491862           speed_x(i) = SQRT(yu_x(i,1)**2+yv_x(i,1)**2)
    17501863        ENDDO
    1751         CALL cdrag(knon, nsrf, &
     1864
     1865
     1866            CALL cdrag(knon, nsrf, &
    17521867            speed_x, yt_x(:,1), yq_x(:,1), zgeo1_x, ypaprs(:,1),&
    1753             yts_x, yqsurf, yz0m, yz0h, &
     1868            yts_x, yqsurf_x, yz0m, yz0h, &
    17541869            ycdragm_x, ycdragh_x, zri1_x, pref_x )
    17551870
     
    17781893        CALL cdrag(knon, nsrf, &
    17791894            speed_w, yt_w(:,1), yq_w(:,1), zgeo1_w, ypaprs(:,1),&
    1780             yts_w, yqsurf, yz0m, yz0h, &
     1895            yts_w, yqsurf_w, yz0m, yz0h, &
    17811896            ycdragm_w, ycdragh_w, zri1_w, pref_w )
    17821897!
     
    18181933      print *,' args coef_diff_turb: ycdragh ', ycdragh
    18191934      print *,' args coef_diff_turb: ytke ', ytke
     1935
    18201936       ENDIF
    18211937        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
     
    18471963      print *,' args coef_diff_turb: ycdragh_x ', ycdragh_x
    18481964      print *,' args coef_diff_turb: ytke_x ', ytke_x
     1965
    18491966       ENDIF
    18501967        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    1851             ypaprs, ypplay, yu_x, yv_x, yq_x, yt_x, yts_x, yqsurf, ycdragm_x, &
     1968            ypaprs, ypplay, yu_x, yv_x, yq_x, yt_x, yts_x, yqsurf_x, ycdragm_x, &
    18521969            ycoefm_x, ycoefh_x, ytke_x,y_treedrg)
    18531970!            ycoefm_x, ycoefh_x, ytke_x)
     
    18771994       ENDIF
    18781995        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    1879             ypaprs, ypplay, yu_w, yv_w, yq_w, yt_w, yts_w, yqsurf, ycdragm_w, &
     1996            ypaprs, ypplay, yu_w, yv_w, yq_w, yt_w, yts_w, yqsurf_w, ycdragm_w, &
    18801997            ycoefm_w, ycoefh_w, ytke_w,y_treedrg)
    18811998!            ycoefm_w, ycoefh_w, ytke_w)
     
    20292146         yxt1(:,:) = yxt(:,:,1)
    20302147#endif
    2031 !!       ELSE IF (iflag_split .eq. 1) THEN
    2032 !!!
    2033 !jyg<
    2034 !!         CALL wx_pbl_fuse_no_dts(knon, dtime, ypplay, ywake_s, &
    2035 !!                                 yt_x, yt_w, yq_x, yq_w, &
    2036 !!                                 yu_x, yu_w, yv_x, yv_w, &
    2037 !!                                 ycdragh_x, ycdragh_w, ycdragm_x, ycdragm_w, &
    2038 !!                                 AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, &
    2039 !!                                 AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, &
    2040 !!                                 BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, &
    2041 !!                                 BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, &
    2042 !!                                 AcoefH, AcoefQ, AcoefU, AcoefV, &
    2043 !!                                 BcoefH, BcoefQ, BcoefU, BcoefV, &
    2044 !!                                 ycdragh, ycdragm, &
    2045 !!                                 yt1, yq1, yu1, yv1 &
    2046 !!                                 )
     2148
    20472149       ELSE IF (iflag_split .ge. 1) THEN
    2048          CALL wx_pbl0_fuse(knon, dtime, ypplay, ywake_s, &
     2150#ifdef ISO
     2151        call abort_gcm('pbl_surface_mod 2149','isos pas encore dans iflag_split=1',1)
     2152#endif
     2153
     2154!
     2155! Cdragq computation
     2156! ------------------
     2157    !******************************************************************************
     2158    ! Cdragq computed from cdrag
     2159    ! The difference comes only from a factor (f_z0qh_oce) on z0, so that
     2160    ! it can be computed inside wx_pbl0_merge
     2161    ! More complicated appraches may require the propagation through
     2162    ! pbl_surface of an independant cdragq variable.
     2163    !******************************************************************************
     2164!
     2165    IF ( f_z0qh_oce .ne. 1. .and. nsrf .eq.is_oce) THEN
     2166       ! Si on suit les formulations par exemple de Tessel, on
     2167       ! a z0h=0.4*nu/u*, z0q=0.62*nu/u*, d'ou f_z0qh_oce=0.62/0.4=1.55
     2168!!       ycdragq_x(1:knon)=ycdragh_x(1:knon)*                                      &
     2169!!            log(z1lay(1:knon)/yz0h(1:knon))/log(z1lay(1:knon)/(f_z0qh_oce*yz0h(1:knon)))
     2170!!       ycdragq_w(1:knon)=ycdragh_w(1:knon)*                                      &
     2171!!            log(z1lay(1:knon)/yz0h(1:knon))/log(z1lay(1:knon)/(f_z0qh_oce*yz0h(1:knon)))
     2172!
     2173       DO j = 1,knon
     2174         z1lay = zgeo1(j)/RG
     2175         fact_cdrag = log(z1lay/yz0h(j))/log(z1lay/(f_z0qh_oce*yz0h(j)))
     2176         ycdragq_x(j)=ycdragh_x(j)*fact_cdrag
     2177         ycdragq_w(j)=ycdragh_w(j)*fact_cdrag
     2178!!     Print *,'YYYYpbl0: fact_cdrag ', fact_cdrag
     2179       ENDDO  ! j = 1,knon
     2180!
     2181!!  Print *,'YYYYpbl0: z1lay, yz0h, f_z0qh_oce, ycdragh_w, ycdragq_w ', &
     2182!!                z1lay, yz0h(1:knon), f_z0qh_oce, ycdragh_w(1:knon), ycdragq_w(1:knon)
     2183    ELSE
     2184       ycdragq_x(1:knon)=ycdragh_x(1:knon)
     2185       ycdragq_w(1:knon)=ycdragh_w(1:knon)
     2186    ENDIF  ! ( f_z0qh_oce .ne. 1. .and. nsrf .eq.is_oce)
     2187!
     2188         CALL wx_pbl_prelim_0(knon, nsrf, dtime, ypplay, ypaprs, ywake_s,  &
     2189                         yts, y_delta_tsurf, ygustiness, &
    20492190                         yt_x, yt_w, yq_x, yq_w, &
    20502191                         yu_x, yu_w, yv_x, yv_w, &
    2051                          ycdragh_x, ycdragh_w, ycdragm_x, ycdragm_w, &
     2192                         ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, &
     2193                         ycdragm_x, ycdragm_w, &
     2194                         AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, &
     2195                         AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, &
     2196                         BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, &
     2197                         BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w  &
     2198                         )
     2199         CALL wx_pbl_prelim_beta(knon, dtime, ywake_s, ybeta,  &
     2200                         BcoefQ_x, BcoefQ_w  &
     2201                         )
     2202         CALL wx_pbl0_merge(knon, ypplay, ypaprs,  &
     2203                         ywake_s, ydTs0, ydqs0, &
     2204                         yt_x, yt_w, yq_x, yq_w, &
     2205                         yu_x, yu_w, yv_x, yv_w, &
     2206                         ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, &
     2207                         ycdragm_x, ycdragm_w, &
    20522208                         AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, &
    20532209                         AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, &
    20542210                         BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, &
    20552211                         BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, &
    2056                          AcoefH, AcoefQ, AcoefU, AcoefV, &
    2057                          BcoefH, BcoefQ, BcoefU, BcoefV, &
    2058                          ycdragh, ycdragm, &
     2212                         AcoefH_0, AcoefQ_0, AcoefU, AcoefV, &
     2213                         BcoefH_0, BcoefQ_0, BcoefU, BcoefV, &
     2214                         ycdragh, ycdragq, ycdragm, &
    20592215                         yt1, yq1, yu1, yv1 &
    2060 #ifdef ISO
    2061                          ,yxt_x,yxt_w,yxt1 &
    2062 #endif
    20632216                         )
    2064 !!       ELSE IF (iflag_split .ge.2) THEN
    2065 !!!    Provisoire
    2066 !!         ah(:) = 0.
    2067 !!         bh(:) = 0.
    2068 !!         IF (nsrf == is_oce) THEN
    2069 !!           ybeta(:) = 1.
    2070 !!         ELSE
    2071 !!           ybeta(:) = beta_land
    2072 !!         ENDIF
    2073 !!         ycdragh(:) = ywake_s(:)*ycdragh_w(:) + (1.-ywake_s(:))*ycdragh_x(:)
    2074 !!         CALL wx_dts(knon, nsrf, ywake_cstar, ywake_s, ywake_dens, &
    2075 !!                     yts, ypplay(:,1), ybeta, ycdragh , ypaprs(:,1), &
    2076 !!                     yq(:,1), yt(:,1), yu(:,1), yv(:,1), ygustiness, &
    2077 !!                     ah, bh &
    2078 !!                     )
    2079 !!!
    2080 !!         CALL wx_pbl_fuse(knon, dtime, ypplay, ywake_s, &
    2081 !!                         yt_x, yt_w, yq_x, yq_w, &
    2082 !!                         yu_x, yu_w, yv_x, yv_w, &
    2083 !!                         ycdragh_x, ycdragh_w, ycdragm_x, ycdragm_w, &
    2084 !!                         AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, &
    2085 !!                         AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, &
    2086 !!                         BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, &
    2087 !!                         BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, &
    2088 !!                         ah, bh, &
    2089 !!                         AcoefH, AcoefQ, AcoefU, AcoefV, &
    2090 !!                         BcoefH, BcoefQ, BcoefU, BcoefV, &
    2091 !!                         ycdragh, ycdragm, &
    2092 !!                         yt1, yq1, yu1, yv1 &
    2093 !!                         )
    2094 !>jyg
    2095 !!!
    2096          ENDIF  ! (iflag_split .eq.0)
     2217         IF (iflag_split .eq. 2 .AND. nsrf .ne. is_oce) THEN
     2218           CALL wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, &
     2219                           ywake_s, ybeta, ywake_cstar, ywake_dens, &
     2220                           AcoefH_x, AcoefH_w, &
     2221                           BcoefH_x, BcoefH_w, &
     2222                           AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0,  &
     2223                           AcoefH, AcoefQ, BcoefH, BcoefQ,  &
     2224                           HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, &
     2225                           phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, &
     2226                           yg_T, yg_Q, &
     2227                           yGamma_dTs_phiT, yGamma_dQs_phiQ, &
     2228                           ydTs_ins, ydqs_ins &
     2229                           )
     2230         ELSE !
     2231           AcoefH(:) = AcoefH_0(:)
     2232           AcoefQ(:) = AcoefQ_0(:)
     2233           BcoefH(:) = BcoefH_0(:)
     2234           BcoefQ(:) = BcoefQ_0(:)
     2235           yg_T(:) = 0.
     2236           yg_Q(:) = 0.
     2237           yGamma_dTs_phiT(:) = 0.
     2238           yGamma_dQs_phiQ(:) = 0.
     2239           ydTs_ins(:) = 0.
     2240           ydqs_ins(:) = 0.
     2241         ENDIF   ! (iflag_split .eq. 2)
     2242       ENDIF  ! (iflag_split .eq.0)
    20972243!!!
    20982244       IF (prt_level >=10) THEN
    2099          PRINT *,'pbl_surface (fuse->): yt(1,:) ',yt(1,:)
    2100          PRINT *,'pbl_surface (fuse->): yq(1,:) ',yq(1,:)
    2101          PRINT *,'pbl_surface (fuse->): yu(1,:) ',yu(1,:)
    2102          PRINT *,'pbl_surface (fuse->): yv(1,:) ',yv(1,:)
    2103          PRINT *,'pbl_surface (fuse->): AcoefH(1) ',AcoefH(1)
    2104          PRINT *,'pbl_surface (fuse->): BcoefH(1) ',BcoefH(1)
     2245         PRINT *,'pbl_surface (merge->): yt(1,:) ',yt(1,:)
     2246         PRINT *,'pbl_surface (merge->): yq(1,:) ',yq(1,:)
     2247         PRINT *,'pbl_surface (merge->): yu(1,:) ',yu(1,:)
     2248         PRINT *,'pbl_surface (merge->): yv(1,:) ',yv(1,:)
     2249         PRINT *,'pbl_surface (merge->): AcoefH(1), AcoefQ(1), AcoefU(1), AcoefV(1) ', &
     2250                                         AcoefH(1), AcoefQ(1), AcoefU(1), AcoefV(1)
     2251         PRINT *,'pbl_surface (merge->): BcoefH(1), BcoefQ(1), BcoefU(1), BcoefV(1) ', &
     2252                                         BcoefH(1), BcoefQ(1), BcoefU(1), BcoefV(1)
     2253
    21052254       ENDIF
    21062255
     2256!  Save initial value of z0h for use in evappot (z0h wiil be computed again in the surface models)
     2257          yz0h_old(1:knon) = yz0h(1:knon)
     2258!
    21072259!****************************************************************************************
    21082260!
     
    21192271
    21202272          ! Calculate the temperature et relative humidity at 2m and the wind at 10m
     2273          IF (iflag_new_t2mq2m==1) THEN
     2274           CALL stdlevvarn(klon, knon, is_ter, zxli, &
     2275               yu(:,1), yv(:,1), yt(:,1), yq(:,1), zgeo1, &
     2276               yts, yqsurf, yz0m, yz0h, ypaprs(:,1), ypplay(:,1), &
     2277               yt2m, yq2m, yt10m, yq10m, yu10m, yustar, &
     2278               yn2mout(:, nsrf, :))
     2279          ELSE
    21212280          CALL stdlevvar(klon, knon, is_ter, zxli, &
    21222281               yu(:,1), yv(:,1), yt(:,1), yq(:,1), zgeo1, &
    21232282               yts, yqsurf, yz0m, yz0h, ypaprs(:,1), ypplay(:,1), &
    21242283               yt2m, yq2m, yt10m, yq10m, yu10m, yustar)
     2284          ENDIF
    21252285         
    21262286       ENDIF
     
    22112371          CALL surf_landice(itap, dtime, knon, ni, &
    22122372               rlon, rlat, debut, lafin, &
    2213                yrmu0, ylwdown, yalb, ypphi(:,1), &
     2373               yrmu0, ylwdown, yalb, zgeo1, &
    22142374               ysolsw, ysollw, yts, ypplay(:,1), &
    22152375!!jyg               ycdragh, ycdragm, yrain_f, ysnow_f, yt(:,1), yq(:,1),&
     
    22212381               ytsoil, yz0m, yz0h, SFRWL, yalb_dir_new, yalb_dif_new, yevap,yfluxsens,yfluxlat, &
    22222382               ytsurf_new, y_dflux_t, y_dflux_q, &
    2223                yzsig, ycldt, &
     2383               yzmea, yzsig, ycldt, &
    22242384               ysnowhgt, yqsnow, ytoice, ysissnow, &
    22252385               yalb3_new, yrunoff, &
     
    22842444               yz0m, yz0h, SFRWL,yalb_dir_new, yalb_dif_new, yevap, yfluxsens,yfluxlat,&
    22852445               ytsurf_new, y_dflux_t, y_dflux_q, slab_wfbils, &
    2286                y_flux_u1, y_flux_v1 &
     2446               y_flux_u1, y_flux_v1, ydelta_sst(:knon), ydelta_sal(:knon), &
     2447               yds_ns(:knon), ydt_ns(:knon), ydter(:knon), ydser(:knon), &
     2448               ytkt(:knon), ytks(:knon), ytaur(:knon), ysss &
    22872449#ifdef ISO
    22882450         &      ,yxtrain_f, yxtsnow_f,yxt1,Roce, &
     
    23922554!
    23932555!****************************************************************************************
    2394 
    2395 !!!
    2396 !!! jyg le 10/04/2013
     2556!!
     2557!!!
     2558!!! jyg le 10/04/2013 et EV 10/2020
     2559
     2560        IF (ok_forc_tsurf) THEN
     2561            DO j=1,knon
     2562                ytsurf_new(j)=tg
     2563                y_d_ts(j) = ytsurf_new(j) - yts(j)
     2564            ENDDO
     2565        ENDIF ! ok_forc_tsurf
     2566
    23972567!!!
    23982568        IF (ok_flux_surf) THEN
     
    24292599#endif
    24302600          ENDDO
    2431         ENDIF
    2432 
    2433        IF (prt_level >=10) THEN
    2434         DO j=1,knon
    2435          print*,'y_flux_t1,yfluxlat,wakes' &
    2436  &             ,  y_flux_t1(j), yfluxlat(j), ywake_s(j)
    2437          print*,'beta,ytsurf_new', ybeta(j), ytsurf_new(j)
    2438          print*,'inertia,facteur,cstar', inertia, facteur,wake_cstar(j)
    2439         ENDDO
    2440        ENDIF
    2441 
    2442 !!! jyg le 07/02/2012 puis le 10/04/2013
    2443 !!       IF (iflag_split .eq.1) THEN
    2444 !!!!!
    2445 !!!jyg<
    2446 !!         CALL wx_pbl_split_no_dts(knon, ywake_s, &
    2447 !!                                AcoefH_x, AcoefH_w, &
    2448 !!                                AcoefQ_x, AcoefQ_w, &
    2449 !!                                AcoefU_x, AcoefU_w, &
    2450 !!                                AcoefV_x, AcoefV_w, &
    2451 !!                                y_flux_t1, y_flux_q1, y_flux_u1, y_flux_v1, &
    2452 !!                                y_flux_t1_x, y_flux_t1_w, &
    2453 !!                                y_flux_q1_x, y_flux_q1_w, &
    2454 !!                                y_flux_u1_x, y_flux_u1_w, &
    2455 !!                                y_flux_v1_x, y_flux_v1_w, &
    2456 !!                                yfluxlat_x, yfluxlat_w &
    2457 !!                                )
    2458 !!       ELSE IF (iflag_split .ge. 2) THEN
     2601        ENDIF ! (ok_flux_surf)
     2602!
     2603! ------------------------------------------------------------------------------
     2604! 12a)  Splitting
     2605! ------------------------------------------------------------------------------
     2606
    24592607       IF (iflag_split .GE. 1) THEN
    2460          CALL wx_pbl0_split(knon, dtime, ywake_s, &
     2608#ifdef ISO
     2609        call abort_gcm('pbl_surface_mod 2607','isos pas encore dans iflag_split=1',1)
     2610#endif
     2611!
     2612         IF (nsrf .ne. is_oce) THEN
     2613!
     2614!         Compute potential evaporation and aridity factor  (jyg, 20200328)
     2615          ybeta_prev(:) = ybeta(:)
     2616             DO j = 1, knon
     2617               yqa(j) = AcoefQ(j) - BcoefQ(j)*yevap(j)*dtime
     2618             ENDDO
     2619!
     2620          CALL wx_evappot(knon, yqa, yTsurf_new, yevap_pot)
     2621!
     2622          ybeta(1:knon) = min(yevap(1:knon)/yevap_pot(1:knon), 1.)
     2623         
     2624          IF (prt_level >=10) THEN
     2625           DO j=1,knon
     2626            print*,'y_flux_t1,yfluxlat,wakes' &
     2627 &                ,  y_flux_t1(j), yfluxlat(j), ywake_s(j)
     2628            print*,'beta_prev, beta, ytsurf_new', ybeta_prev(j), ybeta(j), ytsurf_new(j)
     2629            print*,'inertia,facteur,cstar', inertia, facteur,wake_cstar(j)
     2630           ENDDO
     2631          ENDIF  ! (prt_level >=10)
     2632!
     2633! Second call to wx_pbl0_merge and wx_pbl_dts_merge in order to take into account
     2634! the update of the aridity coeficient beta.
     2635!
     2636        CALL wx_pbl_prelim_beta(knon, dtime, ywake_s, ybeta,  &
     2637                        BcoefQ_x, BcoefQ_w  &
     2638                        )
     2639        CALL wx_pbl0_merge(knon, ypplay, ypaprs,  &
     2640                          ywake_s, ydTs0, ydqs0, &
     2641                          yt_x, yt_w, yq_x, yq_w, &
     2642                          yu_x, yu_w, yv_x, yv_w, &
     2643                          ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, &
     2644                          ycdragm_x, ycdragm_w, &
     2645                          AcoefH_x, AcoefH_w, AcoefQ_x, AcoefQ_w, &
     2646                          AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, &
     2647                          BcoefH_x, BcoefH_w, BcoefQ_x, BcoefQ_w, &
     2648                          BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, &
     2649                          AcoefH_0, AcoefQ_0, AcoefU, AcoefV, &
     2650                          BcoefH_0, BcoefQ_0, BcoefU, BcoefV, &
     2651                          ycdragh, ycdragq, ycdragm, &
     2652                          yt1, yq1, yu1, yv1 &
     2653                          )
     2654          IF (iflag_split .eq. 2) THEN
     2655            CALL wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, &
     2656                            ywake_s, ybeta, ywake_cstar, ywake_dens, &
     2657                            AcoefH_x, AcoefH_w, &
     2658                            BcoefH_x, BcoefH_w, &
     2659                            AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0,  &
     2660                            AcoefH, AcoefQ, BcoefH, BcoefQ,  &
     2661                            HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, &
     2662                            phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, &
     2663                            yg_T, yg_Q, &
     2664                            yGamma_dTs_phiT, yGamma_dQs_phiQ, &
     2665                            ydTs_ins, ydqs_ins &
     2666                            )
     2667          ELSE !
     2668            AcoefH(:) = AcoefH_0(:)
     2669            AcoefQ(:) = AcoefQ_0(:)
     2670            BcoefH(:) = BcoefH_0(:)
     2671            BcoefQ(:) = BcoefQ_0(:)
     2672            yg_T(:) = 0.
     2673            yg_Q(:) = 0.
     2674            yGamma_dTs_phiT(:) = 0.
     2675            yGamma_dQs_phiQ(:) = 0.
     2676            ydTs_ins(:) = 0.
     2677            ydqs_ins(:) = 0.
     2678          ENDIF   ! (iflag_split .eq. 2)
     2679!
     2680        ELSE    ! (nsrf .ne. is_oce)
     2681          ybeta(1:knon) = 1.
     2682          yevap_pot(1:knon) = yevap(1:knon)
     2683          AcoefH(:) = AcoefH_0(:)
     2684          AcoefQ(:) = AcoefQ_0(:)
     2685          BcoefH(:) = BcoefH_0(:)
     2686          BcoefQ(:) = BcoefQ_0(:)
     2687          yg_T(:) = 0.
     2688          yg_Q(:) = 0.
     2689          yGamma_dTs_phiT(:) = 0.
     2690          yGamma_dQs_phiQ(:) = 0.
     2691          ydTs_ins(:) = 0.
     2692          ydqs_ins(:) = 0.
     2693        ENDIF   ! (nsrf .ne. is_oce)
     2694!
     2695        CALL wx_pbl_split(knon, nsrf, dtime, ywake_s, ybeta, iflag_split, &
     2696                       yg_T, yg_Q, &
     2697                       yGamma_dTs_phiT, yGamma_dQs_phiQ, &
     2698                       ydTs_ins, ydqs_ins, &
    24612699                       y_flux_t1, y_flux_q1, y_flux_u1, y_flux_v1, &
     2700!!!!                       HTRn_b, dd_HTRn, HTphiT_b, dd_HTphiT, &
     2701                       phiQ0_b, phiT0_b, &
    24622702                       y_flux_t1_x, y_flux_t1_w, &
    24632703                       y_flux_q1_x, y_flux_q1_w, &
     
    24652705                       y_flux_v1_x, y_flux_v1_w, &
    24662706                       yfluxlat_x, yfluxlat_w, &
    2467                        y_delta_tsurf &
     2707                       y_delta_qsats, &
     2708                       y_delta_tsurf_new, y_delta_qsurf &
    24682709                       )
     2710!
     2711         CALL wx_pbl_check(knon, dtime, ypplay, ypaprs, ywake_s, ybeta, iflag_split, &
     2712                       yTs, y_delta_tsurf,  &
     2713                       yqsurf, yTsurf_new,  &
     2714                       y_delta_tsurf_new, y_delta_qsats,  &
     2715                       AcoefH_x, AcoefH_w, &
     2716                       BcoefH_x, BcoefH_w, &
     2717                       AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0,  &
     2718                       AcoefH, AcoefQ, BcoefH, BcoefQ,  &
     2719                       y_flux_t1, y_flux_q1,  &
     2720                       y_flux_t1_x, y_flux_t1_w, &
     2721                       y_flux_q1_x, y_flux_q1_w)
     2722!
     2723         IF (nsrf .ne. is_oce) THEN
     2724           CALL wx_pbl_dts_check(knon, dtime, ypplay, ypaprs, ywake_s, ybeta, iflag_split, &
     2725                         yTs, y_delta_tsurf,  &
     2726                         yqsurf, yTsurf_new,  &
     2727                         y_delta_qsats, y_delta_tsurf_new, y_delta_qsurf,  &
     2728                         AcoefH_x, AcoefH_w, &
     2729                         BcoefH_x, BcoefH_w, &
     2730                         AcoefH_0, AcoefQ_0, BcoefH_0, BcoefQ_0,  &
     2731                         AcoefH, AcoefQ, BcoefH, BcoefQ,  &
     2732                         HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, &
     2733                         phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, &
     2734                         yg_T, yg_Q, &
     2735                         yGamma_dTs_phiT, yGamma_dQs_phiQ, &
     2736                         ydTs_ins, ydqs_ins, &
     2737                         y_flux_t1, y_flux_q1,  &
     2738                         y_flux_t1_x, y_flux_t1_w, &
     2739                         y_flux_q1_x, y_flux_q1_w )
     2740         ENDIF   ! (nsrf .ne. is_oce)
     2741!
     2742       ELSE  ! (iflag_split .ge. 1)
     2743         ybeta(1:knon) = 1.
     2744         yevap_pot(1:knon) = yevap(1:knon)
    24692745       ENDIF  ! (iflag_split .ge. 1)
     2746!
     2747       IF (prt_level >= 10) THEN
     2748         print *,'pbl_surface, ybeta , yevap, yevap_pot ', &
     2749                               ybeta , yevap, yevap_pot
     2750       ENDIF  ! (prt_level >= 10)
     2751!
    24702752!>jyg
    24712753!
     
    26482930       ENDIF  ! (iflag_split .eq.0)
    26492931!!!
    2650 
    2651         DO j = 1, knon
    2652           y_dflux_t(j) = y_dflux_t(j) * ypct(j)
    2653           y_dflux_q(j) = y_dflux_q(j) * ypct(j)
    2654         ENDDO
    2655 
     2932!!
     2933!!        DO j = 1, knon
     2934!!          y_dflux_t(j) = y_dflux_t(j) * ypct(j)
     2935!!          y_dflux_q(j) = y_dflux_q(j) * ypct(j)
     2936!!        ENDDO
     2937!!
    26562938!****************************************************************************************
    26572939! 13) Transform variables for output format :
     
    28213103          i = ni(j)
    28223104          evap(i,nsrf) = - flux_q(i,1,nsrf)                  !jyg
     3105          beta(i,nsrf) = ybeta(j)                             !jyg
    28233106          d_ts(i,nsrf) = y_d_ts(j)
    28243107!albedo SB >>>
     
    28363119          cdragh(i) = cdragh(i) + ycdragh(j)*ypct(j)
    28373120          cdragm(i) = cdragm(i) + ycdragm(j)*ypct(j)
    2838           dflux_t(i) = dflux_t(i) + y_dflux_t(j)
    2839           dflux_q(i) = dflux_q(i) + y_dflux_q(j)
     3121          dflux_t(i) = dflux_t(i) + y_dflux_t(j)*ypct(j)
     3122          dflux_q(i) = dflux_q(i) + y_dflux_q(j)*ypct(j)
    28403123#ifdef ISO
    28413124        do ixt=1,niso
     
    28443127        do ixt=1,ntraciso
    28453128          xtevap(ixt,i,nsrf) = - flux_xt(ixt,i,1,nsrf)
    2846           dflux_xt(ixt,i) = dflux_xt(ixt,i) + y_dflux_xt(ixt,j) 
     3129          dflux_xt(ixt,i) = dflux_xt(ixt,i) + y_dflux_xt(ixt,j)*ypct(j)
    28473130        enddo 
    28483131        IF (nsrf == is_lic) THEN
     
    28743157!!! nrlmd le 13/06/2011
    28753158!!jyg20170131          delta_tsurf(i,nsrf)=y_delta_tsurf(j)*ypct(j)
    2876           delta_tsurf(i,nsrf)=y_delta_tsurf(j)
     3159!!jyg20210118          delta_tsurf(i,nsrf)=y_delta_tsurf(j)
     3160          delta_tsurf(i,nsrf)=y_delta_tsurf_new(j)
     3161!
     3162          delta_qsurf(i) = delta_qsurf(i) + y_delta_qsurf(j)*ypct(j)
    28773163!
    28783164          cdragh_x(i) = cdragh_x(i) + ycdragh_x(j)*ypct(j)
     
    29183204              tke_x(i,k,nsrf)    = ytke(j,k)
    29193205              tke_x(i,k,is_ave) = tke_x(i,k,is_ave) + ytke(j,k)*ypct(j)
     3206
    29203207!>jyg
    29213208           ENDDO
     
    29313218!!            tke(i,k,is_ave) = tke(i,k,is_ave) + tke(i,k,nsrf)*ypct(j)
    29323219            tke_x(i,k,nsrf)   = ytke_x(j,k)
    2933             tke_x(i,k,is_ave)   = tke_x(i,k,is_ave) + tke_x(i,k,nsrf)*ypct(j)
     3220            tke_x(i,k,is_ave)   = tke_x(i,k,is_ave) + tke_x(i,k,nsrf)*ypct(j)       
    29343221            wake_dltke(i,k,is_ave)   = wake_dltke(i,k,is_ave) + wake_dltke(i,k,nsrf)*ypct(j)
     3222           
    29353223
    29363224!>jyg
     
    30693357          d_t_w(:,1), d_t_x(:,1), d_t(:,1)
    30703358       ENDIF
     3359
     3360       if (nsrf == is_oce .and. activate_ocean_skin >= 1) then
     3361          delta_sal = missing_val
     3362          ds_ns = missing_val
     3363          dt_ns = missing_val
     3364          delta_sst = missing_val
     3365          dter = missing_val
     3366          dser = missing_val
     3367          tkt = missing_val
     3368          tks = missing_val
     3369          taur = missing_val
     3370          sss = missing_val
     3371         
     3372          delta_sal(ni(:knon)) = ydelta_sal(:knon)
     3373          ds_ns(ni(:knon)) = yds_ns(:knon)
     3374          dt_ns(ni(:knon)) = ydt_ns(:knon)
     3375          delta_sst(ni(:knon)) = ydelta_sst(:knon)
     3376          dter(ni(:knon)) = ydter(:knon)
     3377          dser(ni(:knon)) = ydser(:knon)
     3378          tkt(ni(:knon)) = ytkt(:knon)
     3379          tks(ni(:knon)) = ytks(:knon)
     3380          taur(ni(:knon)) = ytaur(:knon)
     3381          sss(ni(:knon)) = ysss(:knon)
     3382       end if
     3383
     3384
    30713385
    30723386!****************************************************************************************
     
    31063420               * (ypaprs(j,1)-ypplay(j,1))
    31073421          tairsol(j) = yts(j) + y_d_ts(j)
    3108           tairsol_x(j) = tairsol(j) - ywake_s(j)*y_delta_tsurf(j)
     3422!!          tairsol_x(j) = tairsol(j) - ywake_s(j)*y_delta_tsurf(j)
     3423          tairsol_x(j) = tairsol(j) - ywake_s(j)*y_delta_tsurf_new(j)
    31093424          qairsol(j) = yqsurf(j)
    31103425        ENDDO
     
    31453460!!! jyg le 07/02/2012
    31463461       IF (iflag_split .eq.0) THEN
     3462        IF (iflag_new_t2mq2m==1) THEN
     3463         CALL stdlevvarn(klon, knon, nsrf, zxli, &
     3464            uzon, vmer, tair1, qair1, zgeo1, &
     3465            tairsol, qairsol, yz0m, yz0h_oupas, psfce, patm, &
     3466            yt2m, yq2m, yt10m, yq10m, yu10m, yustar, &
     3467            yn2mout(:, nsrf, :))
     3468        ELSE
    31473469        CALL stdlevvar(klon, knon, nsrf, zxli, &
    31483470            uzon, vmer, tair1, qair1, zgeo1, &
    31493471            tairsol, qairsol, yz0m, yz0h_oupas, psfce, patm, &
    31503472            yt2m, yq2m, yt10m, yq10m, yu10m, yustar)
     3473        ENDIF
    31513474       ELSE  !(iflag_split .eq.0)
     3475        IF (iflag_new_t2mq2m==1) THEN
     3476         CALL stdlevvarn(klon, knon, nsrf, zxli, &
     3477            uzon_x, vmer_x, tair1_x, qair1_x, zgeo1_x, &
     3478            tairsol_x, qairsol, yz0m, yz0h_oupas, psfce, patm, &
     3479            yt2m_x, yq2m_x, yt10m_x, yq10m_x, yu10m_x, yustar_x, &
     3480            yn2mout_x(:, nsrf, :))
     3481         CALL stdlevvarn(klon, knon, nsrf, zxli, &
     3482            uzon_w, vmer_w, tair1_w, qair1_w, zgeo1_w, &
     3483            tairsol_w, qairsol, yz0m, yz0h_oupas, psfce, patm, &
     3484            yt2m_w, yq2m_w, yt10m_w, yq10m_w, yu10m_w, yustar_w, &
     3485            yn2mout_w(:, nsrf, :))
     3486        ELSE
    31523487        CALL stdlevvar(klon, knon, nsrf, zxli, &
    31533488            uzon_x, vmer_x, tair1_x, qair1_x, zgeo1_x, &
     
    31583493            tairsol_w, qairsol, yz0m, yz0h_oupas, psfce, patm, &
    31593494            yt2m_w, yq2m_w, yt10m_w, yq10m_w, yu10m_w, yustar_w)
     3495        ENDIF
    31603496!!!
    31613497       ENDIF  ! (iflag_split .eq.0)
     
    31713507          u10m(i,nsrf)=(yu10m(j) * uzon(j))/SQRT(uzon(j)**2+vmer(j)**2)
    31723508          v10m(i,nsrf)=(yu10m(j) * vmer(j))/SQRT(uzon(j)**2+vmer(j)**2)
     3509!
     3510          DO k = 1, 6
     3511           n2mout(i,nsrf,k) = yn2mout(j,nsrf,k)
     3512          END DO 
     3513!
    31733514        ENDDO
    31743515       ELSE  !(iflag_split .eq.0)
     
    31813522          u10m_x(i,nsrf)=(yu10m_x(j) * uzon_x(j))/SQRT(uzon_x(j)**2+vmer_x(j)**2)
    31823523          v10m_x(i,nsrf)=(yu10m_x(j) * vmer_x(j))/SQRT(uzon_x(j)**2+vmer_x(j)**2)
     3524!
     3525          DO k = 1, 6
     3526           n2mout_x(i,nsrf,k) = yn2mout_x(j,nsrf,k)
     3527          END DO 
     3528!
    31833529        ENDDO
    31843530        DO j=1, knon
     
    31943540          u10m(i,nsrf) = u10m_x(i,nsrf) + wake_s(i)*(u10m_w(i,nsrf)-u10m_x(i,nsrf))
    31953541          v10m(i,nsrf) = v10m_x(i,nsrf) + wake_s(i)*(v10m_w(i,nsrf)-v10m_x(i,nsrf))
     3542!
     3543          DO k = 1, 6
     3544           n2mout_w(i,nsrf,k) = yn2mout_w(j,nsrf,k)
     3545          END DO 
     3546!
    31963547        ENDDO
    31973548!!!
     
    34793830#endif
    34803831!!!
    3481    
     3832
    34823833!
    34833834! Incrementer la temperature du sol
    34843835!
    34853836    zxtsol(:) = 0.0  ; zxfluxlat(:) = 0.0
    3486     zt2m(:) = 0.0    ; zq2m(:) = 0.0
     3837    zt2m(:) = 0.0    ; zq2m(:) = 0.0 ; zn2mout(:,:) = 0
    34873838    zustar(:)=0.0 ; zu10m(:) = 0.0   ; zv10m(:) = 0.0
    34883839    s_pblh(:) = 0.0  ; s_plcl(:) = 0.0
     
    35373888          zt2m(i)  = zt2m(i)  + t2m(i,nsrf)  * pctsrf(i,nsrf)
    35383889          zq2m(i)  = zq2m(i)  + q2m(i,nsrf)  * pctsrf(i,nsrf)
     3890!
     3891          DO k = 1, 6
     3892           zn2mout(i,k)  = zn2mout(i,k)  + n2mout(i,nsrf,k)  * pctsrf(i,nsrf)
     3893          ENDDO 
     3894!
    35393895          zustar(i) = zustar(i) + ustar(i,nsrf) * pctsrf(i,nsrf)
    35403896          wstar(i,is_ave)=wstar(i,is_ave)+wstar(i,nsrf)*pctsrf(i,nsrf)
     
    35683924          zt2m(i)  = zt2m(i)  + (t2m_x(i,nsrf)+wake_s(i)*(t2m_w(i,nsrf)-t2m_x(i,nsrf))) * pctsrf(i,nsrf)
    35693925          zq2m(i)  = zq2m(i)  + q2m_x(i,nsrf)  * pctsrf(i,nsrf)
     3926!
     3927          DO k = 1, 6
     3928           zn2mout(i,k)  = zn2mout(i,k)  + n2mout_x(i,nsrf,k)  * pctsrf(i,nsrf)
     3929          ENDDO
     3930!
    35703931          zustar(i) = zustar(i) + ustar_x(i,nsrf) * pctsrf(i,nsrf)
    35713932          wstar(i,is_ave)=wstar(i,is_ave)+wstar_x(i,nsrf)*pctsrf(i,nsrf)
     
    36494010    DO nsrf = 1, nbsrf
    36504011       DO i = 1, klon
    3651           zxqsurf(i) = zxqsurf(i) + qsurf(i,nsrf) * pctsrf(i,nsrf)
     4012          zxqsurf(i) = zxqsurf(i) + MAX(qsurf(i,nsrf),0.0) * pctsrf(i,nsrf)
    36524013          zxsnow(i)  = zxsnow(i)  + snow(i,nsrf)  * pctsrf(i,nsrf)
    36534014#ifdef ISO
     
    37184079    IF (ALLOCATED(qsurf)) DEALLOCATE(qsurf)
    37194080    IF (ALLOCATED(ftsoil)) DEALLOCATE(ftsoil)
     4081    IF (ALLOCATED(ydTs0)) DEALLOCATE(ydTs0)
     4082    IF (ALLOCATED(ydqs0)) DEALLOCATE(ydqs0)
    37204083#ifdef ISO
    37214084    IF (ALLOCATED(xtsnow)) DEALLOCATE(xtsnow)
     
    37394102
    37404103!albedo SB >>>
    3741 SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, &
    3742      evap, z0m, z0h, agesno,                                  &
    3743      tsurf,alb_dir,alb_dif, ustar, u10m, v10m, tke &
     4104  SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, &
     4105       evap, z0m, z0h, agesno,                                  &
     4106       tsurf,alb_dir,alb_dif, ustar, u10m, v10m, tke &
    37444107#ifdef ISO
    37454108     ,xtevap  &
     
    37504113
    37514114    USE indice_sol_mod
     4115    use phys_state_var_mod, only: delta_sal, ds_ns, dt_ns, delta_sst
     4116    use config_ocean_skin_m, only: activate_ocean_skin
    37524117#ifdef ISO
    37534118  USE infotrac_phy, ONLY: ntraciso   
     
    38524217                      alb_dif(i,k,nsrf) = 0.06
    38534218                   ENDDO
     4219                   if (activate_ocean_skin >= 1) then
     4220                      if (activate_ocean_skin == 2 &
     4221                           .and. type_ocean == "couple") then
     4222                         delta_sal(i) = 0.
     4223                         delta_sst(i) = 0.
     4224                      end if
     4225                     
     4226                      ds_ns(i) = 0.
     4227                      dt_ns(i) = 0.
     4228                   end if
    38544229                ELSE IF (nsrf.EQ.is_sic) THEN
    38554230                   tsurf(i,nsrf) = 271.35
  • LMDZ6/trunk/libf/phylmdiso/phyetat0.F90

    r3927 r3940  
    1 ! $Id: phyetat0.F90 3581 2019-10-10 12:35:59Z oboucher $
     1! $Id: phyetat0.F90 3890 2021-05-05 15:15:06Z jyg $
    22
    33SUBROUTINE phyetat0 (fichnom, clesphy0, tabcntr0)
     
    1919       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, radpas, radsol, rain_fall, ratqs, &
    2020       rnebcon, rugoro, sig1, snow_fall, solaire_etat0, sollw, sollwdown, &
    21        solsw, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
    22        wake_deltat, wake_delta_pbl_TKE, delta_tsurf, wake_fip, wake_pe, &
     21       solsw, solswfdiff, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
     22       wake_deltat, wake_delta_pbl_TKE, delta_tsurf, beta_aridity, wake_fip, wake_pe, &
    2323       wake_s, wake_dens, zgam, zmax0, zmea, zpic, zsig, &
    2424#ifdef ISO
     
    2727#endif
    2828       zstd, zthe, zval, ale_bl, ale_bl_trig, alp_bl, u10m, v10m, treedrg, &
    29        ale_wake, ale_bl_stat
     29       ale_wake, ale_bl_stat, ds_ns, dt_ns, delta_sst, delta_sal, ratqs_inter
    3030!FC
    3131  USE geometry_mod, ONLY : longitude_deg, latitude_deg
     
    3838  USE ocean_slab_mod, ONLY: nslay, tslab, seaice, tice, ocean_slab_init
    3939  USE time_phylmdz_mod, ONLY: init_iteration, pdtphys, itau_phy
     40#ifdef CPP_XIOS
     41  USE wxios, ONLY: missing_val
     42#else
     43  use netcdf, only: missing_val => nf90_fill_real
     44#endif
     45  use config_ocean_skin_m, only: activate_ocean_skin
    4046#ifdef ISO
    4147  USE infotrac_phy, ONLY: ntraciso,niso,iso_num
     
    5258  ! Objet: Lecture de l'etat initial pour la physique
    5359  !======================================================================
    54   include "netcdf.inc"
    5560  include "dimsoil.h"
    5661  include "clesphys.h"
     
    330335
    331336  found=phyetat0_get(1,solsw,"solsw","net SW radiation surf",0.)
     337  found=phyetat0_get(1,solswfdiff,"solswfdiff","fraction of SW radiation surf that is diffuse",1.)
    332338  found=phyetat0_get(1,sollw,"sollw","net LW radiation surf",0.)
    333339  found=phyetat0_get(1,sollwdown,"sollwdown","down LW radiation surf",0.)
    334340  IF (.NOT. found) THEN
    335      sollwdown = 0. ;  zts=0.
    336      do nsrf=1,nbsrf
     341     sollwdown(:) = 0. ;  zts(:)=0.
     342     DO nsrf=1,nbsrf
    337343        zts(:)=zts(:)+ftsol(:,nsrf)*pctsrf(:,nsrf)
    338      enddo
     344     ENDDO
    339345     sollwdown(:)=sollw(:)+RSIGMA*zts(:)**4
    340346  ENDIF
     
    413419  IF (iflag_pbl>1 .AND. iflag_wake>=1  .AND. iflag_pbl_split >=1 ) then
    414420    found=phyetat0_srf(klev+1,wake_delta_pbl_tke,"DELTATKE","Del TKE wk/env",0.)
    415     found=phyetat0_srf(1,delta_tsurf,"DELTA_TSURF","Delta Ts wk/env ",0.)
     421!!    found=phyetat0_srf(1,delta_tsurf,"DELTA_TSURF","Delta Ts wk/env ",0.)
     422    found=phyetat0_srf(1,delta_tsurf,"DELTATS","Delta Ts wk/env ",0.)
     423!!    found=phyetat0_srf(1,beta_aridity,"BETA_S","Aridity factor ",1.)
     424    found=phyetat0_srf(1,beta_aridity,"BETAS","Aridity factor ",1.)
    416425  ENDIF   !(iflag_pbl>1 .AND. iflag_wake>=1 .AND. iflag_pbl_split >=1 )
    417426
     
    452461  found=phyetat0_get(1,ale_bl_stat,"ALE_BL_STAT","ALE_BL_STAT",0.)
    453462
     463! fisrtilp/Clouds 0.002 could be ratqsbas. But can stay like this as well
     464  found=phyetat0_get(klev,ratqs_inter,"RATQS_INTER","Relative width of the lsc sugrid scale water",0.002)
     465
     466
    454467!===========================================
    455468  ! Read and send field trs to traclmdz
     
    468481  ENDIF
    469482
    470 !--OB now this is for co2i
    471   IF (type_trac == 'co2i') THEN
     483!--OB now this is for co2i - ThL: and therefore also for inco
     484  IF (type_trac == 'co2i' .OR. type_trac == 'inco') THEN
    472485     IF (carbon_cycle_cpl) THEN
    473486        ALLOCATE(co2_send(klon), stat=ierr)
     
    564577  ENDIF ! Slab       
    565578
     579  if (activate_ocean_skin >= 1) then
     580     if (activate_ocean_skin == 2 .and. type_ocean == 'couple') then
     581        found = phyetat0_get(1, delta_sal, "delta_sal", &
     582             "ocean-air interface salinity minus bulk salinity", 0.)
     583        found = phyetat0_get(1, delta_sst, "delta_SST", &
     584             "ocean-air interface temperature minus bulk SST", 0.)
     585     end if
     586     
     587     found = phyetat0_get(1, ds_ns, "dS_ns", "delta salinity near surface", 0.)
     588     found = phyetat0_get(1, dt_ns, "dT_ns", "delta temperature near surface", &
     589          0.)
     590
     591     where (pctsrf(:, is_oce) == 0.)
     592        ds_ns = missing_val
     593        dt_ns = missing_val
     594        delta_sst = missing_val
     595        delta_sal = missing_val
     596     end where
     597  end if
     598 
     599
    566600  ! on ferme le fichier
    567601  CALL close_startphy
  • LMDZ6/trunk/libf/phylmdiso/phyredem.F90

    r3927 r3940  
    1212  USE fonte_neige_mod,  ONLY : fonte_neige_final
    1313  USE pbl_surface_mod,  ONLY : pbl_surface_final
    14   USE phys_state_var_mod, ONLY: radpas, zmasq, pctsrf, ftsol, falb_dir,      &
     14  USE phys_state_var_mod, ONLY: radpas, zmasq, pctsrf,                       &
     15                                ftsol, beta_aridity, delta_tsurf, falb_dir,  &
    1516                                falb_dif, qsol, fevap, radsol, solsw, sollw, &
    1617                                sollwdown, rain_fall, snow_fall, z0m, z0h,   &
     
    2627                                detr_therm, ale_bl, ale_bl_trig, alp_bl,     &
    2728                                ale_wake, ale_bl_stat,                       &
    28                                 du_gwd_rando, du_gwd_front, u10m, v10m,      &
    29                                 treedrg
     29                                du_gwd_rando, du_gwd_front, u10m, v10m, &
     30                                treedrg, solswfdiff, delta_sal, ds_ns, dt_ns, &
     31                                delta_sst, ratqs_inter
    3032#ifdef ISO
    3133  USE phys_state_var_mod, ONLY: xtsol, fxtevap,xtrain_fall, xtsnow_fall,     &
     
    4850  USE ocean_slab_mod, ONLY : nslay, tslab, seaice, tice, fsic
    4951  USE time_phylmdz_mod, ONLY: annee_ref, day_end, itau_phy, pdtphys
     52  use config_ocean_skin_m, only: activate_ocean_skin 
    5053
    5154  IMPLICIT none
     
    179182    END IF
    180183
     184!    Surface variables
    181185    CALL put_field_srf1(pass,"TS","Temperature",ftsol(:,:))
     186
     187!!    CALL put_field_srf1(pass,"DELTA_TS","w-x surface temperature difference", delta_tsurf(:,:))
     188    CALL put_field_srf1(pass,"DELTATS","w-x surface temperature difference", delta_tsurf(:,:))
     189
     190!    CALL put_field_srf1(pass,"BETA_S","Aridity factor", beta_aridity(:,:))
     191    CALL put_field_srf1(pass,"BETAS","Aridity factor", beta_aridity(:,:))
     192!    End surface variables
    182193
    183194! ================== Albedo =======================================
     
    209220
    210221    CALL put_field(pass,"solsw", "Rayonnement solaire a la surface", solsw)
     222
     223    CALL put_field(pass,"solswfdiff", "Fraction du rayonnement solaire a la surface qui est diffus", solswfdiff)
    211224
    212225    CALL put_field(pass,"sollw", "Rayonnement IF a la surface", sollw)
     
    323336
    324337    CALL put_field(pass,"ALE_BL_STAT", "ALE_BL_STAT", ale_bl_stat)
     338
     339
     340    ! fisrtilp/clouds
     341    CALL put_field(pass,"RATQS_INTER","Relative width of the lsc sugrid scale water",ratqs_inter)
    325342
    326343
     
    366383    IF (.not. ok_hines .and. ok_gwd_rando) call put_field(pass,"du_gwd_front", &
    367384         "tendency on zonal wind due to acama gravity waves", du_gwd_front)
     385
     386    if (activate_ocean_skin >= 1) then
     387       if (activate_ocean_skin == 2 .and. type_ocean == 'couple') then
     388          CALL put_field(pass, "delta_sal", &
     389               "ocean-air interface salinity minus bulk salinity", delta_sal)
     390          CALL put_field(pass, "delta_SST", &
     391               "ocean-air interface temperature minus bulk SST", delta_sst)
     392       end if
     393       
     394       CALL put_field(pass, "dS_ns", "delta salinity near surface", ds_ns)
     395       CALL put_field(pass, "dT_ns", "delta temperature near surface", dT_ns)
     396    end if
    368397
    369398#ifdef ISO
  • LMDZ6/trunk/libf/phylmdiso/phys_local_var_mod.F90

    r3927 r3940  
    11!
    2 ! $Id: phys_local_var_mod.F90 3662 2020-04-12 16:41:53Z oboucher $
     2! $Id: phys_local_var_mod.F90 3888 2021-05-05 10:50:37Z jyg $
    33!
    44      MODULE phys_local_var_mod
     
    1616      REAL, SAVE, ALLOCATABLE :: u_seri(:,:), v_seri(:,:)
    1717      !$OMP THREADPRIVATE(u_seri, v_seri)
    18       REAL, SAVE, ALLOCATABLE :: l_mixmin(:,:,:), l_mix(:,:,:)
    19       !$OMP THREADPRIVATE(l_mixmin, l_mix)
    20 
     18      REAL, SAVE, ALLOCATABLE :: l_mixmin(:,:,:), l_mix(:,:,:), tke_dissip(:,:,:)
     19      !$OMP THREADPRIVATE(l_mixmin, l_mix, tke_dissip)
    2120      REAL, SAVE, ALLOCATABLE :: tr_seri(:,:,:)
    2221      !$OMP THREADPRIVATE(tr_seri)
     
    334333!!!OMP THREADPRIVATE(d_s_the, d_dens_the)
    335334      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)           :: d_deltat_ajs_cv, d_deltaq_ajs_cv
    336 !$OMP THREADPRIVATE(d_deltat_ajs_cv, d_deltaq_ajs_cv)
     335!$OMP THREADPRIVATE(d_deltat_ajs_cv, d_deltaq_ajs_cv)                       
    337336#ifdef ISO
    338337    REAL, SAVE, ALLOCATABLE,DIMENSION(:,:,:)          :: d_deltaxt_wk
     
    355354      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: cldh, cldl, cldm, cldq, cldt, qsat2m
    356355!$OMP THREADPRIVATE(cldh, cldl, cldm, cldq, cldt, qsat2m )
    357       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: cldhjn, cldljn, cldmjn,cldtjn
    358 !$OMP THREADPRIVATE(cldhjn, cldljn, cldmjn, cldtjn)
     356!AS: cldhjn, cldljn, cldmjn,cldtjn pas utilisés en tant que variables, juste noms de diagnostics
    359357      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: JrNt
    360358!$OMP THREADPRIVATE(JrNt)
     
    411409      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxfluxlat_x, zxfluxlat_w
    412410!$OMP THREADPRIVATE(zxfluxlat_x, zxfluxlat_w)
     411      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: delta_qsurf
     412!$OMP THREADPRIVATE(delta_qsurf)
    413413!jyg<
    414414!!! Entrees supplementaires couche-limite
     
    455455      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: t2m_min_mon, t2m_max_mon
    456456!$OMP THREADPRIVATE(t2m_min_mon, t2m_max_mon)
    457       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zq2m_cor, zt2m_cor
    458 !$OMP THREADPRIVATE(zq2m_cor, zt2m_cor)
    459       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zu10m_cor, zv10m_cor
    460 !$OMP THREADPRIVATE(zu10m_cor, zv10m_cor)
    461       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zrh2m_cor, zqsat2m_cor
    462 !$OMP THREADPRIVATE(zrh2m_cor, zqsat2m_cor)
    463457      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: weak_inversion
    464458!$OMP THREADPRIVATE(weak_inversion)
     
    567561      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: ref_liq_pi, ref_ice_pi
    568562!$OMP THREADPRIVATE(ref_liq_pi, ref_ice_pi)
    569       REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zx_rh
    570 !$OMP THREADPRIVATE(zx_rh)
     563      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zx_rh, zx_rhl, zx_rhi
     564!$OMP THREADPRIVATE(zx_rh, zx_rhl, zx_rhi)
    571565      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: prfl, psfl, fraca
    572566!$OMP THREADPRIVATE(prfl, psfl, fraca)
     
    604598      REAL, ALLOCATABLE, SAVE, DIMENSION(:) :: p_tropopause, z_tropopause, t_tropopause
    605599!$OMP THREADPRIVATE(p_tropopause, z_tropopause, t_tropopause)
     600
     601      INTEGER,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zn2mout
     602!$OMP THREADPRIVATE(zn2mout)
    606603
    607604#ifdef CPP_StratAer
     
    694691      ALLOCATE(t_seri(klon,klev),q_seri(klon,klev),ql_seri(klon,klev),qs_seri(klon,klev))
    695692      ALLOCATE(u_seri(klon,klev),v_seri(klon,klev))
    696       ALLOCATE(l_mixmin(klon,klev+1,nbsrf), l_mix(klon,klev+1,nbsrf))
    697       l_mix(:,:,:)=0. ; l_mixmin(:,:,:)=0. ! doit etre initialse car pas toujours remplis
     693      ALLOCATE(l_mixmin(klon,klev+1,nbsrf), l_mix(klon,klev+1,nbsrf), tke_dissip(klon,klev+1,nbsrf))
     694      l_mix(:,:,:)=0. ; l_mixmin(:,:,:)=0. ; tke_dissip(:,:,:)=0. ! doit etre initialse car pas toujours remplis
    698695
    699696      ALLOCATE(tr_seri(klon,klev,nbtr))
     
    882879      ALLOCATE(cdragm(klon), cdragh(klon), cldh(klon), cldl(klon))
    883880      ALLOCATE(cldm(klon), cldq(klon), cldt(klon), qsat2m(klon))
    884       ALLOCATE(cldhjn(klon), cldljn(klon), cldmjn(klon), cldtjn(klon))
    885881      ALLOCATE(JrNt(klon))
    886882      ALLOCATE(dthmin(klon), evap(klon), fder(klon), plcl(klon), plfc(klon))
     
    915911      ALLOCATE(sens_x(klon), sens_w(klon))
    916912      ALLOCATE(zxfluxlat_x(klon), zxfluxlat_w(klon))
     913      ALLOCATE(delta_qsurf(klon))
    917914!jyg<
    918915!!      ALLOCATE(t_x(klon,klev), t_w(klon,klev))
     
    940937      ALLOCATE(zt2m_min_mon(klon), zt2m_max_mon(klon))
    941938      ALLOCATE(t2m_min_mon(klon), t2m_max_mon(klon))
    942       ALLOCATE(zq2m_cor(klon), zt2m_cor(klon), zu10m_cor(klon), zv10m_cor(klon))
    943       ALLOCATE(zrh2m_cor(klon), zqsat2m_cor(klon))
    944939      ALLOCATE(sens(klon), flwp(klon), fiwp(klon))
    945940      ALLOCATE(alp_bl_conv(klon), alp_bl_det(klon))
     
    962957      ALLOCATE(ref_liq(klon, klev), ref_ice(klon, klev), theta(klon, klev))
    963958      ALLOCATE(ref_liq_pi(klon, klev), ref_ice_pi(klon, klev))
    964       ALLOCATE(zphi(klon, klev), zx_rh(klon, klev))
     959      ALLOCATE(zphi(klon, klev), zx_rh(klon, klev), zx_rhl(klon,klev), zx_rhi(klon,klev))
    965960      ALLOCATE(pmfd(klon, klev), pmfu(klon, klev))
    966961
     
    10391034      ALLOCATE (z_tropopause(klon))
    10401035      ALLOCATE (t_tropopause(klon))
     1036
     1037      ALLOCATE(zn2mout(klon,6))
    10411038
    10421039#ifdef CPP_StratAer
     
    10891086      DEALLOCATE(t_seri,q_seri,ql_seri,qs_seri)
    10901087      DEALLOCATE(u_seri,v_seri)
    1091       DEALLOCATE(l_mixmin,l_mix)
     1088      DEALLOCATE(l_mixmin,l_mix, tke_dissip)
    10921089
    10931090      DEALLOCATE(tr_seri)
     
    12531250      DEALLOCATE(cdragm, cdragh, cldh, cldl)
    12541251      DEALLOCATE(cldm, cldq, cldt, qsat2m)
    1255       DEALLOCATE(cldljn, cldmjn, cldhjn, cldtjn, JrNt)
     1252      DEALLOCATE(JrNt)
    12561253      DEALLOCATE(dthmin, evap, fder, plcl, plfc)
    12571254      DEALLOCATE(prw, prlw, prsw, zustar, zu10m, zv10m, rh2m, s_lcl)
     
    12741271      DEALLOCATE(sens_x, sens_w)
    12751272      DEALLOCATE(zxfluxlat_x, zxfluxlat_w)
     1273      DEALLOCATE(delta_qsurf)
    12761274!jyg<
    12771275!!      DEALLOCATE(t_x, t_w)
     
    13031301      DEALLOCATE(zt2m_min_mon, zt2m_max_mon)
    13041302      DEALLOCATE(t2m_min_mon, t2m_max_mon)
    1305       DEALLOCATE(zq2m_cor, zt2m_cor, zu10m_cor, zv10m_cor)
    1306       DEALLOCATE(zrh2m_cor, zqsat2m_cor)
    13071303      DEALLOCATE(sens, flwp, fiwp)
    13081304      DEALLOCATE(alp_bl_conv,alp_bl_det)
     
    13221318      DEALLOCATE(ref_liq, ref_ice, theta)
    13231319      DEALLOCATE(ref_liq_pi, ref_ice_pi)
    1324       DEALLOCATE(zphi, zx_rh)
     1320      DEALLOCATE(zphi, zx_rh, zx_rhl, zx_rhi)
    13251321      DEALLOCATE(pmfd, pmfu)
    13261322
     
    13911387      DEALLOCATE (z_tropopause)
    13921388      DEALLOCATE (t_tropopause)
     1389      DEALLOCATE(zn2mout)
    13931390
    13941391#ifdef CPP_StratAer
  • LMDZ6/trunk/libf/phylmdiso/phys_output_ctrlout_mod.F90

    r3927 r3940  
    11!
    2 ! $Id: phys_output_ctrlout_mod.F90 3691 2020-05-30 15:37:19Z oboucher $
     2! $Id: phys_output_ctrlout_mod.F90 3888 2021-05-05 10:50:37Z jyg $
    33!
    44MODULE phys_output_ctrlout_mod
     
    77  USE indice_sol_mod
    88  USE aero_mod
    9 
    109
    1110  IMPLICIT NONE
     
    273272    't2m_sic', "Temp 2m "//clnsurf(4), "K", (/ ('', i=1, 10) /)) /)
    274273
     274  TYPE(ctrl_out), SAVE :: o_nt2mout = ctrl_out((/ 1, 1, 1, 5, 10, 10, 11, 11, 11, 11/), &
     275    'nt2mout', 'Nbt2m out of range complete computation', '-', (/ ('', i=1, 10) /))
     276  TYPE(ctrl_out), SAVE :: o_nq2mout = ctrl_out((/ 1, 1, 1, 5, 10, 10, 11, 11, 11, 11/), &
     277    'nq2mout', 'Nbq2m out of range complete computation', '-', (/ ('', i=1, 10) /))
     278  TYPE(ctrl_out), SAVE :: o_nu2mout = ctrl_out((/ 1, 1, 1, 5, 10, 10, 11, 11, 11, 11/), &
     279    'nu2mout', 'Nbu2m out of range complete computation', '-', (/ ('', i=1, 10) /))
     280
     281  TYPE(ctrl_out), SAVE :: o_nt2moutfg = ctrl_out((/ 1, 1, 1, 5, 10, 10, 11, 11, 11, 11/), &
     282    'nt2moutfg', 'Nbt2m out of range complete/fgRi1 computation', '-', (/ ('', i=1, 10) /))
     283  TYPE(ctrl_out), SAVE :: o_nq2moutfg = ctrl_out((/ 1, 1, 1, 5, 10, 10, 11, 11, 11, 11/), &
     284    'nq2moutfg', 'Nbq2m out of range complete/fgRi1 computation', '-', (/ ('', i=1, 10) /))
     285  TYPE(ctrl_out), SAVE :: o_nu2moutfg = ctrl_out((/ 1, 1, 1, 5, 10, 10, 11, 11, 11, 11/), &
     286    'nu2moutfg', 'Nbu2m out of range complete/fgRi1 computation', '-', (/ ('', i=1, 10) /))
     287
    275288  TYPE(ctrl_out), SAVE :: o_gusts = ctrl_out((/ 1, 1, 1, 10, 10, 10, 11, 11, 11, 11/), &
    276289    'gusts', 'surface gustiness', 'm2/s2', (/ ('', i=1, 10) /))
     
    281294    'wind100m', '100-m wind speed', 'm/s', (/ ('', i=1, 10) /))
    282295  TYPE(ctrl_out), SAVE :: o_loadfactor_wind_onshore = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    283     'load_factor_wind_onshore', 'Load factor for onshore windmill', '-', (/ ('', i=1, 10) /))
     296    'woncfr', 'Onshore Wind Capacity factor', 'kW/kW_installed', (/ ('', i=1, 10) /))
    284297  TYPE(ctrl_out), SAVE :: o_loadfactor_wind_offshore = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    285     'load_factor_wind_offshore', 'Load factor for offshore windmill', '-', (/ ('', i=1, 10) /))
     298    'wofcfr', 'Offshore Wind Capacity factor', 'kW/kW_installed', (/ ('', i=1, 10) /))
    286299  TYPE(ctrl_out), SAVE :: o_wind10max = ctrl_out((/ 10, 1, 10, 10, 10, 10, 11, 11, 11, 11/), &
    287300    'wind10max', '10m wind speed max', 'm/s', &
     
    462475  TYPE(ctrl_out), SAVE :: o_SWupSFCcleanclr = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    463476    'SWupSFCcleanclr', 'SWup clear sky clean (no aerosol) at surface', 'W/m2', (/ ('', i=1, 10) /))
    464   TYPE(ctrl_out), SAVE :: o_SWdnSFC = ctrl_out((/ 1, 1, 10, 10, 5, 10, 11, 11, 11, 11/), &
     477  TYPE(ctrl_out), SAVE :: o_fdiffSWdnSFC = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     478    'fdiffSWdnSFC', 'Fraction of diffuse SWdn at surface', 'W/m2', (/ ('', i=1, 10) /))
     479  TYPE(ctrl_out), SAVE :: o_SWdnSFC = ctrl_out((/ 1, 1, 1, 10, 5, 10, 11, 11, 11, 11/), &
    465480    'SWdnSFC', 'SWdn at surface', 'W/m2', (/ ('', i=1, 10) /))
    466481  TYPE(ctrl_out), SAVE :: o_SWdnSFCclr = ctrl_out((/ 1, 4, 10, 10, 5, 10, 11, 11, 11, 11/), &
     
    567582  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_evappot_srf  = (/ &
    568583      ctrl_out((/ 1, 6, 10, 10, 10, 10, 11, 11, 11, 11/),'evappot_ter',    &
    569       "Temperature"//clnsurf(1),"K", (/ ('', i=1, 10) /)),   &
     584      "Potential evaporation "//clnsurf(1),"kg/(m2*s)", (/ ('', i=1, 10) /)),   &
    570585      ctrl_out((/ 4, 6, 10, 10, 10, 10, 11, 11, 11, 11/),'evappot_lic',    &
    571       "Temperature"//clnsurf(2),"K", (/ ('', i=1, 10) /)),   &
     586      "Potential evaporation "//clnsurf(2),"kg/(m2*s)", (/ ('', i=1, 10) /)),   &
    572587      ctrl_out((/ 4, 6, 10, 10, 10, 10, 11, 11, 11, 11/),'evappot_oce',    &
    573       "Temperature"//clnsurf(3),"K", (/ ('', i=1, 10) /)),   &
     588      "Potential evaporation "//clnsurf(3),"kg/(m2*s)", (/ ('', i=1, 10) /)),   &
    574589      ctrl_out((/ 4, 6, 10, 10, 10, 10, 11, 11, 11, 11/),'evappot_sic',    &
    575       "Temperature"//clnsurf(4),"K", (/ ('', i=1, 10) /)) /)
     590      "Potential evaporation "//clnsurf(4),"kg/(m2*s)", (/ ('', i=1, 10) /)) /)
    576591
    577592  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_sens_srf     = (/          &
     
    695710    'iwp', 'Cloud ice water path', 'kg/m2', (/ ('', i=1, 10) /))
    696711  TYPE(ctrl_out), SAVE :: o_ue = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    697     'ue', 'Zonal dry static energy transport', '-', (/ ('', i=1, 10) /))
     712    'ue', 'Zonal dry static energy transport', 'J/m/s', (/ ('', i=1, 10) /))
    698713  TYPE(ctrl_out), SAVE :: o_ve = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    699     've', 'Merid dry static energy transport', '-', (/ ('', i=1, 10) /))
     714    've', 'Merid dry static energy transport', 'J/m/s', (/ ('', i=1, 10) /))
    700715  TYPE(ctrl_out), SAVE :: o_uq = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    701     'uq', 'Zonal humidity transport', '-', (/ ('', i=1, 10) /))
     716    'uq', 'Zonal humidity transport', 'kg/m/s', (/ ('', i=1, 10) /))
    702717  TYPE(ctrl_out), SAVE :: o_vq = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    703     'vq', 'Merid humidity transport', '-', (/ ('', i=1, 10) /))
     718    'vq', 'Merid humidity transport', 'kg/m/s', (/ ('', i=1, 10) /))
    704719  TYPE(ctrl_out), SAVE :: o_uwat = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    705     'uwat', 'Zonal total water transport', '-', (/ ('', i=1, 10) /))
     720    'uwat', 'Zonal total water transport', 'kg/m/s', (/ ('', i=1, 10) /))
    706721  TYPE(ctrl_out), SAVE :: o_vwat = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    707     'vwat', 'Merid total water transport', '-', (/ ('', i=1, 10) /))
     722    'vwat', 'Merid total water transport', 'kg/m/s', (/ ('', i=1, 10) /))
    708723  TYPE(ctrl_out), SAVE :: o_cape = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    709724    'cape', 'Conv avlbl pot ener', 'J/kg', (/ ('', i=1, 10) /))
     
    799814'flat_w', 'flat within_wake', 'W/m2', (/ ('', i=1, 10) /))
    800815!!
    801   type(ctrl_out),save :: o_delta_tsurf    = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    802 'delta_tsurf', 'Temperature difference (w-x)', 'K', (/ ('', i=1, 10) /))
    803816  type(ctrl_out),save :: o_cdragh_x       = ctrl_out((/ 1, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    804817'cdragh_x', 'cdragh off-wake', '', (/ ('', i=1, 10) /))
     
    10041017  TYPE(ctrl_out), SAVE :: o_tke = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    10051018    'tke ', 'TKE', 'm2/s2', (/ ('', i=1, 10) /))
     1019  TYPE(ctrl_out), SAVE :: o_tke_dissip = ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1020    'tke_dissip ', 'TKE DISSIPATION', 'm2/s3', (/ ('', i=1, 10) /))   
    10061021  TYPE(ctrl_out), SAVE :: o_tke_max = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    10071022    'tke_max', 'TKE max', 'm2/s2',                                  &
    10081023      (/ 't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)', &
    10091024         't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)' /))
    1010 
    10111025  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_tke_srf      = (/             &
    10121026      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'tke_ter',       &
     
    10501064      "min PBL mixing length "//clnsurf(4),"m", (/ ('', i=1, 10) /)) /)
    10511065
     1066
    10521067  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_tke_max_srf  = (/                          &
    10531068      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'tke_max_ter',                &
     
    10771092      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'dltpbltke_sic',       &
    10781093      "TKE difference (w - x) "//clnsurf(4),"-", (/ ('', i=1, 10) /)) /)
     1094
     1095  TYPE(ctrl_out), SAVE :: o_delta_tsurf = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1096    'delta_tsurf ', 'T_surf difference (w - x)', 'K', (/ ('', i=1, 10) /))
     1097  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_delta_tsurf_srf      = (/             &
     1098      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'delta_tsurf_ter',       &
     1099      "T_surf difference (w - x) "//clnsurf(1),"-", (/ ('', i=1, 10) /)), &
     1100      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'delta_tsurf_lic',       &
     1101      "T_surf difference (w - x) "//clnsurf(2),"-", (/ ('', i=1, 10) /)), &
     1102      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'delta_tsurf_oce',       &
     1103      "T_surf difference (w - x) "//clnsurf(3),"-", (/ ('', i=1, 10) /)), &
     1104      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'delta_tsurf_sic',       &
     1105      "T_surf difference (w - x) "//clnsurf(4),"-", (/ ('', i=1, 10) /)) /)
    10791106
    10801107  TYPE(ctrl_out), SAVE :: o_kz = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     
    13041331  TYPE(ctrl_out), SAVE :: o_flx_co2_land = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
    13051332    'flx_co2_land', 'CO2 flux from the land', '1', (/ ('', i=1, 10) /))
     1333  TYPE(ctrl_out), SAVE :: o_flx_co2_ocean_cor = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
     1334    'flx_co2_ocean_cor', 'correction of the CO2 flux from the ocean', 'kg CO2 m-2 s-1', (/ ('', i=1, 10) /))
     1335  TYPE(ctrl_out), SAVE :: o_flx_co2_land_cor = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
     1336    'flx_co2_land_cor', 'correction of the CO2 flux from the land', 'kg CO2 m-2 s-1', (/ ('', i=1, 10) /))
    13061337
    13071338#ifdef CPP_StratAer
     
    14331464  TYPE(ctrl_out), SAVE :: o_rhum = ctrl_out((/ 2, 5, 10, 10, 10, 10, 11, 11, 11, 11/), &
    14341465    'rhum', 'Relative humidity', '-', (/ ('', i=1, 10) /))
     1466  TYPE(ctrl_out), SAVE :: o_rhl = ctrl_out((/ 2, 6, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1467    'rhl', 'Relative humidity wrt liquid', '%', (/ ('', i=1, 10) /))
     1468  TYPE(ctrl_out), SAVE :: o_rhi = ctrl_out((/ 2, 6, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1469    'rhi', 'Relative humidity wrt ice', '%', (/ ('', i=1, 10) /))
    14351470  TYPE(ctrl_out), SAVE :: o_ozone = ctrl_out((/ 2, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    14361471    'ozone', 'Ozone mole fraction', '-', (/ ('', i=1, 10) /))
     
    19702005#endif
    19712006
     2007   type(ctrl_out), save:: o_delta_sst &
     2008        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'delta_SST', &
     2009        "ocean-air interface temperature minus bulk SST", "K", '')
     2010
     2011   type(ctrl_out), save:: o_delta_sal &
     2012        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'delta_sal', &
     2013        "ocean-air interface salinity minus bulk salinity", "ppt", '')
     2014
     2015   type(ctrl_out), save:: o_ds_ns &
     2016        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'dS_ns', &
     2017        "subskin salinity minus foundation salinity", "ppt", '')
     2018
     2019   type(ctrl_out), save:: o_dt_ns &
     2020        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'dT_ns', &
     2021        "subskin temperature minus foundation temperature", "K", '')
     2022
     2023   type(ctrl_out), save:: o_dter &
     2024        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'dTer', &
     2025        "ocean-air interface temperature minus sub-skin temperature", "K", '')
     2026
     2027   type(ctrl_out), save:: o_dser &
     2028        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'dSer', &
     2029        "ocean-air interface salinity minus sub-skin salinity", "ppt", '')
     2030
     2031   type(ctrl_out), save:: o_tkt &
     2032        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'tkt', &
     2033        "thickness of thermal microlayer", "m", '')
     2034
     2035   type(ctrl_out), save:: o_tks &
     2036        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'tks', &
     2037        "thickness of salinity microlayer", "m", '')
     2038
     2039   type(ctrl_out), save:: o_taur &
     2040        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'taur', &
     2041        "momentum flux due to rain", "Pa", '')
     2042
     2043   type(ctrl_out), save:: o_sss &
     2044        = ctrl_out([1, 10, 10, 1, 10, 10, 11, 11, 11, 11], 'SSS', &
     2045        "bulk sea-surface salinity", "ppt", '')
     2046
    19722047END MODULE phys_output_ctrlout_mod
  • LMDZ6/trunk/libf/phylmdiso/phys_output_mod.F90

    r3927 r3940  
    1 ! $Id: phys_output_mod.F90 3666 2020-04-20 10:13:34Z lfalletti $
     1! $Id: phys_output_mod.F90 3792 2021-01-04 17:01:25Z evignon $
    22!
    33
     
    4141    USE mod_phys_lmdz_para
    4242    !Martin
    43     USE surface_data, ONLY : ok_snow
     43    USE surface_data, ONLY : landice_opt
    4444    USE phys_output_ctrlout_mod
    4545    USE mod_grid_phy_lmdz, only: klon_glo,nbp_lon,nbp_lat
     
    379379    ENDIF
    380380
    381         write(*,*) 'phys_output_mid 344'
    382381!!! Declaration des axes verticaux de chaque fichier:
    383382    IF (prt_level >= 10) THEN
     
    452451   ENDIF
    453452#endif
    454         write(*,*) 'phys_output_mid 416'
    455453
    456454        IF (clef_files(iff)) THEN
  • LMDZ6/trunk/libf/phylmdiso/phys_output_var_mod.F90

    r3927 r3940  
    133133 !$OMP THREADPRIVATE(sens_prec_liq_o, sens_prec_sol_o,lat_prec_liq_o,lat_prec_sol_o)
    134134
     135  ! Ocean-atmosphere interface, subskin ocean and near-surface ocean:
     136 
     137  REAL, ALLOCATABLE, SAVE:: dter(:)
     138  ! Temperature variation in the diffusive microlayer, that is
     139  ! ocean-air interface temperature minus subskin temperature. In K.
     140     
     141  REAL, SAVE, ALLOCATABLE:: dser(:)
     142  ! Temperature variation in the diffusive microlayer, that is
     143  ! subskin temperature minus ocean-air interface temperature. In K.
     144
     145  REAL, SAVE, ALLOCATABLE:: tkt(:)
     146  ! épaisseur (m) de la couche de diffusion thermique (microlayer)
     147  ! cool skin thickness
     148
     149  REAL, SAVE, ALLOCATABLE:: tks(:)
     150  ! épaisseur (m) de la couche de diffusion de masse (microlayer)
     151 
     152  REAL, SAVE, ALLOCATABLE:: taur(:) ! momentum flux due to rain, in Pa
     153
     154  REAL, SAVE, ALLOCATABLE:: sss(:)
     155  ! bulk salinity of the surface layer of the ocean, in ppt
     156 
     157  !$OMP THREADPRIVATE(dter, dser, tkt, tks, taur, sss)
     158
    135159CONTAINS
    136160
     
    138162  SUBROUTINE phys_output_var_init
    139163    use dimphy
     164    use config_ocean_skin_m, only: activate_ocean_skin
    140165
    141166    IMPLICIT NONE
     
    191216    IF (ok_gwd_rando) allocate(zustr_gwd_rando(klon), zvstr_gwd_rando(klon))
    192217
     218    if (activate_ocean_skin >= 1) allocate(dter(klon), dser(klon), tkt(klon), &
     219         tks(klon), taur(klon), sss(klon))
     220
    193221  END SUBROUTINE phys_output_var_init
    194222
  • LMDZ6/trunk/libf/phylmdiso/phys_output_write_mod.F90

    r3927 r3940  
    11!
    2 ! $Id: phys_output_write_mod.F90 3692 2020-06-02 14:57:54Z oboucher $
     2! $Id: phys_output_write_mod.F90 3900 2021-05-17 13:35:58Z evignon $
    33!
    44MODULE phys_output_write_mod
     
    3838         o_t2m, o_t2m_min, o_t2m_max, &
    3939         o_t2m_min_mon, o_t2m_max_mon, &
     40         o_nt2mout, o_nt2moutfg, &
     41         o_nq2mout, o_nq2moutfg, &
     42         o_nu2mout, o_nu2moutfg, &
    4043         o_q2m, o_ustar, o_u10m, o_v10m, &
    4144         o_wind10m, o_wind10max, o_wind100m, o_gusts, o_sicf, &
     
    4548         o_snow, o_msnow, o_fsnow, o_evap, o_ep,o_epmax_diag, & ! epmax_cape
    4649         o_tops, o_tops0, o_topl, o_topl0, &
    47          o_SWupTOA, o_SWupTOAclr, o_SWupTOAcleanclr, o_SWdnTOA, &
     50         o_SWupTOA, o_SWupTOAclr, o_SWupTOAcleanclr, o_SWdnTOA, o_fdiffSWdnSFC, &
    4851         o_SWdnTOAclr, o_nettop, o_SWup200, &
    4952         o_SWup200clr, o_SWdn200, o_SWdn200clr, &
     
    8487         o_dtvdf_x    , o_dtvdf_w    , o_dqvdf_x    , o_dqvdf_w    , &
    8588         o_sens_x     , o_sens_w     , o_flat_x     , o_flat_w     , &
    86          o_delta_tsurf, &
     89         o_delta_tsurf, o_delta_tsurf_srf, &
    8790         o_cdragh_x   , o_cdragh_w   , o_cdragm_x   , o_cdragm_w   , &
    8891         o_kh         , o_kh_x       , o_kh_w       , &
     
    133136         o_vitu, o_vitv, o_vitw, o_pres, o_paprs, &
    134137         o_zfull, o_zhalf, o_rneb, o_rnebjn, o_rnebcon, &
    135          o_rnebls, o_rneblsvol, o_rhum, o_ozone, o_ozone_light, &
     138         o_rnebls, o_rneblsvol, o_rhum, o_rhl, o_rhi, o_ozone, o_ozone_light, &
    136139         o_duphy, o_dtphy, o_dqphy, o_dqphy2d, o_dqlphy, o_dqlphy2d, &
    137140         o_dqsphy, o_dqsphy2d, o_albe_srf, o_z0m_srf, o_z0h_srf, &
    138          o_ages_srf, o_snow_srf, o_alb1, o_alb2, o_tke, &
     141         o_ages_srf, o_snow_srf, o_alb1, o_alb2, o_tke, o_tke_dissip, &
    139142         o_tke_max, o_kz, o_kz_max, o_clwcon, &
    140143         o_dtdyn, o_dqdyn, o_dqdyn2d, o_dqldyn, o_dqldyn2d, &
     
    206209         o_col_O3_strato, o_col_O3_tropo,                 &
    207210!--interactive CO2
    208          o_flx_co2_ocean, o_flx_co2_land, o_flx_co2_ff, o_flx_co2_bb
    209 
     211         o_flx_co2_ocean, o_flx_co2_ocean_cor, &
     212         o_flx_co2_land, o_flx_co2_land_cor, &
     213         o_flx_co2_ff, o_flx_co2_bb, &
     214         o_delta_sst, o_delta_sal, o_ds_ns, o_dt_ns, o_dter, o_dser, o_tkt, &
     215         o_tks, o_taur, o_sss
    210216
    211217#ifdef CPP_StratAer
     
    230236         qsol, z0m, z0h, fevap, agesno, &
    231237         nday_rain, rain_con, snow_con, &
    232          topsw, toplw, toplw0, swup, swdn, &
     238         topsw, toplw, toplw0, swup, swdn, solswfdiff, &
    233239         topsw0, swupc0, swdnc0, swup0, swdn0, SWup200, SWup200clr, &
    234240         SWdn200, SWdn200clr, LWup200, LWup200clr, &
     
    252258         T2sumSTD, nlevSTD, du_gwd_rando, du_gwd_front, &
    253259         ulevSTD, vlevSTD, wlevSTD, philevSTD, qlevSTD, tlevSTD, &
    254          rhlevSTD, O3STD, O3daySTD, uvSTD, vqSTD, vTSTD, wqSTD, &
     260         rhlevSTD, O3STD, O3daySTD, uvSTD, vqSTD, vTSTD, wqSTD, vphiSTD, &
     261         wTSTD, u2STD, v2STD, T2STD, missing_val_nf90, delta_sal, ds_ns, &
    255262#ifdef ISO
    256263        xtrain_con, xtsnow_con, xtrain_fall, xtsnow_fall, &
    257264#endif
    258          vphiSTD, wTSTD, u2STD, v2STD, T2STD, missing_val_nf90
     265         dt_ns, delta_sst
    259266
    260267    USE phys_local_var_mod, ONLY: zxfluxlat, slp, ptstar, pt0, zxtsol, zt2m, &
    261          zt2m_cor,zq2m_cor,zu10m_cor,zv10m_cor, zrh2m_cor, zqsat2m_cor, &
    262          t2m_min_mon, t2m_max_mon, evap, &
    263          l_mixmin,l_mix, &
     268         zn2mout, t2m_min_mon, t2m_max_mon, evap, &
     269         l_mixmin,l_mix, tke_dissip, &
    264270         zu10m, zv10m, zq2m, zustar, zxqsurf, &
    265271         rain_lsc, rain_num, snow_lsc, bils, sens, fder, &
     
    268274         sissnow, runoff, albsol3_lic, evap_pot, &
    269275         t2m, fluxt, fluxlat, fsollw, fsolsw, &
    270          wfbils, wfbilo, wfevap, wfrain, wfsnow, & 
     276         wfbils, wfbilo, wfevap, wfrain, wfsnow, &
    271277         cdragm, cdragh, cldl, cldm, &
    272          cldh, cldt, JrNt, cldljn, cldmjn, cldhjn, &
    273          cldtjn, cldq, flwp, fiwp, ue, ve, uq, vq, &
     278         cldh, cldt, JrNt,   & ! only output names: cldljn,cldmjn,cldhjn,cldtjn
     279         cldq, flwp, fiwp, ue, ve, uq, vq, &
    274280         uwat, vwat, &
    275281         plcl, plfc, wbeff, convoccur, upwd, dnwd, dnwd0, prw, prlw, prsw, &
     
    307313         ql_seri, qs_seri, tr_seri, &
    308314         zphi, u_seri, v_seri, omega, cldfra, &
    309          rneb, rnebjn, rneblsvol, zx_rh, d_t_dyn,  &
     315         rneb, rnebjn, rneblsvol, zx_rh, zx_rhl, zx_rhi, d_t_dyn,  &
    310316         d_q_dyn,  d_ql_dyn, d_qs_dyn, &
    311317         d_q_dyn2d,  d_ql_dyn2d, d_qs_dyn2d, &
     
    349355
    350356    USE carbon_cycle_mod, ONLY: fco2_ff, fco2_bb, fco2_land, fco2_ocean
     357    USE carbon_cycle_mod, ONLY: fco2_ocean_cor, fco2_land_cor
    351358
    352359    USE phys_output_var_mod, ONLY: vars_defined, snow_o, zfra_o, bils_diss, &
     
    370377         alt_tropo, &
    371378!Ionela
    372          ok_4xCO2atm
     379         ok_4xCO2atm, dter, dser, tkt, tks, taur, sss
    373380
    374381    USE ocean_slab_mod, ONLY: nslay, tslab, slab_bilg, tice, seaice, &
     
    383390#endif
    384391    USE geometry_mod, ONLY: cell_area, latitude_deg, longitude_deg
    385     USE surface_data, ONLY: type_ocean, version_ocean, ok_veget, ok_snow
     392    USE surface_data, ONLY: type_ocean, version_ocean, ok_veget, landice_opt
    386393    USE aero_mod, ONLY: naero_tot, id_STRAT_phy
    387394    USE ioipsl, ONLY: histend, histsync
     
    402409#endif
    403410    USE tracinca_mod, ONLY: config_inca
     411    use config_ocean_skin_m, only: activate_ocean_skin
    404412
    405413    USE vertical_layers_mod, ONLY: presnivs
     
    471479    REAL,DIMENSION(klon)      :: zrho, zt
    472480
     481    INTEGER :: nqup
     482
    473483    ! On calcul le nouveau tau:
    474484    itau_w = itau_phy + itap
     
    484494#ifdef CPP_XIOS
    485495    CALL wxios_set_context
     496#endif
     497
     498#ifndef CPP_XIOS
     499    missing_val=missing_val_nf90
    486500#endif
    487501
     
    702716       CALL histwrite_phy(o_slp, slp)
    703717       CALL histwrite_phy(o_tsol, zxtsol)
    704        CALL histwrite_phy(o_t2m, zt2m_cor)
    705        CALL histwrite_phy(o_t2m_min, zt2m_cor)
    706        CALL histwrite_phy(o_t2m_max, zt2m_cor)
     718       CALL histwrite_phy(o_t2m, zt2m)
     719       CALL histwrite_phy(o_t2m_min, zt2m)
     720       CALL histwrite_phy(o_t2m_max, zt2m)
    707721       CALL histwrite_phy(o_t2m_max_mon, t2m_max_mon)
    708722       CALL histwrite_phy(o_t2m_min_mon, t2m_min_mon)
     
    710724       IF (vars_defined) THEN
    711725          DO i=1, klon
    712              zx_tmp_fi2d(i)=SQRT(zu10m_cor(i)*zu10m_cor(i)+zv10m_cor(i)*zv10m_cor(i))
     726             zx_tmp_fi2d(i)=real(zn2mout(i,1))
     727          ENDDO
     728       ENDIF
     729       CALL histwrite_phy(o_nt2mout, zx_tmp_fi2d)
     730
     731       IF (vars_defined) THEN
     732          DO i=1, klon
     733             zx_tmp_fi2d(i)=real(zn2mout(i,2))
     734          ENDDO
     735       ENDIF
     736       CALL histwrite_phy(o_nt2moutfg, zx_tmp_fi2d)
     737
     738       IF (vars_defined) THEN
     739          DO i=1, klon
     740             zx_tmp_fi2d(i)=real(zn2mout(i,3))
     741          ENDDO
     742       ENDIF
     743       CALL histwrite_phy(o_nq2mout, zx_tmp_fi2d)
     744
     745       IF (vars_defined) THEN
     746          DO i=1, klon
     747             zx_tmp_fi2d(i)=real(zn2mout(i,4))
     748          ENDDO
     749       ENDIF
     750       CALL histwrite_phy(o_nq2moutfg, zx_tmp_fi2d)
     751
     752       IF (vars_defined) THEN
     753          DO i=1, klon
     754             zx_tmp_fi2d(i)=real(zn2mout(i,5))
     755          ENDDO
     756       ENDIF
     757       CALL histwrite_phy(o_nu2mout, zx_tmp_fi2d)
     758
     759       IF (vars_defined) THEN
     760          DO i=1, klon
     761             zx_tmp_fi2d(i)=real(zn2mout(i,6))
     762          ENDDO
     763       ENDIF
     764       CALL histwrite_phy(o_nu2moutfg, zx_tmp_fi2d)
     765
     766       IF (vars_defined) THEN
     767          DO i=1, klon
     768             zx_tmp_fi2d(i)=SQRT(zu10m(i)*zu10m(i)+zv10m(i)*zv10m(i))
    713769          ENDDO
    714770       ENDIF
     
    717773       IF (vars_defined) THEN
    718774          DO i=1, klon
    719              zx_tmp_fi2d(i)=SQRT(zu10m_cor(i)*zu10m_cor(i)+zv10m_cor(i)*zv10m_cor(i))
     775             zx_tmp_fi2d(i)=SQRT(zu10m(i)*zu10m(i)+zv10m(i)*zv10m(i))
    720776          ENDDO
    721777       ENDIF
     
    725781
    726782       IF (vars_defined) THEN
    727           missing_val=missing_val_nf90
    728783          DO k = 1, kmax_100m                                      !--we could stop much lower
    729784            zrho(:) = pplay(:,k)/t_seri(:,k)/RD                    ! air density in kg/m3
     
    759814               zx_tmp_fi2d(i)=1.0
    760815             ELSE
    761                zx_tmp_fi2d(i)= 1.059e-09*x**10. - 1.351e-07*x**9. + 7.478e-06*x**8. - 0.0002352*x**7. + 0.004627*x**6. &
    762                              - 0.05898*x**5. + 0.4893*x**4. - 2.59*x**3. + 8.339*x**2. - 14.69*x + 10.73
     816               zx_tmp_fi2d(i)= 10.73 + x*(-14.69 + x*(8.339 + x*(-2.59 + x*(0.4893 + x*(-0.05898 + x*(0.004627 + &
     817                               x*(-0.0002352 + x*(7.478e-06 + x*(-1.351e-07 + x*(1.059e-09))))))))))
    763818               zx_tmp_fi2d(i)=MIN(MAX(zx_tmp_fi2d(i),0.0),1.0)
    764819             ENDIF
     
    780835               zx_tmp_fi2d(i)=1.0
    781836             ELSE
    782                zx_tmp_fi2d(i)= 3.352e-10*x**10. - 4.959e-08*x**9. + 3.195e-06*x**8. - 0.0001175*x**7. + 0.002716*x**6. &
    783                              - 0.04099*x**5. + 0.4065*x**4. - 2.601*x**3. + 10.25*x**2. - 22.39*x + 20.59
     837               zx_tmp_fi2d(i)= 20.59 + x*(-22.39 + x*(10.25 + x*(-2.601 + x*(0.4065 + x*(-0.04099 + x*(0.002716 + &
     838                               x*(-0.0001175 + x*(3.195e-06 + x*(-4.959e-08 + x*(3.352e-10))))))))))
    784839               zx_tmp_fi2d(i)=MIN(MAX(zx_tmp_fi2d(i),0.0),1.0)
    785840             ENDIF
     
    797852       ENDIF
    798853       CALL histwrite_phy(o_sicf, zx_tmp_fi2d)
    799        CALL histwrite_phy(o_q2m, zq2m_cor)
    800        CALL histwrite_phy(o_ustar, zustar)
    801        CALL histwrite_phy(o_u10m, zu10m_cor)
    802        CALL histwrite_phy(o_v10m, zv10m_cor)
     854       CALL histwrite_phy(o_q2m, zq2m)
     855       IF (vars_defined) zx_tmp_fi2d = zustar
     856       CALL histwrite_phy(o_ustar, zx_tmp_fi2d)
     857       CALL histwrite_phy(o_u10m, zu10m)
     858       CALL histwrite_phy(o_v10m, zv10m)
     859
    803860
    804861       IF (vars_defined) THEN
     
    9611018       ENDIF
    9621019       CALL histwrite_phy(o_SWdnSFCcleanclr, zx_tmp_fi2d)
     1020
     1021       CALL histwrite_phy(o_fdiffSWdnSFC, solswfdiff)
    9631022
    9641023       IF (vars_defined) THEN
     
    10221081       CALL histwrite_phy(o_tauy, zx_tmp_fi2d)
    10231082
    1024        IF (ok_snow) THEN
    1025           CALL histwrite_phy(o_snowsrf, snow_o)
    1026           CALL histwrite_phy(o_qsnow, qsnow)
    1027           CALL histwrite_phy(o_snowhgt,snowhgt)
    1028           CALL histwrite_phy(o_toice,to_ice)
    1029           CALL histwrite_phy(o_sissnow,sissnow)
    1030           CALL histwrite_phy(o_runoff,runoff)
    1031           CALL histwrite_phy(o_albslw3,albsol3_lic)
    1032        ENDIF
     1083       ! Etienne: test sorties pour compil sur JZ
     1084!       IF (landice_opt .GE. 1) THEN
     1085!          CALL histwrite_phy(o_snowsrf, snow_o)
     1086!          CALL histwrite_phy(o_qsnow, qsnow)
     1087!          CALL histwrite_phy(o_snowhgt,snowhgt)
     1088!          CALL histwrite_phy(o_toice,to_ice)
     1089!          CALL histwrite_phy(o_sissnow,sissnow)
     1090!          CALL histwrite_phy(o_runoff,runoff)
     1091!          CALL histwrite_phy(o_albslw3,albsol3_lic)
     1092!       ENDIF
    10331093
    10341094       DO nsrf = 1, nbsrf
     
    10801140             CALL histwrite_phy(o_l_mixmin(nsrf),  l_mixmin(:,1:klev,nsrf))
    10811141             CALL histwrite_phy(o_tke_max_srf(nsrf),  pbl_tke(:,1:klev,nsrf))
     1142
     1143                         
    10821144          ENDIF
    10831145!jyg<
     
    10901152!            ENDIF
    10911153
    1092 
    10931154       ENDDO
     1155       
     1156               
     1157        IF (iflag_pbl > 1) THEN
     1158          zx_tmp_fi3d=0.
     1159          IF (vars_defined) THEN
     1160             DO nsrf=1,nbsrf
     1161                DO k=1,klev
     1162                   zx_tmp_fi3d(:,k)=zx_tmp_fi3d(:,k) &
     1163                        +pctsrf(:,nsrf)*tke_dissip(:,k,nsrf)
     1164                ENDDO
     1165             ENDDO
     1166          ENDIF
     1167         
     1168          CALL histwrite_phy(o_tke_dissip, zx_tmp_fi3d)   
     1169       ENDIF
    10941170
    10951171       IF (vars_defined) zx_tmp_fi2d(1 : klon) = sens_prec_liq_o(1 : klon, 1)
     
    12151291       ! ATTENTION, LES ANCIENS HISTWRITE ONT ETES CONSERVES EN ATTENDANT MIEUX:
    12161292       ! Champs interpolles sur des niveaux de pression
    1217        missing_val=missing_val_nf90
    12181293       DO iff=1, nfiles
    12191294          ll=0
     
    13071382!
    13081383               CALL histwrite_phy(o_dqvdf_w    ,zx_tmp_fi3d)
    1309                CALL histwrite_phy(o_sens_x     ,sens_x     )
    1310                CALL histwrite_phy(o_sens_w     ,sens_w     )
     1384       IF (vars_defined)  zx_tmp_fi2d(1:klon)=-1*sens_x(1:klon)
     1385               CALL histwrite_phy(o_sens_x     ,zx_tmp_fi2d)
     1386       IF (vars_defined)  zx_tmp_fi2d(1:klon)=-1*sens_w(1:klon)
     1387               CALL histwrite_phy(o_sens_w     ,zx_tmp_fi2d)
    13111388               CALL histwrite_phy(o_flat_x     ,zxfluxlat_x)
    13121389               CALL histwrite_phy(o_flat_w     ,zxfluxlat_w)
    1313                CALL histwrite_phy(o_delta_tsurf,delta_tsurf)
     1390          zx_tmp_fi2d=0.
     1391          IF (vars_defined) THEN
     1392             DO nsrf=1,nbsrf
     1393                   zx_tmp_fi2d(:)=zx_tmp_fi2d(:) &
     1394                        +pctsrf(:,nsrf)*delta_tsurf(:,nsrf)
     1395             ENDDO
     1396          ENDIF
     1397               CALL histwrite_phy(o_delta_tsurf,zx_tmp_fi2d)
    13141398               CALL histwrite_phy(o_cdragh_x   ,cdragh_x   )
    13151399               CALL histwrite_phy(o_cdragh_w   ,cdragh_w   )
     
    14191503       IF (vars_defined) THEN
    14201504          DO i=1, klon
    1421              zx_tmp_fi2d(i)=MIN(100.,rh2m(i)*100.)
     1505             IF (zt2m(i).LE.273.15) then
     1506                zx_tmp_fi2d(i)=MAX(0.,rh2m(i)*100.)
     1507             ELSE
     1508                zx_tmp_fi2d(i)=MAX(0.,MIN(100.,rh2m(i)*100.))
     1509             ENDIF
    14221510          ENDDO
    14231511       ENDIF
     
    14381526!       CALL histwrite_phy(o_rh2m_max, zx_tmp_fi2d)
    14391527
    1440        CALL histwrite_phy(o_qsat2m, zqsat2m_cor)
     1528       CALL histwrite_phy(o_qsat2m, qsat2m)
    14411529       CALL histwrite_phy(o_tpot, tpot)
    14421530       CALL histwrite_phy(o_tpote, tpote)
     
    17171805       CALL histwrite_phy(o_rnebjn, zx_tmp_fi3d)
    17181806       CALL histwrite_phy(o_rhum, zx_rh)
     1807       IF (iflag_ice_thermo .GT. 0) THEN
     1808          IF (vars_defined) zx_tmp_fi3d = zx_rhl * 100.
     1809          CALL histwrite_phy(o_rhl, zx_tmp_fi3d)
     1810          IF (vars_defined) zx_tmp_fi3d = zx_rhi * 100.
     1811          CALL histwrite_phy(o_rhi, zx_tmp_fi3d)
     1812       ENDIF
     1813
    17191814       
    17201815       IF (vars_defined) zx_tmp_fi3d = wo(:, :, 1) * dobson_u * 1e3 / zmasse / rmo3 * rmd
     
    17751870          ENDIF
    17761871          CALL histwrite_phy(o_tke, zx_tmp_fi3d)
    1777 
    1778           CALL histwrite_phy(o_tke_max, zx_tmp_fi3d)
     1872          CALL histwrite_phy(o_tke_max, zx_tmp_fi3d) 
     1873
    17791874       ENDIF
    17801875
     
    21882283          CALL histwrite_phy(o_rldcs4co2, lwdn0p)
    21892284       ENDIF !ok_4xCO2atm
    2190         write(*,*) 'phys_output_write 2188'
    21912285!!!!!!!!!!!! Sorties niveaux de pression NMC !!!!!!!!!!!!!!!!!!!!
    21922286#ifdef CPP_IOIPSL
     
    21952289       ! ATTENTION, LES ANCIENS HISTWRITE ONT ETES CONSERVES EN ATTENDANT MIEUX:
    21962290       ! Champs interpolles sur des niveaux de pression
    2197        missing_val=missing_val_nf90
    21982291       DO iff=7, nfiles-1 !--OB: here we deal with files 7,8,9
    21992292
     
    23292422  ENDIF
    23302423#endif
    2331         write(*,*) 'phys_output_write 2331'
    23322424!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    23332425       IF (iflag_phytrac == 1 ) then
     
    23892481           CALL histwrite_phy(o_flx_co2_land,  fco2_land)
    23902482           CALL histwrite_phy(o_flx_co2_ocean, fco2_ocean)
     2483           CALL histwrite_phy(o_flx_co2_ocean_cor, fco2_ocean_cor)
     2484           CALL histwrite_phy(o_flx_co2_land_cor, fco2_land_cor)
    23912485           CALL histwrite_phy(o_flx_co2_ff,    fco2_ff)
    23922486           CALL histwrite_phy(o_flx_co2_bb,    fco2_bb)
    23932487         ENDIF !--type_trac co2i
    23942488
     2489         IF (type_trac == 'inco') THEN
     2490           nqup = nqo+1
     2491           DO iq=nqo+1, nqup
     2492             !--3D fields
     2493             CALL histwrite_phy(o_trac(iq-nqo), tr_seri(:,:,iq-nqo))
     2494             CALL histwrite_phy(o_dtr_vdf(iq-nqo),d_tr_cl(:,:,iq-nqo))
     2495             CALL histwrite_phy(o_dtr_the(iq-nqo),d_tr_th(:,:,iq-nqo))
     2496             CALL histwrite_phy(o_dtr_con(iq-nqo),d_tr_cv(:,:,iq-nqo))
     2497             !--2D fields
     2498             !--CO2 burden
     2499             zx_tmp_fi2d=0.
     2500             IF (vars_defined) THEN
     2501                DO k=1,klev
     2502                   zx_tmp_fi2d(:)=zx_tmp_fi2d(:)+zmasse(:,k)*tr_seri(:,k,iq-nqo)
     2503                ENDDO
     2504             ENDIF
     2505             CALL histwrite_phy(o_trac_cum(iq-nqo), zx_tmp_fi2d)
     2506           ENDDO !--iq
     2507           !--CO2 net fluxes
     2508           CALL histwrite_phy(o_flx_co2_land,  fco2_land)
     2509           CALL histwrite_phy(o_flx_co2_ocean, fco2_ocean)
     2510           CALL histwrite_phy(o_flx_co2_ocean_cor, fco2_ocean_cor)
     2511           CALL histwrite_phy(o_flx_co2_land_cor, fco2_land_cor)
     2512           CALL histwrite_phy(o_flx_co2_ff,    fco2_ff)
     2513           CALL histwrite_phy(o_flx_co2_bb,    fco2_bb)
     2514         ENDIF !--type_trac inco
     2515
    23952516       ENDIF   !(iflag_phytrac==1)
    23962517
     2518       if (activate_ocean_skin >= 1) then
     2519          CALL histwrite_phy(o_delta_sst, delta_sst)
     2520          CALL histwrite_phy(o_delta_sal, delta_sal)
     2521          CALL histwrite_phy(o_ds_ns, ds_ns)
     2522          CALL histwrite_phy(o_dt_ns, dt_ns)
     2523          CALL histwrite_phy(o_dter, dter)
     2524          CALL histwrite_phy(o_dser, dser)
     2525          CALL histwrite_phy(o_tkt, tkt)
     2526          CALL histwrite_phy(o_tks, tks)
     2527          CALL histwrite_phy(o_taur, taur)
     2528          CALL histwrite_phy(o_sss, sss)
     2529       end if
    23972530
    23982531#ifdef ISO
  • LMDZ6/trunk/libf/phylmdiso/phys_state_var_mod.F90

    r3927 r3940  
    11!
    2 ! $Id: phys_state_var_mod.F90 3496 2019-05-10 10:17:35Z jyg $
     2! $Id: phys_state_var_mod.F90 3888 2021-05-05 10:50:37Z jyg $
    33!
    44      MODULE phys_state_var_mod
     
    3232      REAL, ALLOCATABLE, SAVE :: ftsol(:,:)
    3333!$OMP THREADPRIVATE(ftsol)
     34      REAL, ALLOCATABLE, SAVE :: beta_aridity(:,:)
     35!$OMP THREADPRIVATE(beta_aridity)
    3436      REAL,ALLOCATABLE,SAVE :: qsol(:),fevap(:,:),z0m(:,:),z0h(:,:),agesno(:,:)
    3537!$OMP THREADPRIVATE(qsol,fevap,z0m,z0h,agesno)
     
    4749!albedo SB >>>
    4850      REAL, ALLOCATABLE, SAVE :: falb_dif(:,:,:), falb_dir(:,:,:)
    49       real, allocatable, save :: chl_con(:)
     51      REAL, ALLOCATABLE, SAVE :: chl_con(:)
    5052!$OMP THREADPRIVATE(falb_dir,falb_dif,chl_con)
    5153!albedo SB <<<
     
    106108      REAL, ALLOCATABLE, SAVE :: coefm(:,:,:) ! Kz momentum
    107109!$OMP THREADPRIVATE(pbl_tke, coefh,coefm)
    108 !nrlmd<
    109       REAL, ALLOCATABLE, SAVE :: delta_tsurf(:,:) ! Surface temperature difference inside-outside cold pool
    110 !$OMP THREADPRIVATE(delta_tsurf)
    111 !>nrlmd
    112110      REAL, ALLOCATABLE, SAVE :: zmax0(:), f0(:) !
    113111!$OMP THREADPRIVATE(zmax0,f0)
     
    294292      REAL,ALLOCATABLE,SAVE :: wake_delta_pbl_TKE(:,:,:)
    295293!$OMP THREADPRIVATE(wake_delta_pbl_TKE)
     294!nrlmd<
     295      REAL, ALLOCATABLE, SAVE :: delta_tsurf(:,:) ! Surface temperature difference inside-outside cold pool
     296!$OMP THREADPRIVATE(delta_tsurf)
     297!>nrlmd
    296298!>jyg
    297299!
     
    442444!$OMP THREADPRIVATE(ccm)
    443445
    444 !!! nrlmd le 10/04/2012
    445446      REAL,SAVE,ALLOCATABLE :: ale_bl_trig(:)
    446447!$OMP THREADPRIVATE(ale_bl_trig)
    447 !!! fin nrlmd le 10/04/2012
     448
     449      REAL,SAVE,ALLOCATABLE :: ratqs_inter(:,:)
     450!$OMP THREADPRIVATE(ratqs_inter)
    448451
    449452#ifdef ISO
     
    462465!$OMP THREADPRIVATE(is_initialized)   
    463466
    464 CONTAINS
     467      ! Ocean-atmosphere interface:
     468
     469      REAL, ALLOCATABLE, SAVE:: ds_ns(:) ! (klon)
     470      ! "delta salinity near surface". Salinity variation in the
     471      ! near-surface turbulent layer. That is subskin salinity minus
     472      ! foundation salinity. In ppt.
     473
     474      REAL, ALLOCATABLE, SAVE:: dt_ns(:) ! (klon)
     475      ! "delta temperature near surface". Temperature variation in the
     476      ! near-surface turbulent layer. That is subskin temperature
     477      ! minus foundation temperature. (Can be negative.) In K.
     478     
     479      REAL, ALLOCATABLE, SAVE:: delta_sst(:) ! (klon)
     480      ! Ocean-air interface temperature minus bulk SST, in
     481      ! K. Allocated and defined only if activate_ocean_skin >= 1.
     482
     483      REAL, ALLOCATABLE, SAVE:: delta_sal(:) ! (klon)
     484      ! Ocean-air interface salinity minus bulk salinity, in ppt
     485     
     486      !$OMP THREADPRIVATE(delta_sal, ds_ns, dt_ns, delta_sst)
     487
     488    CONTAINS
    465489
    466490!======================================================================
     
    473497#endif
    474498USE indice_sol_mod
     499use config_ocean_skin_m, only: activate_ocean_skin
    475500IMPLICIT NONE
    476501
     
    485510include "clesphys.h"
    486511
     512      print*, 'is_initialized', is_initialized
    487513      IF (is_initialized) RETURN
    488514      is_initialized=.TRUE.
    489515      ALLOCATE(pctsrf(klon,nbsrf))
    490516      ALLOCATE(ftsol(klon,nbsrf))
     517      ALLOCATE(beta_aridity(klon,nbsrf))
    491518      ALLOCATE(qsol(klon),fevap(klon,nbsrf))
    492519      ALLOCATE(z0m(klon,nbsrf+1),z0h(klon,nbsrf+1),agesno(klon,nbsrf))
     
    496523      ALLOCATE(falb2(klon,nbsrf))
    497524!albedo SB >>>
     525      print*, 'allocate falb'
    498526      ALLOCATE(falb_dir(klon,nsw,nbsrf),falb_dif(klon,nsw,nbsrf))
     527      print*, 'allocate falb good', falb_dir(1,1,1)
    499528      ALLOCATE(chl_con(klon))
    500529!albedo SB <<<
     
    675704#endif     
    676705#endif
    677 !!! nrlmd le 10/04/2012
     706
    678707      ALLOCATE(ale_bl_trig(klon))
    679 !!! fin nrlmd le 10/04/2012
     708      ALLOCATE(ratqs_inter(klon,klev))
    680709      IF (ok_gwd_rando) THEN
    681         allocate(du_gwd_rando(klon, klev))
     710        ALLOCATE(du_gwd_rando(klon, klev))
    682711        du_gwd_rando(:,:)=0.
    683712      ENDIF
     
    686715        du_gwd_front(:,:) = 0 !ym missing init   
    687716      ENDIF
    688 END SUBROUTINE phys_state_var_init
     717      if (activate_ocean_skin >= 1) ALLOCATE(delta_sal(klon), ds_ns(klon), &
     718           dt_ns(klon), delta_sst(klon))
     719
     720    END SUBROUTINE phys_state_var_init
    689721
    690722!======================================================================
    691 SUBROUTINE phys_state_var_end
     723    SUBROUTINE phys_state_var_end
     724      ! Useful only for lmdz1d.
    692725!USE dimphy
    693726USE indice_sol_mod
     727use config_ocean_skin_m, only: activate_ocean_skin
    694728IMPLICIT NONE
    695729include "clesphys.h"
    696730
    697       deallocate(pctsrf, ftsol, falb1, falb2)
    698       deallocate(qsol,fevap,z0m,z0h,agesno)
     731      DEALLOCATE(pctsrf, ftsol, falb1, falb2)
     732      DEALLOCATE(beta_aridity)
     733      DEALLOCATE(qsol,fevap,z0m,z0h,agesno)
    699734!FC
    700       deallocate(treedrg)
    701       deallocate(rain_fall, snow_fall, solsw, solswfdiff, sollw, radsol, swradcorr)
    702       deallocate(zmea, zstd, zsig, zgam)
    703       deallocate(zthe, zpic, zval)
    704       deallocate(rugoro, t_ancien, q_ancien, clwcon, rnebcon)
    705       deallocate(qs_ancien, ql_ancien)
    706       deallocate(prw_ancien, prlw_ancien, prsw_ancien)
    707       deallocate(qtc_cv,sigt_cv)
    708       deallocate(u_ancien, v_ancien)
    709       deallocate(tr_ancien)                           !RomP
    710       deallocate(ratqs, pbl_tke,coefh,coefm)
     735      DEALLOCATE(treedrg)
     736      DEALLOCATE(rain_fall, snow_fall, solsw, solswfdiff, sollw, radsol, swradcorr)
     737      DEALLOCATE(zmea, zstd, zsig, zgam)
     738      DEALLOCATE(zthe, zpic, zval)
     739      DEALLOCATE(rugoro, t_ancien, q_ancien, clwcon, rnebcon)
     740      DEALLOCATE(qs_ancien, ql_ancien)
     741      DEALLOCATE(prw_ancien, prlw_ancien, prsw_ancien)
     742      DEALLOCATE(qtc_cv,sigt_cv)
     743      DEALLOCATE(u_ancien, v_ancien)
     744      DEALLOCATE(tr_ancien)                           !RomP
     745      DEALLOCATE(ratqs, pbl_tke,coefh,coefm)
     746      DEALLOCATE(zmax0, f0)
     747      DEALLOCATE(sig1, w01)
     748      DEALLOCATE(entr_therm, fm_therm)
     749      DEALLOCATE(detr_therm)
     750      DEALLOCATE(clwcon0th, rnebcon0th)
     751! radiation outputs
     752      DEALLOCATE(swdnc0, swdn0, swdn)
     753      DEALLOCATE(swupc0, swup0, swup)
     754      DEALLOCATE(lwdnc0, lwdn0, lwdn)
     755      DEALLOCATE(lwupc0, lwup0, lwup)
     756      DEALLOCATE(SWdn200clr, SWdn200)
     757      DEALLOCATE(SWup200clr, SWup200)
     758      DEALLOCATE(LWdn200clr, LWdn200)
     759      DEALLOCATE(LWup200clr, LWup200)
     760      DEALLOCATE(LWdnTOA, LWdnTOAclr)
     761! pressure level
     762      DEALLOCATE(tsumSTD)
     763      DEALLOCATE(usumSTD, vsumSTD)
     764      DEALLOCATE(wsumSTD, phisumSTD)
     765      DEALLOCATE(tnondef)
     766      DEALLOCATE(qsumSTD, rhsumSTD)
     767      DEALLOCATE(uvsumSTD)
     768      DEALLOCATE(vqsumSTD)
     769      DEALLOCATE(vTsumSTD)
     770      DEALLOCATE(wqsumSTD)
     771      DEALLOCATE(vphisumSTD)
     772      DEALLOCATE(wTsumSTD)
     773      DEALLOCATE(u2sumSTD)
     774      DEALLOCATE(v2sumSTD)
     775      DEALLOCATE(T2sumSTD)
     776      DEALLOCATE(O3sumSTD)
     777      DEALLOCATE(O3daysumSTD)
     778!IM beg
     779      DEALLOCATE(wlevSTD,ulevSTD,vlevSTD,tlevSTD,qlevSTD,rhlevSTD,philevSTD)
     780      DEALLOCATE(uvSTD,vqSTD,vTSTD,wqSTD,vphiSTD,wTSTD,u2STD,v2STD,T2STD,O3STD,O3daySTD)
     781!IM end
     782      DEALLOCATE(seed_old)
     783      DEALLOCATE(zuthe, zvthe)
     784      DEALLOCATE(alb_neig)
     785      DEALLOCATE(ema_cbmf)
     786      DEALLOCATE(ema_pcb, ema_pct)
     787      DEALLOCATE(Mipsh, Ma, qcondc)
     788      DEALLOCATE(wd, sigd)
     789      DEALLOCATE(cin, ALE, ALP)
     790      DEALLOCATE(ftd, fqd)
     791      DEALLOCATE(Ale_bl, Alp_bl)
     792      DEALLOCATE(ale_wake)
     793      DEALLOCATE(ale_bl_stat)
     794      DEALLOCATE(lalim_conv, wght_th)
     795      DEALLOCATE(wake_deltat, wake_deltaq)
     796      DEALLOCATE(wake_s, awake_dens, wake_dens)
     797      DEALLOCATE(wake_Cstar, wake_pe, wake_fip)
     798!jyg<
     799      DEALLOCATE(wake_delta_pbl_TKE)
    711800!nrlmd<
    712       deallocate(delta_tsurf)
     801      DEALLOCATE(delta_tsurf)
    713802!>nrlmd
    714       deallocate(zmax0, f0)
    715       deallocate(sig1, w01)
    716       deallocate(entr_therm, fm_therm)
    717       deallocate(detr_therm)
    718       deallocate(clwcon0th, rnebcon0th)
    719 ! radiation outputs
    720       deallocate(swdnc0, swdn0, swdn)
    721       deallocate(swupc0, swup0, swup)
    722       deallocate(lwdnc0, lwdn0, lwdn)
    723       deallocate(lwupc0, lwup0, lwup)
    724       deallocate(SWdn200clr, SWdn200)
    725       deallocate(SWup200clr, SWup200)
    726       deallocate(LWdn200clr, LWdn200)
    727       deallocate(LWup200clr, LWup200)
    728       deallocate(LWdnTOA, LWdnTOAclr)
    729 ! pressure level
    730       deallocate(tsumSTD)
    731       deallocate(usumSTD, vsumSTD)
    732       deallocate(wsumSTD, phisumSTD)
    733       deallocate(tnondef)
    734       deallocate(qsumSTD, rhsumSTD)
    735       deallocate(uvsumSTD)
    736       deallocate(vqsumSTD)
    737       deallocate(vTsumSTD)
    738       deallocate(wqsumSTD)
    739       deallocate(vphisumSTD)
    740       deallocate(wTsumSTD)
    741       deallocate(u2sumSTD)
    742       deallocate(v2sumSTD)
    743       deallocate(T2sumSTD)
    744       deallocate(O3sumSTD)
    745       deallocate(O3daysumSTD)
    746 !IM beg
    747       deallocate(wlevSTD,ulevSTD,vlevSTD,tlevSTD,qlevSTD,rhlevSTD,philevSTD)
    748       deallocate(uvSTD,vqSTD,vTSTD,wqSTD,vphiSTD,wTSTD,u2STD,v2STD,T2STD,O3STD,O3daySTD)
    749 !IM end
    750       deallocate(seed_old)
    751       deallocate(zuthe, zvthe)
    752       deallocate(alb_neig)
    753       deallocate(ema_cbmf)
    754       deallocate(ema_pcb, ema_pct)
    755       deallocate(Mipsh, Ma, qcondc)
    756       deallocate(wd, sigd)
    757       deallocate(cin, ALE, ALP)
    758       deallocate(ftd, fqd)
    759       deallocate(Ale_bl, Alp_bl)
    760       deallocate(ale_wake)
    761       deallocate(ale_bl_stat)
    762       deallocate(lalim_conv, wght_th)
    763       deallocate(wake_deltat, wake_deltaq)
    764       deallocate(wake_s, awake_dens, wake_dens)
    765       deallocate(wake_Cstar, wake_pe, wake_fip)
    766 !jyg<
    767       deallocate(wake_delta_pbl_TKE)
    768803!>jyg
    769       deallocate(pfrac_impa, pfrac_nucl)
    770       deallocate(pfrac_1nucl)
    771       deallocate(total_rain, nday_rain)
    772       deallocate(paire_ter)
    773       deallocate(albsol1, albsol2)
     804      DEALLOCATE(pfrac_impa, pfrac_nucl)
     805      DEALLOCATE(pfrac_1nucl)
     806      DEALLOCATE(total_rain, nday_rain)
     807      DEALLOCATE(paire_ter)
     808      DEALLOCATE(albsol1, albsol2)
    774809!albedo SB >>>
    775       deallocate(albsol_dir,albsol_dif,falb_dir,falb_dif,chl_con)
     810      DEALLOCATE(albsol_dir,albsol_dif,falb_dir,falb_dif,chl_con)
    776811!albedo SB <<<
    777       deallocate(wo)
    778       deallocate(clwcon0,rnebcon0)
    779       deallocate(heat, heat0)
    780       deallocate(cool, cool0)
    781       deallocate(heat_volc, cool_volc)
    782       deallocate(topsw, toplw)
    783       deallocate(sollwdown, sollwdownclr)
    784       deallocate(gustiness)
    785       deallocate(toplwdown, toplwdownclr)
    786       deallocate(topsw0,toplw0,solsw0,sollw0)
    787       deallocate(albpla)
     812      DEALLOCATE(wo)
     813      DEALLOCATE(clwcon0,rnebcon0)
     814      DEALLOCATE(heat, heat0)
     815      DEALLOCATE(cool, cool0)
     816      DEALLOCATE(heat_volc, cool_volc)
     817      DEALLOCATE(topsw, toplw)
     818      DEALLOCATE(sollwdown, sollwdownclr)
     819      DEALLOCATE(gustiness)
     820      DEALLOCATE(toplwdown, toplwdownclr)
     821      DEALLOCATE(topsw0,toplw0,solsw0,sollw0)
     822      DEALLOCATE(albpla)
    788823!IM ajout variables CFMIP2/CMIP5
    789       deallocate(heatp, coolp)
    790       deallocate(heat0p, cool0p)
    791       deallocate(radsolp, topswp, toplwp)
    792       deallocate(albplap)
    793       deallocate(solswp, solswfdiffp, sollwp)
    794       deallocate(sollwdownp)
    795       deallocate(topsw0p,toplw0p)
    796       deallocate(solsw0p,sollw0p)
    797       deallocate(lwdnc0p, lwdn0p, lwdnp)
    798       deallocate(lwupc0p, lwup0p, lwupp)
    799       deallocate(swdnc0p, swdn0p, swdnp)
    800       deallocate(swupc0p, swup0p, swupp)
    801       deallocate(cape)
    802       deallocate(pbase,bbase)
    803       deallocate(zqasc)
    804       deallocate(ibas_con, itop_con)
    805       deallocate(rain_con, snow_con)
    806       deallocate(rlonPOS)
    807       deallocate(newsst)
    808       deallocate(ustar,u10m, v10m,wstar)
    809       deallocate(topswad, solswad)
    810       deallocate(topswai, solswai)
    811       deallocate(tau_aero,piz_aero,cg_aero)
    812       deallocate(tau_aero_sw_rrtm,piz_aero_sw_rrtm,cg_aero_sw_rrtm)
    813       deallocate(tau_aero_lw_rrtm,piz_aero_lw_rrtm,cg_aero_lw_rrtm)
    814       deallocate(ccm)
    815       if (ok_gwd_rando) deallocate(du_gwd_rando)
    816       if (.not. ok_hines .and. ok_gwd_rando) deallocate(du_gwd_front)
    817        
    818 !!! nrlmd le 10/04/2012
    819       deallocate(ale_bl_trig)
    820 !!! fin nrlmd le 10/04/2012
     824      DEALLOCATE(heatp, coolp)
     825      DEALLOCATE(heat0p, cool0p)
     826      DEALLOCATE(radsolp, topswp, toplwp)
     827      DEALLOCATE(albplap)
     828      DEALLOCATE(solswp, solswfdiffp, sollwp)
     829      DEALLOCATE(sollwdownp)
     830      DEALLOCATE(topsw0p,toplw0p)
     831      DEALLOCATE(solsw0p,sollw0p)
     832      DEALLOCATE(lwdnc0p, lwdn0p, lwdnp)
     833      DEALLOCATE(lwupc0p, lwup0p, lwupp)
     834      DEALLOCATE(swdnc0p, swdn0p, swdnp)
     835      DEALLOCATE(swupc0p, swup0p, swupp)
     836      DEALLOCATE(cape)
     837      DEALLOCATE(pbase,bbase)
     838      DEALLOCATE(zqasc)
     839      DEALLOCATE(ibas_con, itop_con)
     840      DEALLOCATE(rain_con, snow_con)
     841      DEALLOCATE(rlonPOS)
     842      DEALLOCATE(newsst)
     843      DEALLOCATE(ustar,u10m, v10m,wstar)
     844      DEALLOCATE(topswad, solswad)
     845      DEALLOCATE(topswai, solswai)
     846      DEALLOCATE(tau_aero,piz_aero,cg_aero)
     847      DEALLOCATE(tau_aero_sw_rrtm,piz_aero_sw_rrtm,cg_aero_sw_rrtm)
     848      DEALLOCATE(tau_aero_lw_rrtm,piz_aero_lw_rrtm,cg_aero_lw_rrtm)
     849      DEALLOCATE(ccm)
     850      if (ok_gwd_rando) DEALLOCATE(du_gwd_rando)
     851      if (.not. ok_hines .and. ok_gwd_rando) DEALLOCATE(du_gwd_front)
     852      DEALLOCATE(ale_bl_trig)
     853      DEALLOCATE(ratqs_inter)
     854
     855      if (activate_ocean_skin >= 1) deALLOCATE(delta_sal, ds_ns, dt_ns, &
     856           delta_sst)
    821857
    822858#ifdef ISO   
    823       deallocate(xtsol,fxtevap) 
    824       deallocate(xt_ancien,xtl_ancien,xts_ancien, fxtd, wake_deltaxt)
    825       deallocate(xtrain_fall, xtsnow_fall, xtrain_con, xtsnow_con)
     859      DEALLOCATE(xtsol,fxtevap) 
     860      DEALLOCATE(xt_ancien,xtl_ancien,xts_ancien, fxtd, wake_deltaxt)
     861      DEALLOCATE(xtrain_fall, xtsnow_fall, xtrain_con, xtsnow_con)
    826862#ifdef ISOTRAC
    827       deallocate(bassin_map,boite_map)
     863      DEALLOCATE(bassin_map,boite_map)
    828864#endif       
    829865#endif
    830866      is_initialized=.FALSE.
     867     
    831868END SUBROUTINE phys_state_var_end
    832869
  • LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90

    r3927 r3940  
    11!
    2 ! $Id: physiq_mod.F90 3666 2020-04-20 10:13:34Z lfalletti $
     2! $Id: physiq_mod.F90 3908 2021-05-20 07:11:13Z idelkadi $
    33!
    44!#define IO_DEBUG
     
    1616       d_u, d_v, d_t, d_qx, d_ps)
    1717
     18! For clarity, the "USE" section is now arranged in alphabetical order,
     19! with a separate section for CPP keys
     20! PLEASE try to follow this rule
     21
     22    USE ACAMA_GWD_rando_m, only: ACAMA_GWD_rando
     23    USE aero_mod
     24    USE add_phys_tend_mod, only : add_pbl_tend, add_phys_tend, diag_phys_tend, prt_enerbil, &
     25  &      fl_ebil, fl_cor_ebil
    1826    USE assert_m, only: assert
     27    USE change_srf_frac_mod
     28    USE conf_phys_m, only: conf_phys
     29    USE carbon_cycle_mod, ONLY : infocfields_init, RCO2_glo, carbon_cycle_rad
     30    USE CFMIP_point_locations   ! IM stations CFMIP
     31    USE cmp_seri_mod
     32    USE dimphy
     33    USE etat0_limit_unstruct_mod
     34    USE FLOTT_GWD_rando_m, only: FLOTT_GWD_rando
     35    USE fonte_neige_mod, ONLY  : fonte_neige_get_vars
     36    USE geometry_mod, ONLY: cell_area, latitude_deg, longitude_deg
    1937    USE ioipsl, only: histbeg, histvert, histdef, histend, histsync, &
    2038         histwrite, ju2ymds, ymds2ju, getin
    21     USE geometry_mod, ONLY: cell_area, latitude_deg, longitude_deg
     39    USE ioipsl_getin_p_mod, ONLY : getin_p
     40    USE indice_sol_mod
     41    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, type_trac,ok_isotopes, &
     42        nqtottr,itr_indice ! C Risi
     43
     44    USE iophy
     45    USE limit_read_mod, ONLY : init_limit_read
     46    USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, nbp_lev, klon_glo, grid1dTo2d_glo, grid_type, unstructured
     47    USE mod_phys_lmdz_mpi_data, only: is_mpi_root
     48    USE mod_phys_lmdz_para
     49    USE netcdf95, only: nf95_close
     50    USE netcdf, only: nf90_fill_real     ! IM for NMC files
     51    USE open_climoz_m, only: open_climoz ! ozone climatology from a file
     52    USE ozonecm_m, only: ozonecm ! ozone of J.-F. Royer
     53    USE pbl_surface_mod, ONLY : pbl_surface
     54    USE phyaqua_mod, only: zenang_an
     55    USE phystokenc_mod, ONLY: offline, phystokenc
    2256    USE phys_cal_mod, only: year_len, mth_len, days_elapsed, jh_1jan, &
    2357         year_cur, mth_cur,jD_cur, jH_cur, jD_ref, day_cur, hour
     58!!  USE phys_local_var_mod, ONLY : a long list of variables
     59!!              ==> see below, after "CPP Keys" section
     60    USE phys_state_var_mod ! Variables sauvegardees de la physique
     61    USE phys_output_mod
     62    USE phys_output_ctrlout_mod
     63    USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
     64    USE readaerosol_mod, ONLY : init_aero_fromfile
     65    USE readaerosolstrato_m, ONLY : init_readaerosolstrato
     66    USE radlwsw_m, only: radlwsw
     67    USE regr_horiz_time_climoz_m, ONLY: regr_horiz_time_climoz
     68    USE regr_pr_time_av_m, only: regr_pr_time_av
     69    USE surface_data,     ONLY : type_ocean, ok_veget, landice_opt
     70    USE time_phylmdz_mod, only: annee_ref, current_time, day_ini, day_ref, &
     71          day_step_phy, itau_phy, pdtphys, raz_date, start_time, update_time
     72    USE tracinca_mod, ONLY: config_inca
     73    USE tropopause_m,     ONLY: dyn_tropopause
     74    USE vampir
     75    USE VERTICAL_LAYERS_MOD, ONLY: aps,bps, ap, bp
    2476    USE write_field_phy
    25     USE dimphy
    26 USE infotrac_phy, ONLY: nqtot, nbtr, nqo, type_trac,ok_isotopes, &
    27         nqtottr,itr_indice ! C Risi
     77
     78    !USE cmp_seri_mod
     79!    USE add_phys_tend_mod, only : add_pbl_tend, add_phys_tend, diag_phys_tend, prt_enerbil, &
     80!  &      fl_ebil, fl_cor_ebil
     81
     82!!!!!!!!!!!!!!!!!! "USE" section for CPP keys !!!!!!!!!!!!!!!!!!!!!!!!
     83!
     84!
     85#ifdef CPP_Dust
     86    USE phytracr_spl_mod, ONLY: phytracr_spl, phytracr_spl_out_init
     87    USE phys_output_write_spl_mod
     88#else
     89    USE phytrac_mod, ONLY : phytrac_init, phytrac
     90    USE phys_output_write_mod
     91#endif
     92
     93
     94#ifdef REPROBUS
     95    USE CHEM_REP, ONLY : Init_chem_rep_xjour, &
     96         d_q_rep,d_ql_rep,d_qi_rep,ptrop,ttrop, &
     97         ztrop, gravit,itroprep, Z1,Z2,fac,B
     98#endif
     99
     100
     101#ifdef CPP_RRTM
     102    USE YOERAD, ONLY : NRADLP
     103    USE YOESW, ONLY : RSUN
     104#endif
     105
     106
     107#ifdef CPP_StratAer
     108    USE strataer_mod, ONLY: strataer_init
     109#endif
     110
     111
     112#ifdef CPP_XIOS
     113    USE xios, ONLY: xios_update_calendar, xios_context_finalize, &
     114            xios_get_field_attr, xios_field_is_active
     115    USE wxios, ONLY: missing_val, missing_val_omp
     116#endif
     117#ifndef CPP_XIOS
     118    USE paramLMDZ_phy_mod
     119#endif
     120
     121
    28122#ifdef ISO
    29123    USE infotrac_phy, ONLY:  &
     
    64158#endif
    65159#endif
    66     USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, nbp_lev, klon_glo, grid1dTo2d_glo, grid_type, unstructured
    67     USE mod_phys_lmdz_para
    68     USE iophy
    69     USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
    70     USE phystokenc_mod, ONLY: offline, phystokenc
    71     USE time_phylmdz_mod, only: raz_date, day_step_phy, update_time,current_time
    72     USE vampir
    73     USE pbl_surface_mod, ONLY : pbl_surface
    74     USE change_srf_frac_mod
    75     USE surface_data,     ONLY : type_ocean, ok_veget, ok_snow
    76     USE tropopause_m,     ONLY: dyn_tropopause
    77 #ifdef CPP_Dust
    78     USE phytracr_spl_mod, ONLY: phytracr_spl
    79 #endif
    80 #ifdef CPP_StratAer
    81     USE strataer_mod, ONLY: strataer_init
    82 #endif
    83     USE phys_local_var_mod, ONLY: phys_local_var_init, phys_local_var_end, &
     160
     161!
     162!
     163!!!!!!!!!!!!!!!!!!  END "USE" for CPP keys !!!!!!!!!!!!!!!!!!!!!!
     164
     165USE phys_local_var_mod, ONLY: phys_local_var_init, phys_local_var_end, &
    84166       ! [Variables internes non sauvegardees de la physique]
    85167       ! Variables locales pour effectuer les appels en serie
     
    160242       cdragm, cdragh,                   &
    161243       zustar, zu10m, zv10m, rh2m, qsat2m, &
    162        zq2m, zt2m, weak_inversion, &
    163        zq2m_cor,zt2m_cor,zu10m_cor,zv10m_cor, & ! pour corriger d'un bug
    164        zrh2m_cor,zqsat2m_cor, &
     244       zq2m, zt2m, zn2mout, weak_inversion, &
    165245       zt2m_min_mon, zt2m_max_mon,   &         ! pour calcul_divers.h
    166246       t2m_min_mon, t2m_max_mon,  &            ! pour calcul_divers.h
     
    175255       zxrunofflic,                            &
    176256       zxtsol, snow_lsc, zxfqfonte, zxqsurf,   &
     257       delta_qsurf,                            &
    177258       rain_lsc, rain_num,                     &
    178259       !
     
    236317       ref_liq, ref_ice, theta,  &
    237318       ref_liq_pi, ref_ice_pi,  &
    238        zphi, zx_rh, &
     319       zphi, zx_rh, zx_rhl, zx_rhi, &
    239320       pmfd, pmfu,  &
    240321       !
     
    284365#endif
    285366       !
    286     USE phys_state_var_mod ! Variables sauvegardees de la physique
    287 #ifdef CPP_Dust
    288     USE phys_output_write_spl_mod
    289 #else
    290     USE phys_output_var_mod ! Variables pour les ecritures des sorties
    291 #endif
    292 
    293     USE phys_output_write_mod
    294     USE fonte_neige_mod, ONLY  : fonte_neige_get_vars
    295     USE phys_output_mod
    296     USE phys_output_ctrlout_mod
    297     USE open_climoz_m, only: open_climoz ! ozone climatology from a file
    298     USE regr_pr_time_av_m, only: regr_pr_time_av
    299     USE netcdf95, only: nf95_close
    300     !IM for NMC files
    301     USE netcdf, only: nf90_fill_real
    302     USE mod_phys_lmdz_mpi_data, only: is_mpi_root
    303     USE aero_mod
    304     USE ozonecm_m, only: ozonecm ! ozone of J.-F. Royer
    305     USE conf_phys_m, only: conf_phys
    306     USE radlwsw_m, only: radlwsw
    307     USE phyaqua_mod, only: zenang_an
    308     USE time_phylmdz_mod, only: day_step_phy, annee_ref, day_ref, itau_phy, &
    309          start_time, pdtphys, day_ini
    310     USE tracinca_mod, ONLY: config_inca
    311 #ifdef CPP_XIOS
    312     USE wxios, ONLY: missing_val, missing_val_omp
    313     USE xios, ONLY: xios_get_field_attr, xios_field_is_active
    314 #endif
    315 #ifdef REPROBUS
    316     USE CHEM_REP, ONLY : Init_chem_rep_xjour, &
    317          d_q_rep,d_ql_rep,d_qi_rep,ptrop,ttrop, &
    318          ztrop, gravit,itroprep, Z1,Z2,fac,B
    319 #endif
    320     USE indice_sol_mod
    321     USE phytrac_mod, ONLY : phytrac_init, phytrac
    322     USE carbon_cycle_mod, ONLY : infocfields_init, RCO2_glo, carbon_cycle_rad
    323 
    324 #ifdef CPP_RRTM
    325     USE YOERAD, ONLY : NRADLP
    326     USE YOESW, ONLY : RSUN
    327 #endif
    328     USE ioipsl_getin_p_mod, ONLY : getin_p
    329 
    330 #ifndef CPP_XIOS
    331     USE paramLMDZ_phy_mod
    332 #endif
    333 
    334     USE cmp_seri_mod
    335     USE add_phys_tend_mod, only : add_pbl_tend, add_phys_tend, diag_phys_tend, prt_enerbil, &
    336   &      fl_ebil, fl_cor_ebil
    337 
    338     !IM stations CFMIP
    339     USE CFMIP_point_locations
    340     USE FLOTT_GWD_rando_m, only: FLOTT_GWD_rando
    341     USE ACAMA_GWD_rando_m, only: ACAMA_GWD_rando
    342     USE VERTICAL_LAYERS_MOD, ONLY: aps,bps, ap, bp
    343     USE etat0_limit_unstruct_mod
    344 #ifdef CPP_XIOS
    345     USE xios, ONLY: xios_update_calendar, xios_context_finalize
    346 #endif
    347     USE limit_read_mod, ONLY : init_limit_read
    348     USE regr_horiz_time_climoz_m, ONLY: regr_horiz_time_climoz
    349     USE readaerosol_mod, ONLY : init_aero_fromfile
    350     USE readaerosolstrato_m, ONLY : init_readaerosolstrato
     367
    351368
    352369    IMPLICIT NONE
     
    671688    !$OMP THREADPRIVATE(iflag_alp_wk_cond)
    672689
    673     INTEGER,  SAVE               :: iflag_bug_t2m_ipslcm61=1 !
    674     !$OMP THREADPRIVATE(iflag_bug_t2m_ipslcm61)
    675     INTEGER,  SAVE               :: iflag_bug_t2m_stab_ipslcm61=-1 !
    676     !$OMP THREADPRIVATE(iflag_bug_t2m_stab_ipslcm61)
    677 
    678690    REAL t_w(klon,klev),q_w(klon,klev) ! temperature and moisture profiles in the wake region
    679691    REAL t_x(klon,klev),q_x(klon,klev) ! temperature and moisture profiles in the off-wake region
     
    14191431       WRITE(lunout,*) 'debut physiq_mod tau_gl=',tau_gl
    14201432
    1421        iflag_bug_t2m_ipslcm61 = 1
    1422        CALL getin_p('iflag_bug_t2m_ipslcm61', iflag_bug_t2m_ipslcm61)
    1423        iflag_bug_t2m_stab_ipslcm61 = -1
    1424        CALL getin_p('iflag_bug_t2m_stab_ipslcm61', iflag_bug_t2m_stab_ipslcm61)
    1425 
    14261433       CALL getin_p('iflag_alp_wk_cond', iflag_alp_wk_cond)
    14271434       CALL getin_p('random_notrig_max',random_notrig_max)
     
    14511458       iflag_phytrac = 1 ! by default we do want to call phytrac
    14521459       CALL getin_p('iflag_phytrac',iflag_phytrac)
     1460#ifdef CPP_Dust
     1461       IF (iflag_phytrac.EQ.0) THEN
     1462         WRITE(lunout,*) 'In order to run with SPLA, iflag_phytrac will be forced to 1'
     1463         iflag_phytrac = 1
     1464       ENDIF
     1465#endif
    14531466       nvm_lmdz = 13
    14541467       CALL getin_p('NVM',nvm_lmdz)
     
    15081521       tau_overturning_th(:)=0.
    15091522
    1510        IF (type_trac == 'inca') THEN
     1523       IF (type_trac == 'inca' .OR. type_trac == 'inco') THEN
    15111524          ! jg : initialisation jusqu'au ces variables sont dans restart
    15121525          ccm(:,:,:) = 0.
     
    17631776#ifdef CPP_COSP
    17641777      IF (ok_cosp) THEN
    1765            DO k = 1, klev
    1766              DO i = 1, klon
    1767                phicosp(i,k) = pphi(i,k) + pphis(i)
    1768              ENDDO
    1769            ENDDO
     1778!           DO k = 1, klev
     1779!             DO i = 1, klon
     1780!               phicosp(i,k) = pphi(i,k) + pphis(i)
     1781!             ENDDO
     1782!           ENDDO
    17701783        CALL phys_cosp(itap,phys_tstep,freq_cosp, &
    17711784               ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, &
     
    17851798#ifdef CPP_COSP2
    17861799        IF (ok_cosp) THEN
    1787            DO k = 1, klev
    1788              DO i = 1, klon
    1789                phicosp(i,k) = pphi(i,k) + pphis(i)
    1790              ENDDO
    1791            ENDDO
     1800!           DO k = 1, klev
     1801!             DO i = 1, klon
     1802!               phicosp(i,k) = pphi(i,k) + pphis(i)
     1803!             ENDDO
     1804!           ENDDO
    17921805          CALL phys_cosp2(itap,phys_tstep,freq_cosp, &
    17931806               ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, &
     
    18351848
    18361849       CALL iniradia(klon,klev,paprs(1,1:klev+1))
    1837        ! Initialisation des champs dans phytrac qui sont utilisés par phys_output_write
     1850
     1851       ! Initialisation des champs dans phytrac* qui sont utilisés par phys_output_write*
     1852#ifdef CPP_Dust
     1853       ! Quand on utilise SPLA, on force iflag_phytrac=1
     1854       CALL phytracr_spl_out_init()
     1855       CALL phys_output_write_spl(itap, pdtphys, paprs, pphis,                  &
     1856                                pplay, lmax_th, aerosol_couple,                 &
     1857                                ok_ade, ok_aie, ivap, ok_sync,                  &
     1858                                ptconv, read_climoz, clevSTD,                   &
     1859                                ptconvth, d_t, qx, d_qx, d_tr_dyn, zmasse,      &
     1860                                flag_aerosol, flag_aerosol_strat, ok_cdnc)
     1861#else
     1862       ! phys_output_write écrit des variables traceurs seulement si iflag_phytrac == 1
     1863       ! donc seulement dans ce cas on doit appeler phytrac_init()
    18381864       IF (iflag_phytrac == 1 ) THEN
    18391865          CALL phytrac_init()
    1840         ENDIF
    1841 
     1866       ENDIF
    18421867       CALL phys_output_write(itap, pdtphys, paprs, pphis,                    &
    18431868                              pplay, lmax_th, aerosol_couple,                 &
     
    18461871                              ptconvth, d_u, d_t, qx, d_qx, zmasse,           &
    18471872                              flag_aerosol, flag_aerosol_strat, ok_cdnc)
     1873#endif
     1874
    18481875
    18491876#ifdef CPP_XIOS
     
    20702097       !c         ENDDO
    20712098       !
    2072        IF (type_trac == 'inca') THEN
     2099       IF (type_trac == 'inca' .OR. type_trac == 'inco') THEN                   ! ModThL
    20732100#ifdef INCA
    20742101          CALL VTe(VTphysiq)
     
    20872114               klon, &
    20882115               nqtot, &
    2089                nqo, &
     2116               nqo+nqCO2, &
    20902117               pdtphys, &
    20912118               annee_ref, &
     
    28402867    !   s_therm,   s_trmb1,   s_trmb2, s_trmb3,
    28412868    !   zu10m,     zv10m,   fder,
    2842     !   zxqsurf,   rh2m,      zxfluxu, zxfluxv,
     2869    !   zxqsurf,   delta_qsurf,
     2870    !   rh2m,      zxfluxu, zxfluxv,
    28432871    !   frugs,     agesno,    fsollw,  fsolsw,
    28442872    !   d_ts,      fevap,     fluxlat, t2m,
     
    28902918            debut,     lafin, &
    28912919            longitude_deg, latitude_deg, rugoro,  zrmu0,      &
    2892             zsig,      sollwdown, pphi,    cldt,      &
    2893             rain_fall, snow_fall, solsw,   sollw,     &
     2920            sollwdown,    cldt,      &
     2921            rain_fall, snow_fall, solsw,   solswfdiff, sollw,     &
    28942922            gustiness,                                &
    28952923            t_seri,    q_seri,    u_seri,  v_seri,    &
     
    29012929                                !albedo SB <<<
    29022930            cdragh,    cdragm,  u1,    v1,            &
     2931            beta_aridity, &
    29032932                                !albedo SB >>>
    29042933                                ! albsol1,   albsol2,   sens,    evap,      &
     
    29062935                                !albedo SB <<<
    29072936            albsol3_lic,runoff,   snowhgt,   qsnow, to_ice, sissnow, &
    2908             zxtsol,    zxfluxlat, zt2m,    qsat2m,  &
     2937            zxtsol,    zxfluxlat, zt2m,    qsat2m,  zn2mout, &
    29092938            d_t_vdf,   d_q_vdf,   d_u_vdf, d_v_vdf, d_t_diss, &
    29102939                                !nrlmd<
     
    29272956            s_therm,   s_trmb1,   s_trmb2, s_trmb3, &
    29282957            zustar, zu10m,     zv10m,   fder, &
    2929             zxqsurf,   rh2m,      zxfluxu, zxfluxv, &
     2958            zxqsurf, delta_qsurf,   rh2m,      zxfluxu, zxfluxv, &
    29302959            z0m, z0h,     agesno,    fsollw,  fsolsw, &
    29312960            d_ts,      fevap,     fluxlat, t2m, &
     
    30413070#endif
    30423071       ENDIF
    3043 
    3044 !add limitation for t,q at and wind at 10m
    3045         if ( iflag_bug_t2m_ipslcm61 == 0 ) THEN
    3046           CALL borne_var_surf( klon,klev,nbsrf,                 &
    3047             iflag_bug_t2m_stab_ipslcm61,                        &
    3048             t_seri(:,1),q_seri(:,1),u_seri(:,1),v_seri(:,1),    &
    3049             ftsol,zxqsurf,pctsrf,paprs,                         &
    3050             t2m, q2m, u10m, v10m,                               &
    3051             zt2m_cor, zq2m_cor, zu10m_cor, zv10m_cor,           &
    3052             zrh2m_cor, zqsat2m_cor)
    3053         ELSE
    3054           zt2m_cor(:)=zt2m(:)
    3055           zq2m_cor(:)=zq2m(:)
    3056           zu10m_cor(:)=zu10m(:)
    3057           zv10m_cor(:)=zv10m(:)
    3058           zqsat2m_cor=999.999
    3059         ENDIF
    30603072
    30613073       !---------------------------------------------------------------------
     
    39773989             ENDDO
    39783990          ENDIF
    3979        
     3991         
    39803992#ifdef ISOVERIF
    39813993        write(*,*) 'physiq 3977: verif des inputs de calwake'
     
    43124324          !
    43134325          ENDIF  ! (mod(iflag_pbl_split/10,10) .GE. 1)
    4314 
    43154326          !
    43164327          CALL add_phys_tend(d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs,  &
     
    44564467         iflag_ratqs,iflag_con,iflag_cld_th,pdtphys,  &
    44574468         ratqsbas,ratqshaut,ratqsp0, ratqsdp, &
    4458          tau_ratqs,fact_cldcon,   &
     4469         tau_ratqs,fact_cldcon,wake_s, wake_deltaq,   &
    44594470         ptconv,ptconvth,clwcon0th, rnebcon0th,     &
    44604471         paprs,pplay,q_seri,zqsat,fm_therm, &
    4461          ratqs,ratqsc)
    4462 
     4472         ratqs,ratqsc,ratqs_inter)
    44634473
    44644474    !
     
    50135023          ENDIF
    50145024          zx_rh(i,k) = q_seri(i,k)/zx_qs
     5025            IF (iflag_ice_thermo .GT. 0) THEN
     5026          zx_rhl(i,k) = q_seri(i,k)/(qsatl(zx_t)/pplay(i,k))
     5027          zx_rhi(i,k) = q_seri(i,k)/(qsats(zx_t)/pplay(i,k))
     5028            ENDIF
    50155029          zqsat(i,k)=zx_qs
    50165030       ENDDO
     
    50395053    ENDDO
    50405054
    5041     IF (type_trac == 'inca') THEN
     5055    IF (type_trac == 'inca' .OR. type_trac == 'inco') THEN      ! ModThL
    50425056#ifdef INCA
    50435057       CALL VTe(VTphysiq)
     
    50825096            nbp_lon, &
    50835097            nbp_lat-1, &
    5084             tr_seri, &
     5098            tr_seri(:,:,1+nqCO2:nbtr), &
    50855099            ftsol, &
    50865100            paprs, &
     
    50935107       CALL VTe(VTinca)
    50945108       CALL VTb(VTphysiq)
    5095 #endif 
    5096     ENDIF !type_trac = inca
     5109#endif
     5110    ENDIF !type_trac = inca or inco
    50975111    IF (type_trac == 'repr') THEN
    50985112#ifdef REPROBUS
     
    52655279
    52665280       IF (ok_newmicro) then
    5267           IF (iflag_rrtm.NE.0) THEN
     5281! AI          IF (iflag_rrtm.NE.0) THEN
     5282          IF (iflag_rrtm.EQ.1) THEN
    52685283#ifdef CPP_RRTM
    52695284             IF (ok_cdnc.AND.NRADLP.NE.3) THEN
     
    54525467               heat,heat0,cool,cool0,albpla, &
    54535468               heat_volc,cool_volc, &
    5454                topsw,toplw,solsw,sollw, &
     5469               topsw,toplw,solsw,solswfdiff,sollw, &
    54555470               sollwdown, &
    54565471               topsw0,toplw0,solsw0,sollw0, &
     
    55395554                     heatp,heat0p,coolp,cool0p,albplap, &
    55405555                     heat_volc,cool_volc, &
    5541                      topswp,toplwp,solswp,sollwp, &
     5556                     topswp,toplwp,solswp,solswfdiffp,sollwp, &
    55425557                     sollwdownp, &
    55435558                     topsw0p,toplw0p,solsw0p,sollw0p, &
     
    59765991
    59775992    CALL tend_to_tke(pdtphys,paprs,exner,t_seri,u_seri,v_seri,dtadd,duadd,dvadd,pctsrf,pbl_tke)
    5978 
     5993   !
     5994   ! Prevent pbl_tke_w from becoming negative
     5995    wake_delta_pbl_tke(:,:,:) = max(wake_delta_pbl_tke(:,:,:), -pbl_tke(:,:,:))
     5996   !
    59795997
    59805998       ENDIF
     
    60906108#ifdef CPP_COSPV2
    60916109       IF (itap.eq.1.or.MOD(itap,NINT(freq_cosp/phys_tstep)).EQ.0) THEN
     6110!        IF (MOD(itap,NINT(freq_cosp/phys_tstep)).EQ.0) THEN
    60926111
    60936112          IF (prt_level .GE.10) THEN
    60946113             print*,'freq_cosp',freq_cosp
    60956114          ENDIF
     6115           DO k = 1, klev
     6116             DO i = 1, klon
     6117               phicosp(i,k) = pphi(i,k) + pphis(i)
     6118             ENDDO
     6119           ENDDO
    60966120          mr_ozone=wo(:, :, 1) * dobson_u * 1e3 / zmasse
    60976121                 print*,'Dans physiq.F avant appel '
     
    61946218    ELSE
    61956219       sh_in(:,:) = qx(:,:,ivap)
    6196        ch_in(:,:) = qx(:,:,iliq)
     6220       IF (nqo .EQ. 3) THEN
     6221          ch_in(:,:) = qx(:,:,iliq) + qx(:,:,isol)
     6222       ELSE
     6223          ch_in(:,:) = qx(:,:,iliq)
     6224       ENDIF
    61976225    ENDIF
    61986226
    6199     IF (iflag_phytrac == 1 ) THEN
    6200 
    62016227#ifdef CPP_Dust
    6202       CALL       phytracr_spl ( debut,lafin , jD_cur,jH_cur,iflag_con,       &  ! I
     6228    !  Avec SPLA, iflag_phytrac est forcé =1
     6229    CALL       phytracr_spl ( debut,lafin , jD_cur,jH_cur,iflag_con,       &  ! I
    62036230                      pdtphys,ftsol,                                   &  ! I
    62046231                      t,q_seri,paprs,pplay,RHcl,                  &  ! I
     
    62166243
    62176244#else
    6218 
    6219     CALL phytrac ( &
     6245    IF (iflag_phytrac == 1 ) THEN
     6246      CALL phytrac ( &
    62206247         itap,     days_elapsed+1,    jH_cur,   debut, &
    62216248         lafin,    phys_tstep,     u, v,     t, &
     
    62546281
    62556282#endif
    6256 
    6257 #endif
    62586283    ENDIF    ! (iflag_phytrac=1)
     6284
     6285#endif
     6286    !ENDIF    ! (iflag_phytrac=1)
    62596287
    62606288    IF (offline) THEN
     
    63456373#endif
    63466374    !
    6347     IF (type_trac == 'inca') THEN
     6375    IF (type_trac == 'inca' .OR. type_trac == 'inco') THEN
    63486376#ifdef INCA
    63496377       CALL VTe(VTphysiq)
     
    63546382            pplay, &
    63556383            t_seri, &
    6356             tr_seri, &
     6384            tr_seri(:,:,1+nqCO2:nbtr), &
    63576385            nbtr, &
    63586386            paprs, &
     
    66206648  CALL phys_output_write_spl(itap, pdtphys, paprs, pphis,  &
    66216649       pplay, lmax_th, aerosol_couple,                 &
    6622        ok_ade, ok_aie, ivap, ok_sync,         &
     6650       ok_ade, ok_aie, ivap, ok_sync,                  &
    66236651       ptconv, read_climoz, clevSTD,                   &
    66246652       ptconvth, d_t, qx, d_qx, d_tr_dyn, zmasse,      &
     
    66396667#endif
    66406668
    6641 ! On remet des variables a .false. apres un premier appel
     6669! Pour XIOS : On remet des variables a .false. apres un premier appel
    66426670    IF (debut) THEN
    66436671#ifdef CPP_XIOS
  • LMDZ6/trunk/libf/phylmdiso/surf_landice_mod.F90

    r3927 r3940  
    1919       tsoil, z0m, z0h, SFRWL, alb_dir, alb_dif, evap, fluxsens, fluxlat, &
    2020       tsurf_new, dflux_s, dflux_l, &
    21        slope, cloudf, &
     21       alt, slope, cloudf, &
    2222       snowhgt, qsnow, to_ice, sissnow, &
    2323       alb3, runoff, &
     
    3030
    3131    USE dimphy
    32     USE surface_data,     ONLY : type_ocean, calice, calsno, ok_snow
    33     USE fonte_neige_mod,  ONLY : fonte_neige, run_off_lic
     32    USE surface_data,     ONLY : type_ocean, calice, calsno, landice_opt, iflag_albcalc
     33    USE fonte_neige_mod,  ONLY : fonte_neige,run_off_lic,fqcalving_global,ffonte_global,fqfonte_global,runofflic_global
    3434    USE cpl_mod,          ONLY : cpl_send_landice_fields
    3535    USE calcul_fluxs_mod
     
    4747    USE ioipsl_getin_p_mod, ONLY : getin_p
    4848
    49 #ifdef CPP_SISVAT
    50     USE surf_sisvat_mod,  ONLY : surf_sisvat
    51 #endif
     49
     50#ifdef CPP_INLANDSIS
     51    USE surf_inlandsis_mod,  ONLY : surf_inlandsis
     52#endif
     53
    5254    USE indice_sol_mod
    5355
     
    8890    REAL, DIMENSION(klon), INTENT(IN)             :: albedo  !mean albedo
    8991    REAL, DIMENSION(klon), INTENT(IN)             :: pphi1   
     92    REAL, DIMENSION(klon), INTENT(IN)             :: alt   !mean altitude of the grid box 
    9093    REAL, DIMENSION(klon), INTENT(IN)             :: slope   !mean slope in grid box 
    9194    REAL, DIMENSION(klon), INTENT(IN)             :: cloudf  !total cloud fraction
     
    108111!    REAL, DIMENSION(klon), INTENT(OUT)            :: alb1  ! new albedo in visible SW interval
    109112!    REAL, DIMENSION(klon), INTENT(OUT)            :: alb2  ! new albedo in near IR interval
    110     REAL, DIMENSION(6), INTENT(IN)              ::SFRWL
    111     REAL, DIMENSION(klon,nsw), INTENT(OUT)        ::alb_dir,alb_dif
     113    REAL, DIMENSION(6), INTENT(IN)                :: SFRWL
     114    REAL, DIMENSION(klon,nsw), INTENT(OUT)        :: alb_dir,alb_dif
    112115!albedo SB <<<
    113116    REAL, DIMENSION(klon), INTENT(OUT)            :: evap, fluxsens, fluxlat
     
    135138    REAL, DIMENSION(klon)    :: zfra, alb_neig
    136139    REAL, DIMENSION(klon)    :: radsol
    137     REAL, DIMENSION(klon)    :: u0, v0, u1_lay, v1_lay
    138     INTEGER                  :: i,j
     140    REAL, DIMENSION(klon)    :: u0, v0, u1_lay, v1_lay, ustar
     141    INTEGER                  :: i,j,nt
     142    REAL, DIMENSION(klon)    :: fqfonte,ffonte
    139143#ifdef ISO       
    140144      real, parameter :: t_coup = 273.15
     
    155159    REAL, DIMENSION(klon)    :: emis_new                  !Emissivity
    156160    REAL, DIMENSION(klon)    :: swdown,lwdown
    157     REAL, DIMENSION(klon)    :: precip_snow_adv, snow_adv !Snow Drift precip./advection
    158     REAL, DIMENSION(klon)    :: bl_height, wind_velo      !height boundary layer, wind spd
     161    REAL, DIMENSION(klon)    :: precip_snow_adv, snow_adv !Snow Drift precip./advection (not used in inlandsis)
     162    REAL, DIMENSION(klon)    :: erod                      !erosion of surface snow (flux, kg/m2/s like evap)
     163    REAL, DIMENSION(klon)    :: zsl_height, wind_velo     !surface layer height, wind spd
    159164    REAL, DIMENSION(klon)    :: dens_air,  snow_cont_air  !air density; snow content air
    160165    REAL, DIMENSION(klon)    :: alb_soil                  !albedo of underlying ice
    161166    REAL, DIMENSION(klon)    :: pexner                    !Exner potential
    162167    REAL                     :: pref
    163     REAL, DIMENSION(klon,nsoilmx) :: tsoil0 !modfi
     168    REAL, DIMENSION(klon,nsoilmx) :: tsoil0               !modif
     169    REAL                          :: dtis                ! subtimestep
     170    LOGICAL                       :: debut_is, lafin_is  ! debut and lafin for inlandsis
    164171
    165172    CHARACTER (len = 20)                      :: modname = 'surf_landice'
     
    167174
    168175
    169 !albedo SB >>>
    170     real,dimension(klon) :: alb1,alb2
    171 !albedo SB <<<
    172 
     176    REAL,DIMENSION(klon) :: alb1,alb2
     177    REAL, DIMENSION (klon,6) :: alb6
    173178! End definition
    174179!****************************************************************************************
     
    181186  LOGICAL, SAVE :: firstcall = .TRUE.
    182187  !$OMP THREADPRIVATE(firstcall)
    183 !FC
    184 
    185 
     188
     189
     190!FC firtscall initializations
     191!******************************************************************************************
    186192#ifdef ISO
    187193#ifdef ISOVERIF
     
    203209  CALL getin_p('alb_nir_sno_lic',alb_nir_sno_lic)
    204210           PRINT*, 'alb_nir_sno_lic',alb_nir_sno_lic
     211 
     212!  z0m=1.e-3
     213!  z0h = z0m
    205214  firstcall=.false.
    206215  ENDIF
     216!******************************************************************************************
    207217!
    208218! Initialize output variables
     
    220230
    221231!****************************************************************************************
    222 !   ok_snow = TRUE  : prepare and call SISVAT snow model
    223 !   ok_snow = FALSE : soil_model, calcul_flux, fonte_neige, ...
    224 !
    225 !****************************************************************************************
    226     IF (ok_snow) THEN
    227 #ifdef CPP_SISVAT
    228        ! Prepare for calling SISVAT
     232!  landice_opt = 0 : soil_model, calcul_flux, fonte_neige, ... 
     233!  landice_opt = 1  : prepare and call INterace Lmdz SISvat (INLANDSIS)
     234!****************************************************************************************
     235
     236
     237    IF (landice_opt .EQ. 1) THEN
     238
     239!****************************************************************************************   
     240! CALL to INLANDSIS interface
     241!****************************************************************************************
     242#ifdef CPP_INLANDSIS
     243
     244#ifdef ISO
     245        CALL abort_gcm('surf_landice 235','isotopes pas dans INLANDSIS',1)
     246#endif
     247
     248        debut_is=debut
     249        lafin_is=.false.
     250        ! Suppose zero surface speed
     251        u0(:)            = 0.0
     252        v0(:)            = 0.0
     253
     254
     255        CALL calcul_flux_wind(knon, dtime, &
     256         u0, v0, u1, v1, gustiness, cdragm, &
     257         AcoefU, AcoefV, BcoefU, BcoefV, &
     258         p1lay, temp_air, &
     259         flux_u1, flux_v1)
     260
    229261       
    230        ! Calculate incoming flux for SW and LW interval: swdown, lwdown
     262       ! Set constants and compute some input for SISVAT
     263       ! = 1000 hPa
     264       ! and calculate incoming flux for SW and LW interval: swdown, lwdown
    231265       swdown(:)        = 0.0
    232266       lwdown(:)        = 0.0
     267       snow_cont_air(:) = 0.  ! the snow content in air is not a prognostic variable of the model     
     268       alb_soil(:)      = 0.4 ! before albedo(:) but here it is the ice albedo that we have to set
     269       ustar(:)         = 0.
     270       pref             = 100000.       
    233271       DO i = 1, knon
    234272          swdown(i)        = swnet(i)/(1-albedo(i))
    235273          lwdown(i)        = lwdownm(i)
    236        END DO
    237        
    238        ! Set constants and compute some input for SISVAT
    239        snow_adv(:)      = 0.                          ! no snow blown in for now
    240        snow_cont_air(:) = 0.       
    241        alb_soil(:)      = albedo(:)
    242        pref             = 100000.                     ! = 1000 hPa
    243        DO i = 1, knon
    244274          wind_velo(i)     = u1(i)**2 + v1(i)**2
    245275          wind_velo(i)     = wind_velo(i)**0.5
    246276          pexner(i)        = (p1lay(i)/pref)**(RD/RCPD)
    247277          dens_air(i)      = p1lay(i)/RD/temp_air(i)  ! dry air density
    248           bl_height(i)     = pphi1(i)/RG             
     278          zsl_height(i)    = pphi1(i)/RG     
     279          tsoil0(i,:)      = tsoil(i,:) 
     280          ustar(i)= (cdragm(i)*(wind_velo(i)**2))**0.5   
    249281       END DO
    250 
    251 !****************************************************************************************
    252 ! CALL to SISVAT interface
    253 !
    254 !****************************************************************************************
    255        ! config: compute everything with SV but temperatures afterwards with soil/calculfluxs
    256        DO i = 1, knon
    257           tsoil0(i,:)=tsoil(i,:)
    258        END DO
    259            ! Martin
    260            PRINT*, 'on appelle surf_sisvat'
    261            ! Martin
    262        CALL surf_sisvat(knon, rlon, rlat, knindex, itime, dtime, debut, lafin, &
    263             rmu0, swdown, lwdown, pexner, ps, p1lay, &
    264             precip_rain, precip_snow, precip_snow_adv, snow_adv, &
    265             bl_height, wind_velo, temp_air, dens_air, spechum, tsurf, &
    266             rugoro, snow_cont_air, alb_soil, slope, cloudf, &
    267             radsol, qsol, tsoil0, snow, snowhgt, qsnow, to_ice,sissnow, agesno, &
    268             AcoefH, AcoefQ, BcoefH, BcoefQ, cdragh, &
    269             run_off_lic, evap, fluxsens, fluxlat, dflux_s, dflux_l, &       
    270             tsurf_new, alb1, alb2, alb3, &
    271             emis_new, z0m, qsurf)
    272        z0h(1:knon)=z0m(1:knon) ! en attendant mieux
    273282       
    274        ! Suppose zero surface speed
    275        u0(:)            = 0.0
    276        v0(:)            = 0.0
    277        ! The calculation of heat/water fluxes, otherwise done by "CALL calcul_fluxs" is
    278        ! integrated in SISVAT, using the same method. It can be found in "sisvat.f", in the
    279        ! subroutine "SISVAT_TS2".
    280        ! u0, v0=0., dif_grnd=0. and beta=1 are assumed there!
    281        
    282        CALL calcul_flux_wind(knon, dtime, &
    283             u0, v0, u1, v1, gustiness, cdragm, &
    284             AcoefU, AcoefV, BcoefU, BcoefV, &
    285             p1lay, temp_air, &
    286             flux_u1, flux_v1)
     283
     284
     285        dtis=dtime
     286
     287          IF (lafin) THEN
     288            lafin_is=.true.
     289          END IF
     290
     291          CALL surf_inlandsis(knon, rlon, rlat, knindex, itime, dtis, debut_is, lafin_is,&
     292            rmu0, swdown, lwdown, albedo, pexner, ps, p1lay, precip_rain, precip_snow,   &
     293            zsl_height, wind_velo, ustar, temp_air, dens_air, spechum, tsurf,&
     294            rugoro, snow_cont_air, alb_soil, alt, slope, cloudf, &
     295            radsol, qsol, tsoil0, snow, zfra, snowhgt, qsnow, to_ice, sissnow,agesno,   &
     296            AcoefH, AcoefQ, BcoefH, BcoefQ, cdragm, cdragh, &
     297            run_off_lic, fqfonte, ffonte, evap, erod, fluxsens, fluxlat,dflux_s, dflux_l, &
     298            tsurf_new, alb1, alb2, alb3, alb6, &
     299            emis_new, z0m, z0h, qsurf)
     300
     301          debut_is=.false.
     302
     303
     304        ! Treatment of snow melting and calving
     305
     306        ! for consistency with standard LMDZ, add calving to run_off_lic
     307        run_off_lic(:)=run_off_lic(:) + to_ice(:)
     308
     309        DO i = 1, knon
     310           ffonte_global(knindex(i),is_lic)    = ffonte(i)
     311           fqfonte_global(knindex(i),is_lic)   = fqfonte(i)! net melting= melting - refreezing
     312           fqcalving_global(knindex(i),is_lic) = to_ice(i) ! flux
     313           runofflic_global(knindex(i)) = run_off_lic(i)
     314        ENDDO
     315        ! Here, we assume that the calving term is equal to the to_ice term
     316        ! (no ice accumulation)
     317
     318
    287319#else
    288        abort_message='Pb de coherence: ok_snow = .true. mais CPP_SISVAT = .false.'
     320       abort_message='Pb de coherence: landice_opt = 1  mais CPP_INLANDSIS = .false.'
    289321       CALL abort_physic(modname,abort_message,1)
    290322#endif
    291 #ifdef ISO
    292        abort_message='surf_landice 267: isotopes pas prevus ici'
    293        CALL abort_physic(modname,abort_message,1)
    294 #endif
    295     ELSE ! ok_snow=FALSE
     323
     324
     325
     326    ELSE
    296327
    297328!****************************************************************************************
     
    299330!
    300331!****************************************************************************************
     332
     333    ! EV: use calbeta
     334    CALL calbeta(dtime, is_lic, knon, snow, qsol, beta, cal, dif_grnd)
     335
     336
     337    ! use soil model and recalculate properly cal
    301338    IF (soil_model) THEN
    302339       CALL soil(dtime, is_lic, knon, snow, tsurf, tsoil, soilcap, soilflux)
     
    313350!
    314351!****************************************************************************************
    315     beta(:) = 1.0
    316     dif_grnd(:) = 0.0
     352!    beta(:) = 1.0
     353!    dif_grnd(:) = 0.0
    317354
    318355! Suppose zero surface speed
     
    367404!   
    368405!****************************************************************************************
    369     CALL fonte_neige( knon, is_lic, knindex, dtime, &
     406    CALL fonte_neige(knon, is_lic, knindex, dtime, &
    370407         tsurf, precip_rain, precip_snow, &
    371408         snow, qsol, tsurf_new, evap &
     
    408445!****************************************************************************************
    409446    CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) 
     447
    410448    WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
    411449    zfra(1:knon) = MAX(0.0,MIN(1.0,snow(1:knon)/(snow(1:knon)+10.0)))
     
    420458!IM: KstaTER0.77 & LMD_ARMIP6   
    421459
    422 ! Attantion: alb1 and alb2 are the same!
     460! Attantion: alb1 and alb2 are not the same!
    423461    alb1(1:knon)  = alb_vis_sno_lic
    424462    alb2(1:knon)  = alb_nir_sno_lic
     
    433471    z0m = SQRT(z0m**2+rugoro**2)
    434472
    435     END IF ! ok_snow
     473
     474
     475    END IF ! landice_opt
    436476
    437477
    438478!****************************************************************************************
    439479! Send run-off on land-ice to coupler if coupled ocean.
    440 ! run_off_lic has been calculated in fonte_neige or surf_sisvat
     480! run_off_lic has been calculated in fonte_neige or surf_inlandsis
    441481!
    442482!****************************************************************************************
     
    450490 
    451491!****************************************************************************************
    452        snow_o=0.
    453        zfra_o = 0.
    454        DO j = 1, knon
    455            i = knindex(j)
    456            snow_o(i) = snow(j)
    457            zfra_o(i) = zfra(j)
    458        ENDDO
    459 
     492! Etienne: comment these lines because of duplication just below
     493!       snow_o=0.
     494!       zfra_o = 0.
     495!       DO j = 1, knon
     496!           i = knindex(j)
     497!           snow_o(i) = snow(j)
     498!           zfra_o(i) = zfra(j)
     499!       ENDDO
     500!
    460501!****************************************************************************************
    461502       snow_o=0.
     
    485526       alb_dir(1:knon,5)=alb2(1:knon)
    486527       alb_dir(1:knon,6)=alb2(1:knon)
     528
     529       IF ((landice_opt .EQ. 1) .AND. (iflag_albcalc .EQ. 2)) THEN
     530       alb_dir(1:knon,1)=alb6(1:knon,1)
     531       alb_dir(1:knon,2)=alb6(1:knon,2)
     532       alb_dir(1:knon,3)=alb6(1:knon,3)
     533       alb_dir(1:knon,4)=alb6(1:knon,4)
     534       alb_dir(1:knon,5)=alb6(1:knon,5)
     535       alb_dir(1:knon,6)=alb6(1:knon,6)
     536       ENDIF
     537
    487538     end select
    488539alb_dif=alb_dir
    489540!albedo SB <<<
    490541
    491 
    492 
     542 
     543 
    493544
    494545  END SUBROUTINE surf_landice
  • LMDZ6/trunk/libf/phylmdiso/surf_ocean_mod.F90

    r3927 r3940  
    2020       z0m, z0h, SFRWL, alb_dir_new, alb_dif_new, evap, fluxsens, fluxlat, &
    2121       tsurf_new, dflux_s, dflux_l, lmt_bils, &
    22        flux_u1, flux_v1 &
     22       flux_u1, flux_v1, delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, tkt, tks, &
     23       taur, sss &
    2324#ifdef ISO
    2425        &       ,xtprecip_rain, xtprecip_snow,xtspechum,Roce, &
     
    2829
    2930    use albedo, only: alboc, alboc_cd
     31    use bulk_flux_m, only: bulk_flux
    3032    USE dimphy, ONLY: klon, zmasq
    3133    USE surface_data, ONLY     : type_ocean
     
    4244#endif
    4345    USE limit_read_mod
     46    use config_ocean_skin_m, only: activate_ocean_skin
    4447    !
    4548    ! This subroutine will make a call to ocean_XXX_noice according to the ocean mode (force,
     
    6265    REAL, DIMENSION(klon), INTENT(IN)        :: lwnet  ! net longwave radiation at surface 
    6366    REAL, DIMENSION(klon), INTENT(IN)        :: alb1   ! albedo in visible SW interval
    64     REAL, DIMENSION(klon), INTENT(IN)        :: windsp
     67    REAL, DIMENSION(klon), INTENT(IN)        :: windsp ! wind at 10 m, in m s-1
    6568    REAL, DIMENSION(klon), INTENT(IN)        :: rmu0 
    6669    REAL, DIMENSION(klon), INTENT(IN)        :: fder
    67     REAL, DIMENSION(klon), INTENT(IN)        :: tsurf_in
     70    REAL, DIMENSION(klon), INTENT(IN)        :: tsurf_in     ! defined only for subscripts 1:knon
    6871    REAL, DIMENSION(klon), INTENT(IN)        :: p1lay,z1lay ! pression (Pa) et altitude (m) du premier niveau
    6972    REAL, DIMENSION(klon), INTENT(IN)        :: cdragh
     
    8790    REAL, DIMENSION(klon), INTENT(INOUT)     :: qsurf
    8891    REAL, DIMENSION(klon), INTENT(INOUT)     :: agesno
     92    REAL, DIMENSION(klon), INTENT(inOUT)     :: z0h
    8993#ifdef ISO
    9094    REAL, DIMENSION(niso,klon), INTENT(IN)        :: xtsnow
    9195    REAL, DIMENSION(niso,klon), INTENT(INOUT)        :: Roce 
    9296#endif
    93     REAL, DIMENSION(klon), INTENT(inOUT):: z0h
     97
     98    REAL, intent(inout):: delta_sst(:) ! (knon)
     99    ! Ocean-air interface temperature minus bulk SST, in K. Defined
     100    ! only if activate_ocean_skin >= 1.
     101
     102    real, intent(inout):: delta_sal(:) ! (knon)
     103    ! Ocean-air interface salinity minus bulk salinity, in ppt. Defined
     104    ! only if activate_ocean_skin >= 1.
     105
     106    REAL, intent(inout):: ds_ns(:) ! (knon)
     107    ! "delta salinity near surface". Salinity variation in the
     108    ! near-surface turbulent layer. That is subskin salinity minus
     109    ! foundation salinity. In ppt.
     110
     111    REAL, intent(inout):: dt_ns(:) ! (knon)
     112    ! "delta temperature near surface". Temperature variation in the
     113    ! near-surface turbulent layer. That is subskin temperature
     114    ! minus foundation temperature. (Can be negative.) In K.
    94115
    95116    ! Output variables
    96     !******************************************************************************
     117    !**************************************************************************
    97118    REAL, DIMENSION(klon), INTENT(OUT)       :: z0m
    98119    !albedo SB >>>
    99     !    REAL, DIMENSION(klon), INTENT(OUT)       :: alb1_new  ! new albedo in visible SW interval
    100     !    REAL, DIMENSION(klon), INTENT(OUT)       :: alb2_new  ! new albedo in near IR interval
    101     REAL, DIMENSION(6), INTENT(IN)          :: SFRWL
    102     REAL, DIMENSION(klon,nsw), INTENT(OUT)       :: alb_dir_new,alb_dif_new
     120    !    REAL, DIMENSION(klon), INTENT(OUT)  :: alb1_new  ! new albedo in visible SW interval
     121    !    REAL, DIMENSION(klon), INTENT(OUT)  :: alb2_new  ! new albedo in near IR interval
     122    REAL, DIMENSION(6), INTENT(IN)           :: SFRWL
     123    REAL, DIMENSION(klon,nsw), INTENT(OUT)   :: alb_dir_new,alb_dif_new
    103124    !albedo SB <<<     
    104125    REAL, DIMENSION(klon), INTENT(OUT)       :: evap, fluxsens, fluxlat
    105     REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
     126    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new    ! sea surface temperature, in K
    106127    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
    107128    REAL, DIMENSION(klon), INTENT(OUT)       :: lmt_bils
    108129    REAL, DIMENSION(klon), INTENT(OUT)       :: flux_u1, flux_v1
     130
     131    REAL, intent(out):: dter(:) ! (knon)
     132    ! Temperature variation in the diffusive microlayer, that is
     133    ! ocean-air interface temperature minus subskin temperature. In
     134    ! K.
     135
     136    REAL, intent(out):: dser(:) ! (knon)
     137    ! Salinity variation in the diffusive microlayer, that is
     138    ! ocean-air interface salinity minus subskin salinity. In ppt.
     139
     140    REAL, intent(out):: tkt(:) ! (knon)
     141    ! épaisseur (m) de la couche de diffusion thermique (microlayer)
     142    ! cool skin thickness
     143
     144    REAL, intent(out):: tks(:) ! (knon)
     145    ! épaisseur (m) de la couche de diffusion de masse (microlayer)
     146
     147    REAL, intent(out):: taur(:) ! (knon)
     148    ! momentum flux due to rain, in Pa
     149
     150    real, intent(out):: sss(:) ! (klon)
     151    ! Bulk salinity of the surface layer of the ocean, in ppt. (Only
     152    ! defined for subscripts 1:knon, but we have to declare it with
     153    ! size klon because of the coupling machinery.)
    109154#ifdef ISO
    110155    REAL, DIMENSION(ntraciso,klon), INTENT(out)        :: xtevap ! isotopes in surface evaporation flux
     
    113158
    114159    ! Local variables
    115     !******************************************************************************
     160    !*************************************************************************
    116161    INTEGER               :: i, k
    117162    REAL                  :: tmp
     
    121166    REAL, DIMENSION(klon) :: cdragq ! Cdrag pour l'evaporation
    122167    CHARACTER(len=20),PARAMETER :: modname="surf_ocean"
    123 
    124     ! End definition
     168    real rhoa(knon) ! density of moist air  (kg / m3)
     169    REAL sens_prec_liq(knon)
     170
     171    REAL t_int(knon) ! ocean-air interface temperature, in K
     172    real s_int(knon) ! ocean-air interface salinity, in ppt
     173
    125174    !******************************************************************************
    126175
     
    165214
    166215
     216    rhoa = PS(:KNON) / (Rd * temp_air(:knon) * (1. + retv * spechum(:knon)))
     217
    167218    !******************************************************************************
    168219    ! Switch according to type of ocean (couple, slab or forced)
     
    177228            AcoefH, AcoefQ, BcoefH, BcoefQ, &
    178229            AcoefU, AcoefV, BcoefU, BcoefV, &
    179             ps, u1, v1, gustiness, &
     230            ps, u1, v1, gustiness, tsurf_in, &
    180231            radsol, snow, agesno, &
    181232            qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    182             tsurf_new, dflux_s, dflux_l)
     233            tsurf_new, dflux_s, dflux_l, sens_prec_liq, sss, delta_sal, rhoa, &
     234            delta_sst)
    183235
    184236    CASE('slab')
     
    200252            AcoefH, AcoefQ, BcoefH, BcoefQ, &
    201253            AcoefU, AcoefV, BcoefU, BcoefV, &
    202             ps, u1, v1, gustiness, &
     254            ps, u1, v1, gustiness, tsurf_in, &
    203255            radsol, snow, agesno, &
    204256            qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    205             tsurf_new, dflux_s, dflux_l &
     257            tsurf_new, dflux_s, dflux_l, sens_prec_liq, rhoa &
    206258#ifdef ISO
    207259            ,xtprecip_rain, xtprecip_snow, xtspechum,Roce,rlat, &
     
    311363       CALL abort_physic(modname,'version non prevue',1)
    312364    ENDIF
    313     !
    314     !******************************************************************************
     365
     366    if (activate_ocean_skin >= 1) then
     367       if (type_ocean /= 'couple') sss(:knon) = 35.
     368       call bulk_flux(tkt, tks, taur, dter, dser, t_int, s_int, ds_ns, dt_ns, &
     369            u = windsp(:knon), t_ocean_1 = tsurf_new(:knon), s1 = sss(:knon), &
     370            rain = precip_rain(:knon) + precip_snow(:knon), &
     371            hf = - fluxsens(:knon), hlb = - fluxlat(:knon), &
     372            rnl = - lwnet(:knon), &
     373            tau = sqrt(flux_u1(:knon)**2 + flux_v1(:knon)**2), rhoa = rhoa, &
     374            xlv = [(rlvtt, i = 1, knon)], rf = - sens_prec_liq, dtime = dtime, &
     375            rns = swnet(:knon))
     376       delta_sst = t_int - tsurf_new(:knon)
     377       delta_sal = s_int - sss(:knon)
     378       if (activate_ocean_skin >= 2) tsurf_new(:knon) = t_int
     379    end if
     380
    315381  END SUBROUTINE surf_ocean
    316   !******************************************************************************
     382  !****************************************************************************
    317383  !
    318384END MODULE surf_ocean_mod
  • LMDZ6/trunk/libf/phylmdiso/surf_seaice_mod.F90

    r3927 r3940  
    4545    INCLUDE "dimsoil.h"
    4646    INCLUDE "clesphys.h"
     47
     48    INCLUDE "YOMCST.h"
     49    ! for rd and retv
    4750
    4851! Input arguments
     
    113116    REAL, DIMENSION(klon) :: alb1_new,alb2_new
    114117!albedo SB <<<
    115 !
     118
     119    real rhoa(knon) ! density of moist air  (kg / m3)
     120
    116121! End definitions
    117122!****************************************************************************************
     
    124129    radsol(:) = 0.0
    125130    radsol(1:knon) = swnet(1:knon) + lwnet(1:knon)
     131
     132    rhoa = PS(:KNON) / (Rd * temp_air(:knon) * (1. + retv * spechum(:knon)))
    126133
    127134!****************************************************************************************
     
    142149            radsol, snow, qsurf, &
    143150            alb1_new, alb2_new, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    144             tsurf_new, dflux_s, dflux_l)
     151            tsurf_new, dflux_s, dflux_l, rhoa)
    145152       
    146153    ELSE IF (type_ocean == 'slab'.AND.version_ocean=='sicINT') THEN
     
    164171            radsol, snow, qsol, agesno, tsoil, &
    165172            qsurf, alb1_new, alb2_new, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    166             tsurf_new, dflux_s, dflux_l &
     173            tsurf_new, dflux_s, dflux_l, rhoa &
    167174#ifdef ISO
    168175            ,xtprecip_rain, xtprecip_snow, xtspechum,Roce, &
  • LMDZ6/trunk/libf/phylmdiso/thermcell_main.F90

    r3927 r3940  
    267267      REAL xtzo_tmp(klon,klev)
    268268      integer ixt
    269       real xtzl(ntraciso,klon,klev) ! Rq: n'est pas utilisé? Juste diagnostiqué dans thermcell_env. A supprimer?
    270269#endif
    271270!
     
    329328!
    330329      CALL thermcell_env(ngrid,nlay,po,pt,pu,pv,pplay,  &
    331      &           pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out &
    332 #ifdef ISO
    333      &            ,xtpo,xtzo,xtzl &
    334 #endif     
    335      &   )
     330     &           pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out)
    336331       
    337332      if (prt_level.ge.1) print*,'thermcell_main apres thermcell_env'
Note: See TracChangeset for help on using the changeset viewer.