Changeset 3780


Ignore:
Timestamp:
Oct 22, 2020, 2:50:18 PM (4 years ago)
Author:
evignon
Message:

Premiere comission Etienne: changements pour le 1D (forcage en Ts au dessus des continents) et inclusion drag arbres dans yamada4_num=6

Location:
LMDZ6/trunk/libf/phylmd
Files:
23 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/calbeta.F90

    r3102 r3780  
    77  USE dimphy
    88  USE indice_sol_mod
     9
    910  IMPLICIT none
     11
     12#include "flux_arp.h"
     13
    1014!======================================================================
    1115! Auteur(s): Z.X. Li (LMD/CNRS) (adaptation du GCM au LMD)
    … …  
    8286     ENDDO
    8387  ENDIF
     88
     89  ! EV: when beta is prescribed for 1D cases:
     90  IF (knon.EQ.1 .AND. ok_prescr_beta) THEN
     91     DO i = 1, knon
     92          vbeta(i)=betaevap
     93      ENDDO
     94  ENDIF
    8495 
    8596END SUBROUTINE calbeta

  • LMDZ6/trunk/libf/phylmd/change_srf_frac_mod.F90

    r2656 r3780  
    183183           tsurf, alb_dir,alb_dif, ustar, u10m, v10m, pbl_tke)
    184184
    185 
    186185    ELSE
    187186       ! No modifcation should be done

  • LMDZ6/trunk/libf/phylmd/dyn1d/1DUTILS.h

    r3686 r3780  
    233233       CALL getin('ok_flux_surf',ok_flux_surf)
    234234
     235!Config  Key  = ok_forc_tsurf
     236!Config  Desc = forcage ou non par la Ts
     237!Config  Def  = false
     238!Config  Help = forcage ou non par la Ts
     239       ok_forc_tsurf=.false.
     240       CALL getin('ok_forc_tsurf',ok_forc_tsurf)
     241
    235242!Config  Key  = ok_prescr_ust
    236243!Config  Desc = ustar impose ou non
    … …  
    239246       ok_prescr_ust = .false.
    240247       CALL getin('ok_prescr_ust',ok_prescr_ust)
     248
     249
     250!Config  Key  = ok_prescr_beta
     251!Config  Desc = betaevap impose ou non
     252!Config  Def  = false
     253!Config  Help = betaevap impose ou non
     254       ok_prescr_beta = .false.
     255       CALL getin('ok_prescr_beta',ok_prescr_beta)
    241256
    242257!Config  Key  = ok_old_disvert
    … …  
    280295!Config  Desc = surface temperature
    281296!Config  Def  = 290.
    282 !Config  Help = not used if type_ts_forcing=1 in lmdz1d.F
     297!Config  Help = surface temperature
    283298       tsurf = 290.
    284299       CALL getin('tsurf',tsurf)
    … …  
    297312       zsurf = 0.
    298313       CALL getin('zsurf',zsurf)
     314! EV pour accord avec format standard       
     315       CALL getin('zorog',zsurf)
     316
    299317
    300318!Config  Key  = rugos
    … …  
    359377       qsolinp = 1.
    360378       CALL getin('qsolinp',qsolinp)
     379
     380
     381
     382!Config  Key  = betaevap
     383!Config  Desc = beta for actual evaporation when prescribed
     384!Config  Def  = 1.0
     385!Config  Help =
     386       betaevap = 1.
     387       CALL getin('betaevap',betaevap)     
    361388
    362389!Config  Key  = zpicinp
    … …  
    520547       CALL getin('forc_ustar',forc_ustar)
    521548       IF (forc_ustar .EQ. 1) ok_prescr_ust=.true.
     549
    522550
    523551!Config  Key  = nudging_u
    … …  
    12481276      END
    12491277
    1250 !======================================================================
    1251        SUBROUTINE read_tsurf1d(knon,sst_out)
    1252 
    1253 ! This subroutine specifies the surface temperature to be used in 1D simulations
    1254 
    1255       USE dimphy, ONLY : klon
    1256 
    1257       INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
    1258       REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out  ! tsurf used to force the single-column model
    1259 
    1260        INTEGER :: i
    1261 ! COMMON defined in lmdz1d.F:
    1262        real ts_cur
    1263        common /sst_forcing/ts_cur
    1264 
    1265        DO i = 1, knon
    1266         sst_out(i) = ts_cur
    1267        ENDDO
    1268 
    1269       END SUBROUTINE read_tsurf1d
    1270 
     1278!!======================================================================
     1279!       SUBROUTINE read_tsurf1d(knon,sst_out)
     1280!
     1281!! This subroutine specifies the surface temperature to be used in 1D simulations
     1282!
     1283!      USE dimphy, ONLY : klon
     1284!
     1285!      INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
     1286!      REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out  ! tsurf used to force the single-column model
     1287!
     1288!       INTEGER :: i
     1289!! COMMON defined in lmdz1d.F:
     1290!       real ts_cur
     1291!       common /sst_forcing/ts_cur
     1292
     1293!       DO i = 1, knon
     1294!        sst_out(i) = ts_cur
     1295!       ENDDO
     1296!
     1297!      END SUBROUTINE read_tsurf1d
     1298!
    12711299!===============================================================
    12721300      subroutine advect_vert(llm,w,dt,q,plev)

  • LMDZ6/trunk/libf/phylmd/dyn1d/1D_decl_cases.h

    r3686 r3780  
    3434        real w_mod(llm), t_mod(llm),q_mod(llm)
    3535        real u_mod(llm),v_mod(llm), ht_mod(llm),vt_mod(llm),ug_mod(llm),vg_mod(llm)
    36         real temp_nudg_mod(llm),qv_nudg_mod(llm),u_nudg_mod(llm),v_nudg_mod(llm)
     36            real temp_nudg_mod(llm),qv_nudg_mod(llm),u_nudg_mod(llm),v_nudg_mod(llm)
    3737        real hq_mod(llm),vq_mod(llm),qv_mod(llm),ql_mod(llm),qt_mod(llm)
    3838        real th_mod(llm)
    3939
    40         real ts_cur
    41         common /sst_forcing/ts_cur ! also in read_tsurf1d.F
     40! EV comment these lines
     41!        real ts_cur
     42!        common /sst_forcing/ts_cur ! also in read_tsurf1d.F
    4243!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    4344! Declarations specifiques au cas RICO

  • LMDZ6/trunk/libf/phylmd/dyn1d/1D_interp_cases.h

    r3686 r3780  
    11
    2          print*,'FORCING CASE forcing_case2'
     2        print*,'FORCING CASE forcing_case2'
    33!       print*,                                                             &
    44!    & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/pdt_cas=',     &
    … …  
    2828     &       ,dth_prof_cas,hth_prof_cas,vth_prof_cas,lat_prof_cas                           &
    2929     &       ,sens_prof_cas,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas,tke_prof_cas)
    30 
    31              ts_cur = ts_prof_cas
     30! EV tg instead of ts_cur
     31             tg = ts_prof_cas
    3232!            psurf=plev_prof_cas(1)
    3333             psurf=ps_prof_cas

  • LMDZ6/trunk/libf/phylmd/dyn1d/1D_read_forc_cases.h

    r3686 r3780  
    7070
    7171! initial and boundary conditions :
    72 !      tsurf = ts_prof_cas
     72!     tsurf = ts_prof_cas
    7373      psurf = ps_prof_cas
    74       ts_cur = ts_prof_cas
     74      !EV tg instead of ts_cur
     75      tg = ts_prof_cas
     76      print*, 'tg=', tg
     77
    7578      do l = 1, llm
    7679       temp(l) = t_mod_cas(l)
    … …  
    108111       IF (ok_prescr_ust) THEN
    109112       ust=ustar_prof_cas
    110        print *,'ust=',ust
    111113       ENDIF
    112114
     115
    113116      endif !forcing_SCM

  • LMDZ6/trunk/libf/phylmd/dyn1d/mod_1D_cases_read_std.F90

    r3688 r3780  
    9292      REAL, ALLOCATABLE :: time_val(:)
    9393
    94       print*,'ON EST VRAIMENT LA'
     94      print*,'ON EST VRAIMENT DASN MOD_1D_CASES_READ_STD'
    9595      fich_cas='cas.nc'
    9696      print*,'fich_cas ',fich_cas
    … …  
    924924!      enddo
    925925
     926!       print*, 'plev_prof_cas', plev_prof_cas
     927!       print*, 'play', play
    926928       do l = 1, llm
    927929
    … …  
    951953
    952954         frac = (plev_prof_cas(k2)-play(l))/(plev_prof_cas(k2)-plev_prof_cas(k1))
     955         
    953956         t_mod_cas(l)= t_prof_cas(k2) - frac*(t_prof_cas(k2)-t_prof_cas(k1))
    954957         theta_mod_cas(l)= th_prof_cas(k2) - frac*(th_prof_cas(k2)-th_prof_cas(k1))
    … …  
    10751078       enddo ! l
    10761079
     1080
    10771081          return
    10781082          end SUBROUTINE interp2_case_vertical_std

  • LMDZ6/trunk/libf/phylmd/dyn1d/old_1D_decl_cases.h

    r3593 r3780  
    3737        real th_mod(llm)
    3838
    39         real ts_cur
    40         common /sst_forcing/ts_cur ! also in read_tsurf1d.F
     39        !real ts_cur
     40        !common /sst_forcing/ts_cur ! also in read_tsurf1d.F
    4141!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    4242! Declarations specifiques au cas RICO

  • LMDZ6/trunk/libf/phylmd/dyn1d/old_1D_interp_cases.h

    r3593 r3780  
    6262     &             ,ts_prof,plev_prof,t_prof,q_prof,u_prof,v_prof,w_prof    &
    6363     &             ,ht_prof,vt_prof,hq_prof,vq_prof)
    64 
    65         if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
     64! EV: tg instead of ts_cur
     65        if (type_ts_forcing.eq.1) tg = ts_prof !
    6666
    6767! vertical interpolation:
    … …  
    113113!     print *,'llm l omega_profd',llm,l,omega_profd(l)
    114114!     enddo
    115 
    116         if (type_ts_forcing.eq.1) ts_cur = tg_prof ! SST used in read_tsurf1d
     115! EV tg instead of ts_cur
     116        if (type_ts_forcing.eq.1) tg = tg_prof ! SST used
    117117
    118118! vertical interpolation:
    … …  
    206206     &             ,ug_gabls4,vg_gabls4,ht_gabls4,hq_gabls4,tg_gabls4                            &
    207207     &             ,ug_profg,vg_profg,ht_profg,hq_profg,tg_profg)
    208 
    209         if (type_ts_forcing.eq.1) ts_cur = tg_prof ! SST used in read_tsurf1d
     208!EV tg instead of ts_cur
     209        if (type_ts_forcing.eq.1) tg = tg_prof ! SST used
    210210
    211211! vertical interpolation:
    … …  
    499499     &             ,nlev_sandu                                              &
    500500     &             ,ts_sandu,ts_prof)
    501 
    502         if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
     501! EV tg instead of ts_cur
     502        if (type_ts_forcing.eq.1) tg = ts_prof ! SST used in read_tsurf1d
    503503
    504504! vertical interpolation:
    … …  
    582582     &             ,ufa_astex,vfa_astex,div_prof,ts_prof,ug_prof,vg_prof    &
    583583     &             ,ufa_prof,vfa_prof)
    584 
    585         if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
    586 
     584! EV tg instead of ts_cur
     585        if (type_ts_forcing.eq.1) tg = ts_prof ! SST used
    587586! vertical interpolation:
    588587      CALL interp_astex_vertical(play,nlev_astex,plev_profa                 &
    … …  
    675674     &       ,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas               &
    676675     &       ,sens_prof_cas,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas)
    677 
    678              ts_cur = ts_prof_cas
     676! EV tg instead of ts_cur
     677
     678             tg = ts_prof_cas
    679679             psurf=plev_prof_cas(1)
    680680
    … …  
    850850     &       ,dth_prof_cas,hth_prof_cas,vth_prof_cas,lat_prof_cas                           &
    851851     &       ,sens_prof_cas,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas,tke_prof_cas)
    852 
    853              ts_cur = ts_prof_cas
     852! EV tg instead of ts_cur
     853
     854             tg = ts_prof_cas
    854855!            psurf=plev_prof_cas(1)
    855856             psurf=ps_prof_cas

  • LMDZ6/trunk/libf/phylmd/dyn1d/old_1D_read_forc_cases.h

    r3679 r3780  
    875875
    876876! initial and boundary conditions :
    877 !      tsurf = ts_prof_cas
    878       ts_cur = ts_prof_cas
     877!     tsurf = ts_prof_cas
     878! EV tg instead of ts_cur
     879      tg= ts_prof_cas
    879880      psurf=plev_prof_cas(1)
    880881      write(*,*) 'SST initiale: ',tsurf
    … …  
    965966! initial and boundary conditions :
    966967!      tsurf = ts_prof_cas
    967       ts_cur = ts_prof_cas
     968! EV tg instead of ts_cur
     969      tg = ts_prof_cas
    968970      psurf=plev_prof_cas(1)
    969971      write(*,*) 'SST initiale: ',tsurf
    … …  
    10631065! initial and boundary conditions :
    10641066!      tsurf = ts_prof_cas
    1065       ts_cur = ts_prof_cas
     1067! EV tg instead of ts_cur
     1068
     1069      tg = ts_prof_cas
    10661070      psurf=plev_prof_cas(1)
    10671071      write(*,*) 'SST initiale: ',tsurf

  • LMDZ6/trunk/libf/phylmd/dyn1d/old_lmdz1d.F90

    r3594 r3780  
    728728
    729729!Al1 pour SST forced, appell?? depuis ocean_forced_noice
    730       ts_cur = tsurf ! SST used in read_tsurf1d
     730! EV tg instead of ts_cur
     731
     732      tg = tsurf ! SST used in read_tsurf1d
    731733!=====================================================================
    732734! Initialisation de la physique :
    … …  
    791793
    792794        fder=0.
     795        print *, 'snsrf', snsrf
    793796        snsrf(1,:)=snowmass ! masse de neige des sous surface
    794797        qsurfsrf(1,:)=qsurf ! humidite de l'air des sous surface
    … …  
    841844     end if
    842845
    843 
    844846        print*,'nat_surf,pctsrf(1,is_oce),pctsrf(1,is_ter)',nat_surf         &
    845847     &        ,pctsrf(1,is_oce),pctsrf(1,is_ter)
    … …  
    851853              ! 6 albedo, mais on peut quand meme tourner avec
    852854              ! moins. Seules les 2 ou 4 premiers seront lus
     855              print*, 'les albedos sont sont', albedo, falb_dir
    853856        falb_dir=albedo
    854857        falb_dif=albedo
     858             print*, falb_dir
    855859        rugoro=rugos
    856860        t_ancien(1,:)=temp(:)
    … …  
    913917        v_ancien(1,:)=v(:)
    914918
    915 u10m=0.
    916 v10m=0.
    917 ale_wake=0.
    918 ale_bl_stat=0.
     919        u10m=0.
     920        v10m=0.
     921        ale_wake=0.
     922        ale_bl_stat=0.
    919923 
    920924!------------------------------------------------------------------------

  • LMDZ6/trunk/libf/phylmd/dyn1d/scm.F90

    r3693 r3780  
    7575      real :: zcufi    = 1.
    7676      real :: zcvfi    = 1.
    77 
    78 !-      real :: nat_surf
    79 !-      logical :: ok_flux_surf
    80 !-      real :: fsens
    81 !-      real :: flat
    82 !-      real :: tsurf
    83 !-      real :: rugos
    84 !-      real :: qsol(1:2)
    85 !-      real :: qsurf
    86 !-      real :: psurf
    87 !-      real :: zsurf
    88 !-      real :: albedo
    89 !-
    90 !-      real :: time     = 0.
    91 !-      real :: time_ini
    92 !-      real :: xlat
    93 !-      real :: xlon
    94 !-      real :: wtsurf
    95 !-      real :: wqsurf
    96 !-      real :: restart_runoff
    97 !-      real :: xagesno
    98 !-      real :: qsolinp
    99 !-      real :: zpicinp
    100 !-
    10177      real :: fnday
    10278      real :: day, daytime
    … …  
    141117        logical :: forcing_case2   = .false.
    142118        logical :: forcing_SCM   = .false.
    143         integer :: type_ts_forcing ! 0 = SST constant; 1 = SST read from a file
    144 !                                                            (cf read_tsurf1d.F)
    145119
    146120!flag forcings
    … …  
    148122        logical :: nudge_thermo=.false.
    149123        logical :: cptadvw=.true.
     124
     125
    150126!=====================================================================
    151127! DECLARATIONS FOR EACH CASE
    … …  
    248224!
    249225      integer :: it_end ! iteration number of the last call
    250 !Al1
     226!Al1,plev,play,phi,phis,presnivs,
    251227      integer ecrit_slab_oc !1=ecrit,-1=lit,0=no file
    252228      data ecrit_slab_oc/-1/
    … …  
    278254      d_v_age(:)=0.
    279255     
     256     
    280257! Initialization of Common turb_forcing
    281258       dtime_frcg = 0.
    … …  
    290267! OPTIONS OF THE 1D SIMULATION (lmdz1d.def => unicol.def)
    291268!---------------------------------------------------------------------
    292 !Al1
    293269        call conf_unicol
    294270!Al1 moves this gcssold var from common fcg_gcssold to
    … …  
    296272! --------------------------------------------------------------------
    297273        close(1)
    298 !Al1
    299274        write(*,*) 'lmdz1d.def lu => unicol.def'
    300275
    … …  
    302277       year_ini_cas=1997
    303278       ! It is possible that those parameters are run twice.
    304 
    305279       ! A REVOIR : LIRE PEUT ETRE AN MOIS JOUR DIRECETEMENT
     280
     281
    306282       call getin('anneeref',year_ini_cas)
    307283       call getin('dayref',day_deb)
    … …  
    309285       call getin('time_ini',heure_ini_cas)
    310286
    311         type_ts_forcing = 0
    312         IF (nat_surf==0) type_ts_forcing=1 ! SST forcee sur OCEAN
    313         print*,'NATURE DE LA SURFACE ',nat_surf
     287      print*,'NATURE DE LA SURFACE ',nat_surf
    314288!
    315289! Initialization of the logical switch for nudging
     290
    316291     jcode = iflag_nudge
    317292     do i = 1,nudge_max
    … …  
    319294       jcode = jcode/10
    320295     enddo
    321 !---------------------------------------------------------------------
     296!-----------------------------------------------------------------------
    322297!  Definition of the run
    323 !---------------------------------------------------------------------
     298!-----------------------------------------------------------------------
    324299
    325300      call conf_gcm( 99, .TRUE. )
    … …  
    343318      allocate( phy_flic(year_len)) ! Fraction de glace
    344319      phy_flic(:)=0.0
     320
     321
    345322!-----------------------------------------------------------------------
    346323!   Choix du calendrier
    … …  
    373350!      Le numero du jour est dans "day". L heure est traitee separement.
    374351!      La date complete est dans "daytime" (l'unite est le jour).
     352
     353
    375354      if (nday>0) then
    376355         fnday=nday
    … …  
    409388! Initialization of dimensions, geometry and initial state
    410389!---------------------------------------------------------------------
    411 !      call init_phys_lmdz(1,1,llm,1,(/1/)) ! job now done via iniphysiq
     390!     call init_phys_lmdz(1,1,llm,1,(/1/)) ! job now done via iniphysiq
    412391!     but we still need to initialize dimphy module (klon,klev,etc.)  here.
    413392      call init_dimphy1D(1,llm)
    … …  
    446425!!! Feedback forcing values for Gateaux differentiation (al1)
    447426!!!=====================================================================
    448 !!! Surface Planck forcing bracketing call radiation
    449 !!      surf_Planck = 0.
    450 !!      surf_Conv   = 0.
    451 !!      write(*,*) 'Gateaux-dif Planck,Conv:',surf_Planck,surf_Conv
    452 !!! a mettre dans le lmdz1d.def ou autre
    453 !!
    454427!!
    455428      qsol = qsolinp
    … …  
    469442      ENDIF
    470443      print*,'Flux sol ',fsens,flat
    471 !!      ok_flux_surf=.false.
    472 !!      fsens=-wtsurf*rcpd*rho(1)
    473 !!      flat=-wqsurf*rlvtt*rho(1)
    474 !!!!
    475444
    476445! Vertical discretization and pressure levels at half and mid levels:
    … …  
    496465      plev =ap+bp*psurf
    497466      play = 0.5*(plev(1:llm)+plev(2:llm+1))
    498       zlay=-rd*300.*log(play/psurf)/rg ! moved after reading profiles
     467      zlay=-rd*300.*log(play/psurf)/rg ! moved after reading profiles.
    499468
    500469      IF (forcing_type .eq. 59) THEN
    … …  
    527496      print*,'mxcalc=',mxcalc
    528497!     print*,'zlay=',zlay(mxcalc)
    529       print*,'play=',play(mxcalc)
    530 
    531 !Al1 pour SST forced, appell?? depuis ocean_forced_noice
    532       ts_cur = tsurf ! SST used in read_tsurf1d
     498!      print*,'play=',play(mxcalc)
     499
     500!! When surface temperature is forced
     501      tg= tsurf ! surface T used in read_tsurf1d
     502
     503
    533504!=====================================================================
    534505! Initialisation de la physique :
    … …  
    546517! airefi,zcufi,zcvfi initialises au debut de ce programme
    547518! rday,ra,rg,rd,rcpd declares dans YOMCST.h et calcules dans suphel.F
     519
     520
    548521      day_step = float(nsplit_phys)*day_step/float(iphysiq)
    549522      write (*,*) 'Time step divided by nsplit_phys (=',nsplit_phys,')'
    … …  
    563536     ! e.g. for cell boundaries, which are meaningless in 1D; so pad these
    564537     ! with '0.' when necessary
     538
    565539      call iniphysiq(iim,jjm,llm, &
    566540           1,comm_lmdz, &
    … …  
    653627              ! 6 albedo, mais on peut quand meme tourner avec
    654628              ! moins. Seules les 2 ou 4 premiers seront lus
     629              print*, 'les albedos sont', albedo
    655630        falb_dir=albedo
    656631        falb_dif=albedo
    … …  
    664639        prw_ancien = 0.
    665640!jyg<
    666 !!        pbl_tke(:,:,:)=1.e-8
     641!!      pbl_tke(:,:,:)=1.e-8
    667642        pbl_tke(:,:,:)=0.
    668         pbl_tke(:,2,:)=1.e-2
     643! EV: pourquoi????
     644!        pbl_tke(:,2,:)=1.e-2
    669645        PRINT *, ' pbl_tke dans lmdz1d '
    670646        if (prt_level .ge. 5) then
    … …  
    675651
    676652!>jyg
    677 
    678653        rain_fall=0.
    679654        snow_fall=0.
    … …  
    715690        v_ancien(1,:)=v(:)
    716691 
    717 u10m=0.
    718 v10m=0.
    719 ale_wake=0.
    720 ale_bl_stat=0.
     692        u10m=0.
     693        v10m=0.
     694        ale_wake=0.
     695        ale_bl_stat=0.
    721696
    722697!------------------------------------------------------------------------
    … …  
    738713! to be set at some arbitratry convenient values.
    739714!------------------------------------------------------------------------
    740 !Al1 =============== restart option ==========================
     715!Al1 =============== restart option ======================================
    741716        if (.not.restart) then
    742717          iflag_pbl = 5
    … …  
    803778       print*,'plev,play,phi,phis,presnivs,u,v,temp,q,omega2'
    804779       print*,'temp(1),q(1,1),u(1),v(1),plev(1),phis :'
    805        print*,temp(1),q(1,1),u(1),v(1),plev(1),phis
     780       print*,temp(1),q(1,1),u(1),v(1),plev(1),phis(1)
    806781! raz for safety
    807782       do l=1,llm
    … …  
    809784       enddo
    810785      endif
    811 !Al1 ================  end restart =================================
     786!======================  end restart =================================
    812787      IF (ecrit_slab_oc.eq.1) then
    813788         open(97,file='div_slab.dat',STATUS='UNKNOWN')
    … …  
    820795       CALL iophys_ini
    821796#endif
     797
     798!=====================================================================
    822799! START OF THE TEMPORAL LOOP :
    823800!=====================================================================
    824801           
    825802      it_end = nint(fnday*day_step)
    826 !test JLD     it_end = 10
    827803      do while(it.le.it_end)
    828804
    … …  
    832808         print*,'PAS DE TEMPS ',timestep
    833809       endif
    834 !Al1 demande de restartphy.nc
    835810       if (it.eq.it_end) lastcall=.True.
    836811
    … …  
    844819!  Geopotential :
    845820!---------------------------------------------------------------------
    846 
     821!        phis(1)=zsurf*RG
     822!        phi(1)=phis(1)+RD*temp(1)*(plev(1)-play(1))/(.5*(plev(1)+play(1)))
    847823        phi(1)=RD*temp(1)*(plev(1)-play(1))/(.5*(plev(1)+play(1)))
     824
    848825        do l = 1, llm-1
    849826          phi(l+1)=phi(l)+RD*(temp(l)+temp(l+1))*                           &
    850827     &    (play(l)-play(l+1))/(play(l)+play(l+1))
    851828        enddo
     829
    852830
    853831!---------------------------------------------------------------------
    … …  
    950928       sfdt = sin(0.5*fcoriolis*timestep)
    951929       cfdt = cos(0.5*fcoriolis*timestep)
    952 !       print *,'fcoriolis,sfdt,cfdt,timestep',fcoriolis,sfdt,cfdt,timestep
    953 !
     930
    954931        d_u_age(1:mxcalc)= -2.*sfdt/timestep*                                &
    955932     &          (sfdt*(u(1:mxcalc)-ug(1:mxcalc)) -                          &
    … …  
    10301007        temp(1:mxcalc)=temp(1:mxcalc)+timestep*(                            &
    10311008     &              dt_phys(1:mxcalc)                                       &
    1032      &             +d_t_adv(1:mxcalc)                                      &
    1033      &             +d_t_nudge(1:mxcalc)                                      &
     1009     &             +d_t_adv(1:mxcalc)                                       &
     1010     &             +d_t_nudge(1:mxcalc)                                     &
    10341011     &             +dt_cooling(1:mxcalc))  ! Taux de chauffage ou refroid.
    10351012
    10361013
    1037 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1014!=======================================================================
    10381015!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
    1039 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    1040 !     endif  ! forcing_sandu or forcing_astex
    1041 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1016!=======================================================================
    10421017
    10431018        teta=temp*(pzero/play)**rkappa
    1044 !
     1019
    10451020!---------------------------------------------------------------------
    10461021!   Nudge soil temperature if requested
    … …  
    10801055
    10811056!  incremente day time
    1082 !        print*,'daytime bef',daytime,1./day_step
    10831057        daytime = daytime+1./day_step
    1084 !Al1dbg
    10851058        day = int(daytime+0.1/day_step)
    10861059!        time = max(daytime-day,0.0)
    … …  
    10881061!cc        time = real(mod(it,day_step))/day_step
    10891062        time = time_ini/24.+real(mod(it,day_step))/day_step
    1090 !        print*,'daytime nxt time',daytime,time
    10911063        it=it+1
    10921064
    10931065      enddo
    10941066
    1095 !Al1
    10961067      if (ecrit_slab_oc.ne.-1) close(97)
    10971068
    10981069!Al1 Call to 1D equivalent of dynredem (an,mois,jour,heure ?)
    1099 ! -------------------------------------
     1070! ---------------------------------------------------------------------------
    11001071       call dyn1dredem("restart1dyn.nc",                                    &
    11011072     &              plev,play,phi,phis,presnivs,                            &

  • LMDZ6/trunk/libf/phylmd/flux_arp.h

    r2329 r3780  
    11!
    22! $Id: flux_arp.h 2010-08-04 17:02:56Z lahellec $
     3! Modif EV, 10/2020
    34!
    45      logical :: ok_flux_surf
    56      logical :: ok_prescr_ust !for prescribed ustar
     7      logical :: ok_prescr_beta
     8      logical :: ok_forc_tsurf
     9
     10
    611      real :: fsens
    712      real :: flat
     13      real :: betaevap
    814      real :: ust
    915      real :: tg
    1016
    11       common /flux_arp/fsens,flat,ust,tg,ok_flux_surf,ok_prescr_ust
     17      common /flux_arp/fsens,flat,ust,tg,ok_flux_surf,ok_prescr_ust,ok_prescr_beta,betaevap,ok_forc_tsurf
    1218
    1319!$OMP THREADPRIVATE(/flux_arp/)
    1420
    1521
     22     
    1623
    17 

  • LMDZ6/trunk/libf/phylmd/ocean_forced_mod.F90

    r3327 r3780  
    3838    INCLUDE "YOMCST.h"
    3939    INCLUDE "clesphys.h"
    40 
     40    INCLUDE "flux_arp.h"
    4141
    4242! Input arguments
    … …  
    7676    REAL, DIMENSION(klon)       :: u1_lay, v1_lay
    7777    LOGICAL                     :: check=.FALSE.
    78     REAL, DIMENSION(klon) :: sens_prec_liq, sens_prec_sol   
    79     REAL, DIMENSION(klon) :: lat_prec_liq, lat_prec_sol   
     78    REAL, DIMENSION(klon)       :: sens_prec_liq, sens_prec_sol   
     79    REAL, DIMENSION(klon)       :: lat_prec_liq, lat_prec_sol   
    8080
    8181!****************************************************************************************
    … …  
    9292!!jyg    if (knon.eq.1) then ! single-column model
    9393    if (klon_glo.eq.1) then ! single-column model
    94       CALL read_tsurf1d(knon,tsurf_lim) ! new
     94      ! EV: now surface Tin flux_arp.h
     95      !CALL read_tsurf1d(knon,tsurf_lim) ! new
     96       DO i = 1, knon
     97        tsurf_lim(i) = tg
     98       ENDDO
     99
    95100    else ! GCM
    96101      CALL limit_read_sst(knon,knindex,tsurf_lim)
    … …  
    104109!****************************************************************************************
    105110! Set some variables for calcul_fluxs
    106     cal = 0.
    107     beta = 1.
    108     dif_grnd = 0.
     111    !cal = 0.
     112    !beta = 1.
     113    !dif_grnd = 0.
     114    ! EV: use calbeta to calculate beta
     115   
     116    CALL calbeta(dtime, is_oce, knon, snow, beta*0., beta, cal, dif_grnd)
     117
     118
    109119    alb_neig(:) = 0.
    110120    agesno(:) = 0.
    … …  
    172182    INCLUDE "YOMCST.h"
    173183    INCLUDE "clesphys.h"
     184    INCLUDE "flux_arp.h"
    174185
    175186! Input arguments
    … …  
    233244    tsurf_tmp(:) = tsurf_in(:)
    234245
    235 ! calculate the parameters cal, beta, capsol and dif_grnd
    236     CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, capsol, dif_grnd)
     246! calculate the parameters cal, beta, capsol and dif_grnd and then recalculate cal
     247    CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, cal, dif_grnd)
    237248
    238249   
    … …  
    249260    ENDIF
    250261
    251     beta = 1.0
     262    !beta = 1.0
    252263    sens_prec_liq = 0.; sens_prec_sol = 0.; lat_prec_liq = 0.; lat_prec_sol = 0.
    253264
    … …  
    304315  END SUBROUTINE ocean_forced_ice
    305316
     317!
    306318!************************************************************************
    307 ! 1D case
    308 !************************************************************************
    309   SUBROUTINE read_tsurf1d(knon,sst_out)
    310 
    311 ! This subroutine specifies the surface temperature to be used in 1D simulations
    312 
    313       USE dimphy, ONLY : klon
    314 
    315       INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
    316       REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out  ! tsurf used to force the single-column model
    317 
    318        INTEGER :: i
    319 ! COMMON defined in lmdz1d.F:
    320        real ts_cur
    321        common /sst_forcing/ts_cur
    322 
    323        DO i = 1, knon
    324         sst_out(i) = ts_cur
    325        ENDDO
    326 
    327       END SUBROUTINE read_tsurf1d
    328 
    329 !
    330 !************************************************************************
    331319!
    332320END MODULE ocean_forced_mod

  • LMDZ6/trunk/libf/phylmd/ocean_slab_mod.F90

    r3102 r3780  
    421421!   
    422422!****************************************************************************************
    423     cal(:)      = 0. ! infinite thermal inertia
    424     beta(:)     = 1. ! wet surface
    425     dif_grnd(:) = 0. ! no diffusion into ground
     423    !cal(:)      = 0. ! infinite thermal inertia
     424    !beta(:)     = 1. ! wet surface
     425    !dif_grnd(:) = 0. ! no diffusion into ground
     426    ! EV: use calbeta
     427    CALL calbeta(dtime, is_oce, knon, snow,qsurf, beta, cal, dif_grnd)
     428
     429
    426430   
    427431! Suppose zero surface speed
    … …  
    742746! set beta, cal, compute conduction fluxes inside ice/snow
    743747    slab_bilg(:)=0.
    744     dif_grnd(:)=0.
    745     beta(:) = 1.
     748    !dif_grnd(:)=0.
     749    !beta(:) = 1.
     750    ! EV: use calbeta to calculate beta and then recalculate properly cal
     751    CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, cal, dif_grnd)
     752
     753
    746754    DO i=1,knon
    747755    ki=knindex(i)

  • LMDZ6/trunk/libf/phylmd/pbl_surface_mod.F90

    r3774 r3780  
    206206!jyg<
    207207!!       zxfluxt,   zxfluxq,   q2m,      flux_q, tke,   &
    208        zxfluxt,   zxfluxq,   q2m,      flux_q, tke_x,   &
     208       zxfluxt,   zxfluxq,   q2m,      flux_q, tke_x,  &
    209209!>jyg
    210210!!! nrlmd+jyg le 02/05/2011 et le 20/02/2012
    … …  
    504504    CHARACTER(len=8), DIMENSION(nbsrf), SAVE :: cl_surf
    505505!$OMP THREADPRIVATE(cl_surf)
    506     REAL, SAVE                               :: beta_land         ! beta for wx_dts
    507 !$OMP THREADPRIVATE(beta_land)
     506 ! EV Ne sert plus:
     507 !   REAL, SAVE                               :: beta_land         ! beta for wx_dts
     508!!$OMP THREADPRIVATE(beta_land)
    508509
    509510! Other local variables
    … …  
    845846       ! Initialize ok_flux_surf (for 1D model)
    846847       if (klon_glo>1) ok_flux_surf=.FALSE.
     848       if (klon_glo>1) ok_forc_tsurf=.FALSE.
     849
    847850
    848851       ! intialize beta_land
    849        beta_land = 0.5
    850        call getin_p('beta_land', beta_land)
     852       !beta_land = 0.5
     853       !call getin_p('beta_land', beta_land)
    851854       
    852855       ! Initilize debug IO
    … …  
    947950!!    tke(:,:,is_ave)=0.
    948951    tke_x(:,:,is_ave)=0.
     952
    949953    wake_dltke(:,:,is_ave)=0.
    950954!>jyg
    … …  
    978982!!    d_t_diss= 0.0 ;d_u = 0.0     ; d_v = 0.0
    979983    yqsol = 0.0   
     984
    980985    ytke=0.
    981986!FC
    … …  
    13881393             ytke_w(j,k)      = tke_x(i,k,nsrf)+wake_dltke(i,k,nsrf)
    13891394             ywake_dltke(j,k) = wake_dltke(i,k,nsrf)
     1395           
    13901396!>jyg
    13911397          ENDDO
    … …  
    14621468      ENDDO
    14631469     ENDIF
     1470
    14641471        IF (prt_level >=10) print *,'clcdrag -> ycdragh ', ycdragh
    14651472       ELSE  !(iflag_split .eq.0)
    … …  
    15451552      print *,' args coef_diff_turb: ycdragh ', ycdragh
    15461553      print *,' args coef_diff_turb: ytke ', ytke
     1554
    15471555       ENDIF
    15481556        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    … …  
    15741582      print *,' args coef_diff_turb: ycdragh_x ', ycdragh_x
    15751583      print *,' args coef_diff_turb: ytke_x ', ytke_x
     1584
    15761585       ENDIF
    15771586        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    … …  
    20202029!
    20212030!****************************************************************************************
    2022 
    2023 !!!
    2024 !!! jyg le 10/04/2013
     2031!!
     2032!!!
     2033!!! jyg le 10/04/2013 et EV 10/2020
     2034
     2035        IF (ok_forc_tsurf) THEN
     2036            DO j=1,knon
     2037                ytsurf_new(j)=tg
     2038                y_d_ts(j) = ytsurf_new(j) - yts(j)
     2039            ENDDO
     2040        ENDIF ! ok_forc_tsurf
     2041
    20252042!!!
    20262043        IF (ok_flux_surf) THEN
    … …  
    24512468              tke_x(i,k,nsrf)    = ytke(j,k)
    24522469              tke_x(i,k,is_ave) = tke_x(i,k,is_ave) + ytke(j,k)*ypct(j)
     2470
    24532471!>jyg
    24542472           ENDDO
    … …  
    24642482!!            tke(i,k,is_ave) = tke(i,k,is_ave) + tke(i,k,nsrf)*ypct(j)
    24652483            tke_x(i,k,nsrf)   = ytke_x(j,k)
    2466             tke_x(i,k,is_ave)   = tke_x(i,k,is_ave) + tke_x(i,k,nsrf)*ypct(j)
     2484            tke_x(i,k,is_ave)   = tke_x(i,k,is_ave) + tke_x(i,k,nsrf)*ypct(j)       
    24672485            wake_dltke(i,k,is_ave)   = wake_dltke(i,k,is_ave) + wake_dltke(i,k,nsrf)*ypct(j)
     2486           
    24682487
    24692488!>jyg

  • LMDZ6/trunk/libf/phylmd/phys_local_var_mod.F90

    r3779 r3780  
    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)
    … …  
    445444      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: ref_liq_pi, ref_ice_pi
    446445!$OMP THREADPRIVATE(ref_liq_pi, ref_ice_pi)
    447       REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zx_rh
    448 !$OMP THREADPRIVATE(zx_rh)
     446      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zx_rh, zx_rhl, zx_rhi
     447!$OMP THREADPRIVATE(zx_rh, zx_rhl, zx_rhi)
    449448      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: prfl, psfl, fraca
    450449!$OMP THREADPRIVATE(prfl, psfl, fraca)
    … …  
    561560      ALLOCATE(t_seri(klon,klev),q_seri(klon,klev),ql_seri(klon,klev),qs_seri(klon,klev))
    562561      ALLOCATE(u_seri(klon,klev),v_seri(klon,klev))
    563       ALLOCATE(l_mixmin(klon,klev+1,nbsrf), l_mix(klon,klev+1,nbsrf))
    564       l_mix(:,:,:)=0. ; l_mixmin(:,:,:)=0. ! doit etre initialse car pas toujours remplis
     562      ALLOCATE(l_mixmin(klon,klev+1,nbsrf), l_mix(klon,klev+1,nbsrf), tke_dissip(klon,klev+1,nbsrf))
     563      l_mix(:,:,:)=0. ; l_mixmin(:,:,:)=0. ; tke_dissip(:,:,:)=0. ! doit etre initialse car pas toujours remplis
    565564
    566565      ALLOCATE(tr_seri(klon,klev,nbtr))
    … …  
    780779      ALLOCATE(ref_liq(klon, klev), ref_ice(klon, klev), theta(klon, klev))
    781780      ALLOCATE(ref_liq_pi(klon, klev), ref_ice_pi(klon, klev))
    782       ALLOCATE(zphi(klon, klev), zx_rh(klon, klev))
     781      ALLOCATE(zphi(klon, klev), zx_rh(klon, klev), zx_rhl(klon,klev), zx_rhi(klon,klev))
    783782      ALLOCATE(pmfd(klon, klev), pmfu(klon, klev))
    784783
    … …  
    879878      DEALLOCATE(t_seri,q_seri,ql_seri,qs_seri)
    880879      DEALLOCATE(u_seri,v_seri)
    881       DEALLOCATE(l_mixmin,l_mix)
     880      DEALLOCATE(l_mixmin,l_mix, tke_dissip)
    882881
    883882      DEALLOCATE(tr_seri)
    … …  
    10741073      DEALLOCATE(ref_liq, ref_ice, theta)
    10751074      DEALLOCATE(ref_liq_pi, ref_ice_pi)
    1076       DEALLOCATE(zphi, zx_rh)
     1075      DEALLOCATE(zphi, zx_rh, zx_rhl, zx_rhi)
    10771076      DEALLOCATE(pmfd, pmfu)
    10781077

  • LMDZ6/trunk/libf/phylmd/phys_output_ctrlout_mod.F90

    r3775 r3780  
    10051005  TYPE(ctrl_out), SAVE :: o_tke = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    10061006    'tke ', 'TKE', 'm2/s2', (/ ('', i=1, 10) /))
     1007  TYPE(ctrl_out), SAVE :: o_tke_dissip = ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1008    'tke_dissip ', 'TKE DISSIPATION', 'm2/s3', (/ ('', i=1, 10) /))   
    10071009  TYPE(ctrl_out), SAVE :: o_tke_max = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    10081010    'tke_max', 'TKE max', 'm2/s2',                                  &
    10091011      (/ 't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)', &
    10101012         't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)', 't_max(X)' /))
    1011 
    10121013  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_tke_srf      = (/             &
    10131014      ctrl_out((/ 10, 4, 10, 10, 10, 10, 11, 11, 11, 11/),'tke_ter',       &
    … …  
    10501051      ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/),'l_mix_sic',       &
    10511052      "min PBL mixing length "//clnsurf(4),"m", (/ ('', i=1, 10) /)) /)
     1053
    10521054
    10531055  TYPE(ctrl_out), SAVE, DIMENSION(4) :: o_tke_max_srf  = (/                          &
    … …  
    14341436  TYPE(ctrl_out), SAVE :: o_rhum = ctrl_out((/ 2, 5, 10, 10, 10, 10, 11, 11, 11, 11/), &
    14351437    'rhum', 'Relative humidity', '-', (/ ('', i=1, 10) /))
     1438  TYPE(ctrl_out), SAVE :: o_rhl = ctrl_out((/ 2, 6, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1439    'rhl', 'Relative humidity wrt liquid', '%', (/ ('', i=1, 10) /))
     1440  TYPE(ctrl_out), SAVE :: o_rhi = ctrl_out((/ 2, 6, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1441    'rhi', 'Relative humidity wrt ice', '%', (/ ('', i=1, 10) /))
    14361442  TYPE(ctrl_out), SAVE :: o_ozone = ctrl_out((/ 2, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    14371443    'ozone', 'Ozone mole fraction', '-', (/ ('', i=1, 10) /))

  • LMDZ6/trunk/libf/phylmd/phys_output_write_mod.F90

    r3779 r3780  
    133133         o_vitu, o_vitv, o_vitw, o_pres, o_paprs, &
    134134         o_zfull, o_zhalf, o_rneb, o_rnebjn, o_rnebcon, &
    135          o_rnebls, o_rneblsvol, o_rhum, o_ozone, o_ozone_light, &
     135         o_rnebls, o_rneblsvol, o_rhum, o_rhl, o_rhi, o_ozone, o_ozone_light, &
    136136         o_duphy, o_dtphy, o_dqphy, o_dqphy2d, o_dqlphy, o_dqlphy2d, &
    137137         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, &
     138         o_ages_srf, o_snow_srf, o_alb1, o_alb2, o_tke, o_tke_dissip, &
    139139         o_tke_max, o_kz, o_kz_max, o_clwcon, &
    140140         o_dtdyn, o_dqdyn, o_dqdyn2d, o_dqldyn, o_dqldyn2d, &
    … …  
    251251         zt2m_cor,zq2m_cor,zu10m_cor,zv10m_cor, zrh2m_cor, zqsat2m_cor, &
    252252         t2m_min_mon, t2m_max_mon, evap, &
    253          l_mixmin,l_mix, &
     253         l_mixmin,l_mix, tke_dissip, &
    254254         zu10m, zv10m, zq2m, zustar, zxqsurf, &
    255255         rain_lsc, rain_num, snow_lsc, bils, sens, fder, &
    … …  
    297297         ql_seri, qs_seri, tr_seri, &
    298298         zphi, u_seri, v_seri, omega, cldfra, &
    299          rneb, rnebjn, rneblsvol, zx_rh, d_t_dyn,  &
     299         rneb, rnebjn, rneblsvol, zx_rh, zx_rhl, zx_rhi, d_t_dyn,  &
    300300         d_q_dyn,  d_ql_dyn, d_qs_dyn, &
    301301         d_q_dyn2d,  d_ql_dyn2d, d_qs_dyn2d, &
    … …  
    10591059             CALL histwrite_phy(o_l_mixmin(nsrf),  l_mixmin(:,1:klev,nsrf))
    10601060             CALL histwrite_phy(o_tke_max_srf(nsrf),  pbl_tke(:,1:klev,nsrf))
     1061
     1062                         
    10611063          ENDIF
    10621064!jyg<
    … …  
    10691071!            ENDIF
    10701072
    1071 
    10721073       ENDDO
     1074       
     1075               
     1076        IF (iflag_pbl > 1) THEN
     1077          zx_tmp_fi3d=0.
     1078          IF (vars_defined) THEN
     1079             DO nsrf=1,nbsrf
     1080                DO k=1,klev
     1081                   zx_tmp_fi3d(:,k)=zx_tmp_fi3d(:,k) &
     1082                        +pctsrf(:,nsrf)*tke_dissip(:,k,nsrf)
     1083                ENDDO
     1084             ENDDO
     1085          ENDIF
     1086         
     1087          CALL histwrite_phy(o_tke_dissip, zx_tmp_fi3d)   
     1088       ENDIF
    10731089
    10741090       IF (vars_defined) zx_tmp_fi2d(1 : klon) = sens_prec_liq_o(1 : klon, 1)
    … …  
    16951711       CALL histwrite_phy(o_rnebjn, zx_tmp_fi3d)
    16961712       CALL histwrite_phy(o_rhum, zx_rh)
     1713       CALL histwrite_phy(o_rhl, zx_rhl*100.)
     1714       CALL histwrite_phy(o_rhi, zx_rhi*100.)
     1715
    16971716       
    16981717       IF (vars_defined) zx_tmp_fi3d = wo(:, :, 1) * dobson_u * 1e3 / zmasse / rmo3 * rmd
    … …  
    17531772          ENDIF
    17541773          CALL histwrite_phy(o_tke, zx_tmp_fi3d)
    1755 
    1756           CALL histwrite_phy(o_tke_max, zx_tmp_fi3d)
     1774          CALL histwrite_phy(o_tke_max, zx_tmp_fi3d) 
     1775
    17571776       ENDIF
    17581777

  • LMDZ6/trunk/libf/phylmd/phys_state_var_mod.F90

    r3756 r3780  
    449449
    450450include "clesphys.h"
    451 
    452451      IF (is_initialized) RETURN
    453452      is_initialized=.TRUE.

  • LMDZ6/trunk/libf/phylmd/physiq_mod.F90

    r3778 r3780  
    273273       ref_liq, ref_ice, theta,  &
    274274       ref_liq_pi, ref_ice_pi,  &
    275        zphi, zx_rh, &
     275       zphi, zx_rh, zx_rhl, zx_rhi, &
    276276       pmfd, pmfu,  &
    277277       !
    … …  
    37383738          ENDIF
    37393739          zx_rh(i,k) = q_seri(i,k)/zx_qs
     3740          zx_rhl(i,k) = q_seri(i,k)/(qsatl(zx_t)/pplay(i,k))
     3741          zx_rhi(i,k) = q_seri(i,k)/(qsats(zx_t)/pplay(i,k))
    37403742          zqsat(i,k)=zx_qs
    37413743       ENDDO

  • LMDZ6/trunk/libf/phylmd/surf_landice_mod.F90

    r3102 r3780  
    241241!
    242242!****************************************************************************************
     243
     244    ! EV: use calbeta
     245    CALL calbeta(dtime, is_lic, knon, snow, qsol, beta, cal, dif_grnd)
     246
     247
     248    ! use soil model and recalculate properly cal
    243249    IF (soil_model) THEN
    244250       CALL soil(dtime, is_lic, knon, snow, tsurf, tsoil, soilcap, soilflux)
    … …  
    255261!
    256262!****************************************************************************************
    257     beta(:) = 1.0
    258     dif_grnd(:) = 0.0
     263 !   beta(:) = 1.0
     264 !   dif_grnd(:) = 0.0
    259265
    260266! Suppose zero surface speed
    … …  
    291297!****************************************************************************************
    292298    CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) 
     299
    293300    WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
    294301    zfra(1:knon) = MAX(0.0,MIN(1.0,snow(1:knon)/(snow(1:knon)+10.0)))
    … …  
    303310!IM: KstaTER0.77 & LMD_ARMIP6   
    304311
    305 ! Attantion: alb1 and alb2 are the same!
     312! Attantion: alb1 and alb2 are not the same!
    306313    alb1(1:knon)  = alb_vis_sno_lic
    307314    alb2(1:knon)  = alb_nir_sno_lic

  • LMDZ6/trunk/libf/phylmd/yamada4.F90

    r3531 r3780  
    66  USE dimphy
    77  USE ioipsl_getin_p_mod, ONLY : getin_p
    8 
     8  USE phys_local_var_mod, only: tke_dissip
     9 
    910  IMPLICIT NONE
    1011  include "iniprint.h"
    … …  
    5657  !             iflag_pbl=11 -> the model starts with NP from start files created by ce0l
    5758  !                          -> the model can run with longer time-steps.
    58   !             2016/11/30 (EV etienne.vignon@univ-grenoble-alpes.fr)
     59  !             2016/11/30 (EV etienne.vignon@lmd.ipsl.fr)
    5960  !               On met tke (=q2/2) en entr??e plut??t que q2
    6061  !               On corrige l'update de la tke
    61   !
     62  !             2020/10/01 (EV)
     63  !               On ajoute la dissipation de la TKE en diagnostique de sortie
     64  !                 
    6265  ! Inpout/Output :
    6366  !==============
    … …  
    121124  REAL,SAVE :: viscom,viscoh
    122125  !$OMP THREADPRIVATE( hboville,viscom,viscoh)
    123   INTEGER ig, k
     126  INTEGER ig, jg, k
    124127  REAL ri, zrif, zalpha, zsm, zsn
    125128  REAL rif(klon, klev+1), sm(klon, klev+1), alpha(klon, klev)
    … …  
    416419  ELSE IF (iflag_pbl>=10) THEN
    417420
     421    shear(:,:)=0.
     422    buoy(:,:)=0.
     423    dissip(:,:)=0.
     424    km(:,:)=0.
     425       
    418426    IF (yamada4_num>=1) THEN
    419427 
    … …  
    424432      shear(ig,k)=km(ig, k)*m2(ig, k)
    425433      buoy(ig,k)=km(ig, k)*m2(ig, k)*(-1.*rif(ig,k))
    426       dissip(ig,k)=((sqrt(q2(ig,k)))**3)/(b1*l(ig,k))
     434      dissip(ig,k)=min(max(((sqrt(q2(ig,k)))**3)/(b1*l(ig,k)),1.E-12),1.E4)
    427435     ENDDO
    428436    ENDDO
    429 
     437   
    430438    IF (yamada4_num==1) THEN ! Schema du MAR tel quel
    431439       DO k = 2, klev - 1
    … …  
    478486        ENDDO
    479487       ENDDO
     488       
     489      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     490      !! Attention, yamada4_num=5 est inexacte car néglige les termes de flottabilité
     491      !! en conditions instables
     492      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    480493    ELSE IF (yamada4_num==5) THEN ! version modifiee avec integration exacte pour la dissipation
    481494       DO k = 2, klev - 1
    … …  
    507520       DO k = 2, klev - 1
    508521         DO ig=1,ngrid
     522         !tkeprov=q2(ig,k)/ydeux
     523         !tkeprov=tkeprov+max(buoy(ig,k)+shear(ig,k),0.)*dt
     524         !disseff=dissip(ig,k)-min(0.,buoy(ig,k)+shear(ig,k))
     525         !tkeexp=exp(-dt*disseff/tkeprov)
     526         !tkeprov= tkeprov*tkeexp
     527         !q2(ig,k)=tkeprov*ydeux
     528         ! En cas stable, on traite la flotabilite comme la
     529         ! dissipation, en supposant que dissipeff/TKE est constant.
     530         ! Puis on prend la solution exacte
     531         !
     532         ! With drag and dissipation from high vegetation (EV & FC, 05/10/2020)
     533         winds(ig,k)=sqrt(u(ig,k)**2+v(ig,k)**2)
    509534         tkeprov=q2(ig,k)/ydeux
    510          tkeprov=tkeprov+max(buoy(ig,k)+shear(ig,k),0.)*dt
    511          disseff=dissip(ig,k)-min(0.,buoy(ig,k)+shear(ig,k))
     535         tkeprov=tkeprov+max(buoy(ig,k)+shear(ig,k)+drgpro(ig,k)*(winds(ig,k))**3,0.)*dt
     536         disseff=dissip(ig,k)-min(0.,buoy(ig,k)+shear(ig,k)+drgpro(ig,k)*(winds(ig,k))**3) + drgpro(ig,k)*tkeprov
    512537         tkeexp=exp(-dt*disseff/tkeprov)
    513538         tkeprov= tkeprov*tkeexp
    514          q2(ig,k)=tkeprov*ydeux
    515          ! En cas stable, on traite la flotabilite comme la
    516          ! dissipation, en supposant que buoy/q2^3 est constant.
    517          ! Puis on prend la solution exacte
     539         q2(ig,k)=tkeprov*ydeux       
     540         
    518541        ENDDO
    519542       ENDDO
    … …  
    533556!     q2(1:ngrid, k) = q2(1:ngrid, k) + dt*km(1:ngrid, k)*m2(1:ngrid, k)*(1.-rif(1:ngrid,k))
    534557      q2(1:ngrid, k) = min(max(q2(1:ngrid,k),1.E-10), 1.E4)
    535        q2(1:ngrid, k) = 1./(1./sqrt(q2(1:ngrid,k))+dt/(yun*l(1:ngrid,k)*b1))
     558      q2(1:ngrid, k) = 1./(1./sqrt(q2(1:ngrid,k))+dt/(yun*l(1:ngrid,k)*b1))
    536559!     q2(1:ngrid, k) = 1./(1./sqrt(q2(1:ngrid,k))+dt/(2*l(1:ngrid,k)*b1))
    537560      q2(1:ngrid, k) = q2(1:ngrid, k)*q2(1:ngrid, k)
    … …  
    725748
    726749!============================================================================
     750! Diagnostique de la dissipation
     751!============================================================================
     752
     753! Diagnostics
     754 tke_dissip(1:ngrid,:,nsrf)=0.
     755 DO k=2,klev
     756    DO ig=1,ngrid
     757       jg=ni(ig)
     758       tke_dissip(jg,k,nsrf)=dissip(ig,k)
     759    ENDDO
     760 ENDDO
     761 
     762!=============================================================================
    727763
    728764  RETURN
    … …  
    10171053!=====================================================================
    10181054
    1019 
     1055  l1(1:ngrid,:)=0.
    10201056  IF (iflag_pbl==8 .OR. iflag_pbl==10) THEN
    10211057
    … …  
    10691105   l2(1:ngrid,:)=0.0
    10701106   l_mixmin(1:ngrid,:,nsrf)=0.
    1071    l_mix(1:ngrid,:,nsrf)=0.
     1107   l_mix(1:ngrid,:,nsrf)=1.E-5
    10721108
    10731109   IF (nsrf .EQ. 1) THEN
    … …  
    11351171
    11361172
    1137  DO k=2,klev
     1173 DO k=1,klev+1
    11381174    DO ig=1,ngrid
    11391175       lmix(ig,k)=MAX(MAX(l1(ig,k), l2(ig,k)),lmixmin)
     1176       lmix(ig,k)=MAX(lmix(ig,k),1.E-5)
    11401177   ENDDO
    11411178 ENDDO
    … …  
    11431180! Diagnostics
    11441181
    1145  DO k=2,klev
     1182 DO k=2,klev+1
    11461183    DO ig=1,ngrid
    11471184       jg=ni(ig)
Note: See TracChangeset for help on using the changeset viewer.