Ignore:
Timestamp:
Jun 30, 2017, 12:17:42 PM (8 years ago)
Author:
fcheruy
Message:

Update tree branch to trunk version

Location:
LMDZ5/branches/LMDZ_tree_FC
Files:
1 deleted
26 edited
2 copied

Legend:

Unmodified
Added
Removed
  • LMDZ5/branches/LMDZ_tree_FC

  • LMDZ5/branches/LMDZ_tree_FC/libf/dynphy_lonlat/phylmd/etat0phys_netcdf.F90

    r2788 r2925  
    4343    sig1, ftsol, clwcon, fm_therm, wake_Cstar,  pctsrf,  entr_therm,radpas, f0,&
    4444    zmax0,fevap, rnebcon,falb_dir, wake_fip,    agesno,  detr_therm, pbl_tke,  &
    45     phys_state_var_init
     45    phys_state_var_init, ql_ancien, qs_ancien, prlw_ancien, prsw_ancien, &
     46    prw_ancien
    4647  USE comconst_mod, ONLY: pi, dtvr
    4748
     
    201202  t_ancien   = 273.15
    202203  q_ancien   = 0.
     204  ql_ancien = 0.
     205  qs_ancien = 0.
     206  prlw_ancien = 0.
     207  prsw_ancien = 0.
     208  prw_ancien = 0.
    203209  agesno     = 0.
    204210
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/cloudth_mod.F90

    r2686 r2925  
    311311      ! (J Jouhaud, JL Dufresne, JB Madeleine)
    312312      REAL,SAVE :: vert_alpha
     313      !$OMP THREADPRIVATE(vert_alpha)
    313314      LOGICAL, SAVE :: firstcall = .TRUE.
     315      !$OMP THREADPRIVATE(firstcall)
    314316     
    315317      REAL zpdf_sig(ngrid),zpdf_k(ngrid),zpdf_delta(ngrid)
     
    858860      ! Change the width of the PDF used for vertical subgrid scale heterogeneity
    859861      ! (J Jouhaud, JL Dufresne, JB Madeleine)
    860       REAL,SAVE :: vert_alpha
     862      REAL,SAVE :: vert_alpha, vert_alpha_th
     863      !$OMP THREADPRIVATE(vert_alpha, vert_alpha_th)
    861864      LOGICAL, SAVE :: firstcall = .TRUE.
     865      !$OMP THREADPRIVATE(firstcall)
    862866
    863867      REAL zpdf_sig(ngrid),zpdf_k(ngrid),zpdf_delta(ngrid)
     
    893897        CALL getin_p('cloudth_vert_alpha',vert_alpha)
    894898        WRITE(*,*) 'cloudth_vert_alpha = ', vert_alpha
     899        ! The factor used for the thermal is equal to that of the environment
     900        !   if nothing is explicitly specified in the def file
     901        vert_alpha_th=vert_alpha
     902        CALL getin_p('cloudth_vert_alpha_th',vert_alpha_th)
     903        WRITE(*,*) 'cloudth_vert_alpha_th = ', vert_alpha_th
    895904        firstcall=.FALSE.
    896905      ENDIF
     
    9971006      qltot(ind1,ind2)=fraca(ind1,ind2)*qlth(ind1,ind2)+(1.-1.*fraca(ind1,ind2))*qlenv(ind1,ind2)
    9981007
    999       ELSE IF (iflag_cloudth_vert == 3) THEN
     1008      ELSE IF (iflag_cloudth_vert >= 3) THEN
    10001009
    10011010!-------------------------------------------------------------------------------
     
    10061015!      deltasenv=aenv*ratqs(ind1,ind2)*zqsatenv(ind1,ind2)
    10071016!      deltasth=ath*ratqs(ind1,ind2)*zqsatth(ind1,ind2)
    1008       deltasenv=aenv*vert_alpha*sigma1s
    1009       deltasth=ath*vert_alpha*sigma2s
     1017      IF (iflag_cloudth_vert == 3) THEN
     1018        deltasenv=aenv*vert_alpha*sigma1s
     1019        deltasth=ath*vert_alpha_th*sigma2s
     1020      ELSE IF (iflag_cloudth_vert == 4) THEN
     1021        deltasenv=vert_alpha*sigma1s
     1022        deltasth=vert_alpha_th*sigma2s
     1023      ENDIF
    10101024     
    10111025      xenv1=-(senv+deltasenv)/(sqrt(2.)*sigma1s)
     
    10491063
    10501064
    1051       ENDIF ! of if (iflag_cloudth_vert==1 or 3)
     1065      ENDIF ! of if (iflag_cloudth_vert==1 or 3 or 4)
    10521066
    10531067      if (cenv(ind1,ind2).lt.1.e-10.or.cth(ind1,ind2).lt.1.e-10) then
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/comsoil.h

    r1907 r2925  
    33!
    44
    5       common /comsoil/inertie_sol,inertie_sno,inertie_ice
    6       real inertie_sol,inertie_sno,inertie_ice
     5      common /comsoil/inertie_sol,inertie_sno,inertie_sic,inertie_lic,  &
     6     &                iflag_sic
     7      real inertie_sol,inertie_sno,inertie_sic,inertie_lic
     8      integer iflag_sic
    79!$OMP THREADPRIVATE(/comsoil/)
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/conf_phys_m.F90

    r2924 r2925  
    173173    REAL,SAVE :: exposant_glace_omp
    174174    REAL,SAVE :: rei_min_omp, rei_max_omp
    175     REAL,SAVE :: inertie_sol_omp,inertie_sno_omp,inertie_ice_omp
     175    INTEGER,SAVE :: iflag_sic_omp
     176    REAL,SAVE :: inertie_sol_omp,inertie_sno_omp,inertie_sic_omp
     177    REAL,SAVE :: inertie_lic_omp
    176178    REAL,SAVE :: qsol0_omp
    177179    REAL,SAVE :: evap0_omp
     
    11371139    !-----------------------------------------------------------------------
    11381140    !
    1139     !Config Key  = inertie_ice
     1141    !Config Key  = iflag_sic
     1142    !Config Desc = 
     1143    !Config Def  = 0
     1144    !Config Help =
     1145    !
     1146    iflag_sic_omp = 0
     1147    CALL getin('iflag_sic',iflag_sic_omp)
     1148    !
     1149    !Config Key  = inertie_sic
    11401150    !Config Desc = 
    11411151    !Config Def  = 2000.
    11421152    !Config Help =
    11431153    !
    1144     inertie_ice_omp = 2000.
    1145     CALL getin('inertie_ice',inertie_ice_omp)
     1154    inertie_sic_omp = 2000.
     1155    CALL getin('inertie_sic',inertie_sic_omp)
     1156    !
     1157    !Config Key  = inertie_lic
     1158    !Config Desc = 
     1159    !Config Def  = 2000.
     1160    !Config Help =
     1161    !
     1162    inertie_lic_omp = 2000.
     1163    CALL getin('inertie_lic',inertie_lic_omp)
    11461164    !
    11471165    !Config Key  = inertie_sno
     
    21542172    evap0 = evap0_omp
    21552173    albsno0 = albsno0_omp
     2174    iflag_sic = iflag_sic_omp
    21562175    inertie_sol = inertie_sol_omp
    2157     inertie_ice = inertie_ice_omp
     2176    inertie_sic = inertie_sic_omp
     2177    inertie_lic = inertie_lic_omp
    21582178    inertie_sno = inertie_sno_omp
    21592179    rad_froid = rad_froid_omp
     
    25722592    write(lunout,*)' evap0 = ', evap0
    25732593    write(lunout,*)' albsno0 = ', albsno0
     2594    write(lunout,*)' iflag_sic = ', iflag_sic
    25742595    write(lunout,*)' inertie_sol = ', inertie_sol
    2575     write(lunout,*)' inertie_ice = ', inertie_ice
     2596    write(lunout,*)' inertie_sic = ', inertie_sic
     2597    write(lunout,*)' inertie_lic = ', inertie_lic
    25762598    write(lunout,*)' inertie_sno = ', inertie_sno
    25772599    write(lunout,*)' f_cdrag_ter = ',f_cdrag_ter
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/cv3_routines.F90

    r2818 r2925  
    111111     ok_optim_yield=.False.
    112112     CALL getin_p('ok_optim_yield',ok_optim_yield)
     113     ok_homo_tend=.TRUE.
     114     CALL getin_p('ok_homo_tend',ok_homo_tend)
     115     ok_entrain=.TRUE.
     116     CALL getin_p('ok_entrain',ok_entrain)
     117
    113118     coef_peel=0.25
    114119     CALL getin_p('coef_peel',coef_peel)
     
    277282
    278283SUBROUTINE cv3_feed(len, nd, ok_conserv_q, &
    279                     t, q, u, v, p, ph, hm, gz, &
     284                    t, q, u, v, p, ph, h, gz, &
    280285                    p1feed, p2feed, wght, &
    281286                    wghti, tnk, thnk, qnk, qsnk, unk, vnk, &
     
    283288
    284289  USE mod_phys_lmdz_transfert_para, ONLY : bcast
     290  USE add_phys_tend_mod, ONLY: fl_cor_ebil
    285291  IMPLICIT NONE
    286292
     
    292298! - here, nk(i)=minorig
    293299! - icb defined differently (plcl compared with ph instead of p)
     300! - dry static energy as argument instead of moist static energy
    294301
    295302! Main differences with convect3:
     
    307314  REAL, DIMENSION (len, nd), INTENT (IN)             :: t, q, p
    308315  REAL, DIMENSION (len, nd), INTENT (IN)             :: u, v
    309   REAL, DIMENSION (len, nd), INTENT (IN)             :: hm, gz
     316  REAL, DIMENSION (len, nd), INTENT (IN)             :: h, gz
    310317  REAL, DIMENSION (len, nd+1), INTENT (IN)           :: ph
    311318  REAL, DIMENSION (len), INTENT (IN)                 :: p1feed
     
    378385    pup(i) = p2feed(i)
    379386  END DO
    380   CALL cv3_vertmix(len, nd, iflag, p1feed, pup, p, ph, &
    381                    t, q, u, v, wght, &
    382                    wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plclup)
     387  IF (fl_cor_ebil >=2 ) THEN
     388    CALL cv3_estatmix(len, nd, iflag, p1feed, pup, p, ph, &
     389                     t, q, u, v, h, gz, wght, &
     390                     wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plclup)
     391  ELSE
     392    CALL cv3_enthalpmix(len, nd, iflag, p1feed, pup, p, ph, &
     393                       t, q, u, v, wght, &
     394                       wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plclup)
     395  ENDIF  ! (fl_cor_ebil >=2 )
    383396! 1.b- LCL associated with ph(nk+1)
    384397  DO i = 1, len
    385398    plo(i) = ph(i, nk(i)+1)
    386399  END DO
    387   CALL cv3_vertmix(len, nd, iflag, p1feed, plo, p, ph, &
    388                    t, q, u, v, wght, &
    389                    wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plcllo)
     400  IF (fl_cor_ebil >=2 ) THEN
     401    CALL cv3_estatmix(len, nd, iflag, p1feed, plo, p, ph, &
     402                     t, q, u, v, h, gz, wght, &
     403                     wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plcllo)
     404  ELSE
     405    CALL cv3_enthalpmix(len, nd, iflag, p1feed, plo, p, ph, &
     406                       t, q, u, v, wght, &
     407                       wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plcllo)
     408  ENDIF  ! (fl_cor_ebil >=2 )
    390409! 2- Iterations
    391410  niter = 5
     
    434453!jyg>
    435454
    436     CALL cv3_vertmix(len, nd, iflag, p1feed, pfeed, p, ph, &
    437                    t, q, u, v, wght, &
    438                    wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plclfeed)
     455    IF (fl_cor_ebil >=2 ) THEN
     456      CALL cv3_estatmix(len, nd, iflag, p1feed, pfeed, p, ph, &
     457                       t, q, u, v, h, gz, wght, &
     458                       wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plclfeed)
     459    ELSE
     460      CALL cv3_enthalpmix(len, nd, iflag, p1feed, pfeed, p, ph, &
     461                         t, q, u, v, wght, &
     462                         wghti, nk, tnk, thnk, qnk, qsnk, unk, vnk, plclfeed)
     463    ENDIF  ! (fl_cor_ebil >=2 )
    439464!jyg20140217<
    440465    IF (ok_new_feed) THEN
     
    16471672      DO i = 1, ncum
    16481673        IF ((k>=icb(i)) .AND. (k<=inb(i))) THEN
     1674!jyg<   (energy conservation tests)
     1675!!          hp(i, k) = hnk(i) + (lv(i,k)+(cpd-cpv)*tp(i,k))*ep(i, k)*clw(i, k)
     1676!!          hp(i, k) = ( hnk(i) + (lv(i,k)+(cpd-cpv)*t(i,k))*ep(i, k)*clw(i, k) ) / &
     1677!!                     (1. - ep(i,k)*clw(i,k))
     1678!!          hp(i, k) = ( hnk(i) + (lv(i,k)+(cpd-cl)*t(i,k))*ep(i, k)*clw(i, k) ) / &
     1679!!                     (1. - ep(i,k)*clw(i,k))
    16491680          hp(i, k) = hnk(i) + (lv(i,k)+(cpd-cpv)*t(i,k))*ep(i, k)*clw(i, k)
    16501681        END IF
     
    29843015                     qcondc, wd, &
    29853016                     ftd, fqd, qnk, qtc, sigt, tau_cld_cv, coefw_cld_cv)
     3017
     3018    USE print_control_mod, ONLY: lunout, prt_level
     3019    USE add_phys_tend_mod, only : fl_cor_ebil
    29863020
    29873021  IMPLICIT NONE
     
    32203254
    32213255      IF ((0.01*grav*work(il)*am(il))>=delti) iflag(il) = 1 !consist vect
    3222       ft(il, 1) = ft(il, 1) + 0.01*grav*work(il)*am(il) * &
    3223                                    (t(il,2)-t(il,1)+(gz(il,2)-gz(il,1))/cpn(il,1))
     3256!jyg<
     3257        IF (fl_cor_ebil >= 2) THEN
     3258          ft(il, 1) = ft(il, 1) + 0.01*grav*work(il)*am(il) * &
     3259                    ((t(il,2)-t(il,1))*cpn(il,2)+gz(il,2)-gz(il,1))/cpn(il,1)
     3260        ELSE
     3261          ft(il, 1) = ft(il, 1) + 0.01*grav*work(il)*am(il) * &
     3262                    (t(il,2)-t(il,1)+(gz(il,2)-gz(il,1))/cpn(il,1))
     3263        ENDIF
     3264!>jyg
    32243265    END IF ! iflag
    32253266  END DO
     
    35293570
    35303571! sature
    3531         ft(il, i) = ft(il, i) + 0.01*grav*dpinv * &
     3572!jyg<
     3573        IF (fl_cor_ebil >= 2) THEN
     3574          ft(il, i) = ft(il, i) + 0.01*grav*dpinv * &
     3575              ( amp1(il)*( (t(il,i+1)-t(il,i))*cpn(il,i+1) + gz(il,i+1)-gz(il,i))*cpinv - &
     3576                ad(il)*( (t(il,i)-t(il,i-1))*cpn(il,i-1) + gz(il,i)-gz(il,i-1))*cpinv)
     3577        ELSE
     3578          ft(il, i) = ft(il, i) + 0.01*grav*dpinv * &
    35323579                     (amp1(il)*(t(il,i+1)-t(il,i) + (gz(il,i+1)-gz(il,i))*cpinv) - &
    35333580                      ad(il)*(t(il,i)-t(il,i-1)+(gz(il,i)-gz(il,i-1))*cpinv))
     3581        ENDIF
     3582!>jyg
    35343583
    35353584
     
    35383587                                    t(il,i)*(cpv-cpd)*(rr(il,i)-qent(il,i,i)))*cpinv
    35393588        END IF
    3540 
    3541 
    3542 
     3589!
    35433590! sb: on ne fait pas encore la correction permettant de mieux
    35443591! conserver l'eau:
     
    38733920
    38743921!!!!      do 700 i=1,icb(il)-1
    3875   DO i = 1, nl
    3876     DO il = 1, ncum
    3877       IF (i<=(icb(il)-1) .AND. iflag(il)<=1) THEN
    3878         ftd(il, i) = esum(il)*t_wake(il, i)/(th_wake(il,i)*hsum(il))
    3879         fqd(il, i) = fsum(il)/gsum(il)
    3880         ft(il, i) = ftd(il, i) + asum(il)*t(il, i)/(th(il,i)*dsum(il))
    3881         fr(il, i) = fqd(il, i) + bsum(il)/csum(il)
    3882       END IF
    3883     END DO
    3884   END DO
     3922  IF (ok_homo_tend) THEN
     3923    DO i = 1, nl
     3924      DO il = 1, ncum
     3925        IF (i<=(icb(il)-1) .AND. iflag(il)<=1) THEN
     3926          ftd(il, i) = esum(il)*t_wake(il, i)/(th_wake(il,i)*hsum(il))
     3927          fqd(il, i) = fsum(il)/gsum(il)
     3928          ft(il, i) = ftd(il, i) + asum(il)*t(il, i)/(th(il,i)*dsum(il))
     3929          fr(il, i) = fqd(il, i) + bsum(il)/csum(il)
     3930        END IF
     3931      END DO
     3932    END DO
     3933  ENDIF
    38853934
    38863935!jyg<
     
    39203969  END DO
    39213970!
    3922 !    print *,' YIELD : alpha_qpos ',alpha_qpos(1)
     3971    IF (prt_level .GE. 5) THEN
     3972      print *,' CV3_YIELD : alpha_qpos ',alpha_qpos(1)
     3973    ENDIF
    39233974!
    39243975  DO il = 1, ncum
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/cv3p_mixing.F90

    r2478 r2925  
    1414
    1515  USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
     16  USE ioipsl_getin_p_mod, ONLY: getin_p
     17  USE add_phys_tend_mod, ONLY: fl_cor_ebil
    1618
    1719  IMPLICIT NONE
     
    110112            fmax, gammas, qqa1, qqa2, Qcoef1max, Qcoef2max
    111113!>jyg
    112 
     114!
    113115  END IF
    114116
     
    165167  DO i = minorig + 1, nl
    166168
    167     DO j = minorig, nl
    168       DO il = 1, ncum
    169         IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. (j>=(icb(il)-1)) &
    170                          .AND. (j<=inb(il))) THEN
    171 
    172           rti = qnk(il) - ep(il, i)*clw(il, i)
    173           bf2 = 1. + lv(il, j)*lv(il, j)*rs(il, j)/(rrv*t(il,j)*t(il,j)*cpd)
    174 !jyg(from aj)<
    175           IF (cvflag_ice) THEN
    176 ! print*,cvflag_ice,'cvflag_ice dans do 700'
    177             IF (t(il,j)<=263.15) THEN
    178               bf2 = 1. + (lf(il,j)+lv(il,j))*(lv(il,j)+frac(il,j)* &
    179                    lf(il,j))*rs(il, j)/(rrv*t(il,j)*t(il,j)*cpd)
    180             END IF
    181           END IF
    182 !>jyg
    183           anum = h(il, j) - hp(il, i) + (cpv-cpd)*t(il, j)*(rti-rr(il,j))
    184           denom = h(il, i) - hp(il, i) + (cpd-cpv)*(rr(il,i)-rti)*t(il, j)
    185           dei = denom
    186           IF (abs(dei)<0.01) dei = 0.01
    187           Sij(il, i, j) = anum/dei
    188           Sij(il, i, i) = 1.0
    189           altem = Sij(il, i, j)*rr(il, i) + (1.-Sij(il,i,j))*rti - rs(il, j)
    190           altem = altem/bf2
    191           cwat = clw(il, j)*(1.-ep(il,j))
    192           stemp = Sij(il, i, j)
    193           IF ((stemp<0.0 .OR. stemp>1.0 .OR. altem>cwat) .AND. j>i) THEN
     169    IF (ok_entrain) THEN
     170      DO j = minorig, nl
     171        DO il = 1, ncum
     172          IF ((i>=icb(il)) .AND. (i<=inb(il)) .AND. (j>=(icb(il)-1)) &
     173                           .AND. (j<=inb(il))) THEN
     174
     175            rti = qnk(il) - ep(il, i)*clw(il, i)
     176            bf2 = 1. + lv(il, j)*lv(il, j)*rs(il, j)/(rrv*t(il,j)*t(il,j)*cpd)
    194177!jyg(from aj)<
    195178            IF (cvflag_ice) THEN
    196               anum = anum - (lv(il,j)+frac(il,j)*lf(il,j))*(rti-rs(il,j)-cwat*bf2)
    197               denom = denom + (lv(il,j)+frac(il,j)*lf(il,j))*(rr(il,i)-rti)
    198             ELSE
    199               anum = anum - lv(il, j)*(rti-rs(il,j)-cwat*bf2)
    200               denom = denom + lv(il, j)*(rr(il,i)-rti)
     179! print*,cvflag_ice,'cvflag_ice dans do 700'
     180              IF (t(il,j)<=263.15) THEN
     181                bf2 = 1. + (lf(il,j)+lv(il,j))*(lv(il,j)+frac(il,j)* &
     182                     lf(il,j))*rs(il, j)/(rrv*t(il,j)*t(il,j)*cpd)
     183              END IF
    201184            END IF
    202185!>jyg
    203             IF (abs(denom)<0.01) denom = 0.01
    204             Sij(il, i, j) = anum/denom
     186            anum = h(il, j) - hp(il, i) + (cpv-cpd)*t(il, j)*(rti-rr(il,j))
     187            denom = h(il, i) - hp(il, i) + (cpd-cpv)*(rr(il,i)-rti)*t(il, j)
     188            dei = denom
     189            IF (abs(dei)<0.01) dei = 0.01
     190            Sij(il, i, j) = anum/dei
     191            Sij(il, i, i) = 1.0
    205192            altem = Sij(il, i, j)*rr(il, i) + (1.-Sij(il,i,j))*rti - rs(il, j)
    206             altem = altem - (bf2-1.)*cwat
    207           END IF
    208           IF (Sij(il,i,j)>0.0) THEN
     193            altem = altem/bf2
     194            cwat = clw(il, j)*(1.-ep(il,j))
     195            stemp = Sij(il, i, j)
     196            IF ((stemp<0.0 .OR. stemp>1.0 .OR. altem>cwat) .AND. j>i) THEN
     197!jyg(from aj)<
     198              IF (cvflag_ice) THEN
     199                anum = anum - (lv(il,j)+frac(il,j)*lf(il,j))*(rti-rs(il,j)-cwat*bf2)
     200                denom = denom + (lv(il,j)+frac(il,j)*lf(il,j))*(rr(il,i)-rti)
     201              ELSE
     202                anum = anum - lv(il, j)*(rti-rs(il,j)-cwat*bf2)
     203                denom = denom + lv(il, j)*(rr(il,i)-rti)
     204              END IF
     205!>jyg
     206              IF (abs(denom)<0.01) denom = 0.01
     207              Sij(il, i, j) = anum/denom
     208              altem = Sij(il, i, j)*rr(il, i) + (1.-Sij(il,i,j))*rti - rs(il, j)
     209              altem = altem - (bf2-1.)*cwat
     210            END IF
     211            IF (Sij(il,i,j)>0.0) THEN
    209212!!!                 Ment(il,i,j)=m(il,i)
    210             Ment(il, i, j) = 1.
    211             elij(il, i, j) = altem
    212             elij(il, i, j) = amax1(0.0, elij(il,i,j))
    213             nent(il, i) = nent(il, i) + 1
    214           END IF
    215 
    216           Sij(il, i, j) = amax1(0.0, Sij(il,i,j))
    217           Sij(il, i, j) = amin1(1.0, Sij(il,i,j))
    218         END IF ! new
    219       END DO
    220     END DO
     213              Ment(il, i, j) = 1.
     214              elij(il, i, j) = altem
     215              elij(il, i, j) = amax1(0.0, elij(il,i,j))
     216              nent(il, i) = nent(il, i) + 1
     217            END IF
     218
     219            Sij(il, i, j) = amax1(0.0, Sij(il,i,j))
     220            Sij(il, i, j) = amin1(1.0, Sij(il,i,j))
     221          END IF ! new
     222        END DO
     223      END DO
     224    ELSE  ! (ok_entrain)
     225      DO il = 1,ncum
     226        nent(il,i) = 0
     227      ENDDO
     228    ENDIF ! (ok_entrain)
    221229
    222230!jygdebug<
     
    243251        uent(il, i, i) = unk(il)
    244252        vent(il, i, i) = vnk(il)
     253        IF (fl_cor_ebil .GE. 2) THEN
     254          hent(il, i, i) = hp(il,i)
     255        ENDIF
    245256        elij(il, i, i) = clw(il, i)*(1.-ep(il,i))
    246257        Sij(il, i, i) = 0.0
    247258      END IF
    248259    END DO
    249   END DO
     260  END DO ! i = minorig + 1, nl
    250261
    251262!jyg!  DO j = 1, ntra
     
    652663        vent(il, i, i) = vnk(il)
    653664        elij(il, i, i) = clw(il, i)*(1.-ep(il,i))
    654         Sij(il, i, i) = 0.0
     665        IF (fl_cor_ebil .GE. 2) THEN
     666          hent(il, i, i) = hp(il,i)
     667          Sigij(il, i, i) = 0.0
     668        ELSE
     669          Sij(il, i, i) = 0.0
     670        ENDIF
    655671      END IF
    656672    END DO ! il
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/cv3param.h

    r2761 r2925  
    77!------------------------------------------------------------
    88
     9      logical ok_homo_tend
    910      logical ok_optim_yield
     11      logical ok_entrain
    1012      logical ok_convstop
    1113      logical ok_intermittent
     
    4143                      ,noff, minorig, nl, nlp, nlm  &
    4244                      ,ok_convstop, ok_intermittent &
    43                       ,ok_optim_yield
     45                      ,ok_optim_yield &
     46                      ,ok_entrain &
     47                      ,ok_homo_tend
    4448!$OMP THREADPRIVATE(/cv3param/)
    4549
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/cva_driver.F90

    r2853 r2925  
    4040
    4141  USE print_control_mod, ONLY: prt_level, lunout
     42  USE add_phys_tend_mod, ONLY: fl_cor_ebil
    4243  IMPLICIT NONE
    4344
     
    730731             PRINT *, 'cva_driver -> cv3_feed'
    731732    CALL cv3_feed(len, nd, ok_conserv_q, &                 ! nd->na
    732                   t1, q1, u1, v1, p1, ph1, hm1, gz1, &
     733                  t1, q1, u1, v1, p1, ph1, h1, gz1, &
    733734                  p1feed1, p2feed1, wght1, &
    734735                  wghti1, tnk1, thnk1, qnk1, qsnk1, unk1, vnk1, &
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/dyn1d/1DUTILS.h

    r2857 r2925  
    482482       forc_omega =0
    483483       CALL getin('forc_omega',forc_omega)
     484
     485!Config  Key  = forc_u
     486!Config  Desc = forcage ou non par u
     487!Config  Def  = false
     488!Config  Help = forcage ou non par u
     489       forc_u =0
     490       CALL getin('forc_u',forc_u)
     491
     492!Config  Key  = forc_v
     493!Config  Desc = forcage ou non par v
     494!Config  Def  = false
     495!Config  Help = forcage ou non par v
     496       forc_v =0
     497       CALL getin('forc_v',forc_v)
    484498
    485499!Config  Key  = forc_w
     
    653667
    654668      modname = 'dyn1deta0 : '
    655 !!      nmq(1)="vap"
    656 !!      nmq(2)="cond"
    657 !!      do iq=3,nqtot
    658 !!        write(nmq(iq),'("tra",i1)') iq-2
    659 !!      enddo
    660       DO iq = 1,nqtot
    661         nmq(iq) = trim(tname(iq))
    662       ENDDO
     669      nmq(1)="vap"
     670      nmq(2)="cond"
     671      do iq=3,nqtot
     672        write(nmq(iq),'("tra",i1)') iq-2
     673      enddo
    663674      print*,'in dyn1deta0 ',fichnom,klon,klev,nqtot
    664675      CALL open_startphy(fichnom)
     
    807818      CALL open_restartphy(fichnom)
    808819      print*,'redm1 ',fichnom,klon,klev,nqtot
    809 !!      nmq(1)="vap"
    810 !!      nmq(2)="cond"
    811 !!      nmq(3)="tra1"
    812 !!      nmq(4)="tra2"
    813       DO iq = 1,nqtot
    814         nmq(iq) = trim(tname(iq))
    815       ENDDO
     820      nmq(1)="vap"
     821      nmq(2)="cond"
     822      nmq(3)="tra1"
     823      nmq(4)="tra2"
    816824
    817825      modname = 'dyn1dredem'
     
    13961404      cor(:) = rkappa*temp*(1.+q(:,1)*rv/rd)/(play*(1.+q(:,1)))
    13971405
    1398 
    13991406      do k=2,llm-1
    14001407
     
    27892796         hq_mod_cas(l)= hq_prof_cas(k2) - frac*(hq_prof_cas(k2)-hq_prof_cas(k1))
    27902797         vq_mod_cas(l)= vq_prof_cas(k2) - frac*(vq_prof_cas(k2)-vq_prof_cas(k1))
    2791          dtrad_mod_cas(l)= dtrad_prof_cas(k2) - frac*(dtrad_prof_cas(k2)-dtrad_prof_cas(k1))
    27922798     
    27932799         else !play>plev_prof_cas(1)
     
    28162822         hq_mod_cas(l)= frac1*hq_prof_cas(k1) - frac2*hq_prof_cas(k2)
    28172823         vq_mod_cas(l)= frac1*vq_prof_cas(k1) - frac2*vq_prof_cas(k2)
    2818          dtrad_mod_cas(l)= frac1*dtrad_prof_cas(k1) - frac2*dtrad_prof_cas(k2)
    28192824
    28202825         endif ! play.le.plev_prof_cas(1)
     
    28452850         hq_mod_cas(l)= hq_prof_cas(nlev_cas)*fact                            !jyg
    28462851         vq_mod_cas(l)= vq_prof_cas(nlev_cas)*fact                            !jyg
    2847          dtrad_mod_cas(l)= dtrad_prof_cas(nlev_cas)*fact                      !jyg
    28482852 
    28492853        endif ! play
     
    51715175         hq_mod_cas(l)= hq_prof_cas(k2) - frac*(hq_prof_cas(k2)-hq_prof_cas(k1))
    51725176         vq_mod_cas(l)= vq_prof_cas(k2) - frac*(vq_prof_cas(k2)-vq_prof_cas(k1))
    5173          dtrad_mod_cas(l)= dtrad_prof_cas(k2) - frac*(dtrad_prof_cas(k2)-dtrad_prof_cas(k1))
    51745177     
    51755178         else !play>plev_prof_cas(1)
     
    52085211         hq_mod_cas(l)= frac1*hq_prof_cas(k1) - frac2*hq_prof_cas(k2)
    52095212         vq_mod_cas(l)= frac1*vq_prof_cas(k1) - frac2*vq_prof_cas(k2)
    5210          dtrad_mod_cas(l)= frac1*dtrad_prof_cas(k1) - frac2*dtrad_prof_cas(k2)
    52115213
    52125214         endif ! play.le.plev_prof_cas(1)
     
    52455247         hq_mod_cas(l)= hq_prof_cas(nlev_cas)*fact                     !jyg
    52465248         vq_mod_cas(l)= vq_prof_cas(nlev_cas)*fact                     !jyg
    5247          dtrad_mod_cas(l)= dtrad_prof_cas(nlev_cas)*fact               !jyg
    52485249 
    52495250        endif ! play
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/dyn1d/1D_interp_cases.h

    r2716 r2925  
    3737
    3838       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    39        d_th_adv(l) = ht_gcssold(l)
     39       d_t_adv(l) = ht_gcssold(l)
    4040       d_q_adv(l,1) = hq_gcssold(l)
    4141       dt_cooling(l) = 0.0
     
    8484
    8585       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    86        d_th_adv(l) = alpha*omega(l)/rcpd-(ht_mod(l)+vt_mod(l))
     86       d_t_adv(l) = alpha*omega(l)/rcpd-(ht_mod(l)+vt_mod(l))
    8787       d_q_adv(l,1) = -(hq_mod(l)+vq_mod(l))
    8888       dt_cooling(l) = 0.0
     
    183183
    184184       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    185        d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-d_t_dyn_z(l)
     185       d_t_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-d_t_dyn_z(l)
    186186       d_q_adv(l,1) = hq_mod(l)-d_q_dyn_z(l)
    187187       d_u_adv(l) = hu_mod(l)-d_u_dyn_z(l)
     
    224224       ug(l)= ug_mod(l)
    225225       vg(l)= vg_mod(l)
    226        d_th_adv(l)=ht_mod(l)
     226       d_t_adv(l)=ht_mod(l)
    227227       d_q_adv(l,1)=hq_mod(l)
    228228      enddo
     
    315315!calcul de l'advection totale
    316316        if (cptadvw) then
    317         d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-d_t_dyn_z(l)
     317        d_t_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-d_t_dyn_z(l)
    318318!        print*,'temp vert adv',l,ht_mod(l),vt_mod(l),-d_t_dyn_z(l)
    319319        d_q_adv(l,1) = hq_mod(l)-d_q_dyn_z(l)
    320320!        print*,'q vert adv',l,hq_mod(l),vq_mod(l),-d_q_dyn_z(l)
    321321        else
    322         d_th_adv(l) = alpha*omega(l)/rcpd+(ht_mod(l)+vt_mod(l))
     322        d_t_adv(l) = alpha*omega(l)/rcpd+(ht_mod(l)+vt_mod(l))
    323323        d_q_adv(l,1) = (hq_mod(l)+vq_mod(l))
    324324        endif
     
    392392       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    393393!calcul de l'advection totale
    394 !        d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-omega(l)*d_t_z(l)
     394!        d_t_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-omega(l)*d_t_z(l)
    395395!attention: on impose dth
    396         d_th_adv(l) = alpha*omega(l)/rcpd+                                  &
     396        d_t_adv(l) = alpha*omega(l)/rcpd+                                  &
    397397     &         ht_mod(l)*(play(l)/pzero)**rkappa-omega(l)*d_t_z(l)
    398 !        d_th_adv(l) = 0.
     398!        d_t_adv(l) = 0.
    399399!        print*,'temp vert adv',l,ht_mod(l),vt_mod(l),-d_t_dyn_z(l)
    400400        d_q_adv(l,1) = hq_mod(l)-omega(l)*d_q_z(l)
     
    417417!---------------------------------------------------------------------
    418418      if (forcing_rico) then
    419 !       call lstendH(llm,omega,dt_dyn,dq_dyn,du_dyn, dv_dyn,
    420 !     :  q,temp,u,v,play)
     419!      call lstendH(llm,omega,dt_dyn,dq_dyn,du_dyn, dv_dyn,q,temp,u,v,play)
    421420       call lstendH(llm,nqtot,omega,dt_dyn,dq_dyn,q,temp,u,v,play)
    422421
    423422        do l=1,llm
    424        d_th_adv(l) =  (dth_rico(l) +  dt_dyn(l))
     423       d_t_adv(l) =  (dth_rico(l) +  dt_dyn(l))
    425424       d_q_adv(l,1) = (dqh_rico(l) +  dq_dyn(l,1))
    426425       d_q_adv(l,2) = 0.
     
    458457       vg(l)= v_mod(l)
    459458       IF((phi(l)/RG).LT.1000) THEN
    460          d_th_adv(l) = (adv_theta_prof + rad_theta_prof)/3600.
     459         d_t_adv(l) = (adv_theta_prof + rad_theta_prof)/3600.
    461460         d_q_adv(l,1) = adv_qt_prof/1000./3600.
    462461         d_q_adv(l,2) = 0.0
    463462!        print *,'INF1000: phi dth dq1 dq2',
    464 !    :  phi(l)/RG,d_th_adv(l),d_q_adv(l,1),d_q_adv(l,2)
     463!    :  phi(l)/RG,d_t_adv(l),d_q_adv(l,1),d_q_adv(l,2)
    465464       ELSEIF ((phi(l)/RG).GE.1000.AND.(phi(l)/RG).lt.3000) THEN
    466465         fact=((phi(l)/RG)-1000.)/2000.
    467466         fact=1-fact
    468          d_th_adv(l) = (adv_theta_prof + rad_theta_prof)*fact/3600.
     467         d_t_adv(l) = (adv_theta_prof + rad_theta_prof)*fact/3600.
    469468         d_q_adv(l,1) = adv_qt_prof*fact/1000./3600.
    470469         d_q_adv(l,2) = 0.0
    471470!        print *,'SUP1000: phi fact dth dq1 dq2',
    472 !    :  phi(l)/RG,fact,d_th_adv(l),d_q_adv(l,1),d_q_adv(l,2)
     471!    :  phi(l)/RG,fact,d_t_adv(l),d_q_adv(l,1),d_q_adv(l,2)
    473472       ELSE
    474          d_th_adv(l) = 0.0
     473         d_t_adv(l) = 0.0
    475474         d_q_adv(l,1) = 0.0
    476475         d_q_adv(l,2) = 0.0
    477476!        print *,'SUP3000: phi dth dq1 dq2',
    478 !    :  phi(l)/RG,d_th_adv(l),d_q_adv(l,1),d_q_adv(l,2)
     477!    :  phi(l)/RG,d_t_adv(l),d_q_adv(l,1),d_q_adv(l,2)
    479478       ENDIF
    480479      dt_cooling(l) = 0.0 
    481480!     print *,'Interp armcu: phi dth dq1 dq2',
    482 !    :  l,phi(l),d_th_adv(l),d_q_adv(l,1),d_q_adv(l,2)
     481!    :  l,phi(l),d_t_adv(l),d_q_adv(l,1),d_q_adv(l,2)
    483482      enddo
    484483      endif ! forcing_armcu
     
    554553       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    555554!
    556 !      d_th_adv(l) = 0.0
     555!      d_t_adv(l) = 0.0
    557556!      d_q_adv(l,1) = 0.0
    558557!CR:test advection=0
    559558!calcul de l'advection verticale
    560         d_th_adv(l) = alpha*omega(l)/rcpd-d_t_dyn_z(l)
     559        d_t_adv(l) = alpha*omega(l)/rcpd-d_t_dyn_z(l)
    561560!        print*,'temp adv',l,-d_t_dyn_z(l)
    562561        d_q_adv(l,1) = -d_q_dyn_z(l)
     
    637636       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    638637!
    639 !      d_th_adv(l) = 0.0
     638!      d_t_adv(l) = 0.0
    640639!      d_q_adv(l,1) = 0.0
    641640!CR:test advection=0
    642641!calcul de l'advection verticale
    643         d_th_adv(l) = alpha*omega(l)/rcpd-d_t_dyn_z(l)
     642        d_t_adv(l) = alpha*omega(l)/rcpd-d_t_dyn_z(l)
    644643!        print*,'temp adv',l,-d_t_dyn_z(l)
    645644        d_q_adv(l,1) = -d_q_dyn_z(l)
     
    715714
    716715!Calcul de l advection verticale
     716
    717717      d_t_dyn_z(:)=w_mod_cas(:)*d_t_z(:)
     718
    718719      d_q_dyn_z(:)=w_mod_cas(:)*d_q_z(:)
    719720      d_u_dyn_z(:)=w_mod_cas(:)*d_u_z(:)
     
    764765
    765766      do l = 1, llm
    766        omega(l) = w_mod_cas(l)
     767       omega(l) = w_mod_cas(l)  ! juste car w_mod_cas en Pa/s (MPL 20170310)
    767768       omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
    768769       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
     
    782783
    783784        if ((tend_t.eq.1).and.(tend_w.eq.0)) then
    784 !           d_th_adv(l)=alpha*omega(l)/rcpd+dt_mod_cas(l)
    785            d_th_adv(l)=alpha*omega(l)/rcpd-dt_mod_cas(l)
     785!           d_t_adv(l)=alpha*omega(l)/rcpd+dt_mod_cas(l)
     786           d_t_adv(l)=alpha*omega(l)/rcpd-dt_mod_cas(l)
    786787        else if ((tend_t.eq.1).and.(tend_w.eq.1)) then
    787 !           d_th_adv(l)=alpha*omega(l)/rcpd+ht_mod_cas(l)-d_t_dyn_z(l)
    788            d_th_adv(l)=alpha*omega(l)/rcpd-ht_mod_cas(l)-d_t_dyn_z(l)
     788!           d_t_adv(l)=alpha*omega(l)/rcpd+ht_mod_cas(l)-d_t_dyn_z(l)
     789           d_t_adv(l)=alpha*omega(l)/rcpd-ht_mod_cas(l)-d_t_dyn_z(l)
    789790        endif
    790791
     
    816817      ENDIF
    817818      endif ! forcing_case
    818 
    819819
    820820!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     
    877877      d_th_z(:)=0.
    878878      d_q_z(:)=0.
     879      d_u_z(:)=0.
     880      d_v_z(:)=0.
    879881      d_t_dyn_z(:)=0.
    880882      d_th_dyn_z(:)=0.
    881883      d_q_dyn_z(:)=0.
     884      d_u_dyn_z(:)=0.
     885      d_v_dyn_z(:)=0.
    882886      DO l=2,llm-1
    883887       d_t_z(l)=(temp(l+1)-temp(l-1))/(play(l+1)-play(l-1))
    884888       d_th_z(l)=(teta(l+1)-teta(l-1))/(play(l+1)-play(l-1))
    885889       d_q_z(l)=(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1))
     890       d_u_z(l)=(u(l+1)-u(l-1))/(play(l+1)-play(l-1))
     891       d_v_z(l)=(v(l+1)-v(l-1))/(play(l+1)-play(l-1))
    886892      ENDDO
    887893      d_t_z(1)=d_t_z(2)
    888894      d_th_z(1)=d_th_z(2)
    889895      d_q_z(1)=d_q_z(2)
     896      d_u_z(1)=d_u_z(2)
     897      d_v_z(1)=d_v_z(2)
    890898      d_t_z(llm)=d_t_z(llm-1)
    891899      d_th_z(llm)=d_th_z(llm-1)
    892900      d_q_z(llm)=d_q_z(llm-1)
     901      d_u_z(llm)=d_u_z(llm-1)
     902      d_v_z(llm)=d_v_z(llm-1)
    893903
    894904!Calcul de l advection verticale
    895       d_t_dyn_z(:)=w_mod_cas(:)*d_t_z(:)
    896       d_th_dyn_z(:)=w_mod_cas(:)*d_th_z(:)
    897       d_q_dyn_z(:)=w_mod_cas(:)*d_q_z(:)
    898 
     905! Modif w_mod_cas -> omega_mod_cas (MM+MPL 20170310)
     906      d_t_dyn_z(:)=omega_mod_cas(:)*d_t_z(:)
     907      d_th_dyn_z(:)=omega_mod_cas(:)*d_th_z(:)
     908      d_q_dyn_z(:)=omega_mod_cas(:)*d_q_z(:)
     909      d_u_dyn_z(:)=omega_mod_cas(:)*d_u_z(:)
     910      d_v_dyn_z(:)=omega_mod_cas(:)*d_v_z(:)
     911
     912!geostrophic wind
     913      if (forc_geo.eq.1) then
     914        do l=1,llm
     915        ug(l) = ug_mod_cas(l)
     916        vg(l) = vg_mod_cas(l)
     917        enddo
     918      endif
    899919!wind nudging
    900920      if (nudging_u.gt.0.) then
     
    902922           u(l)=u(l)+timestep*(u_mod_cas(l)-u(l))/(nudge_u)
    903923        enddo
    904       else
    905         do l=1,llm
    906         ug(l) = u_mod_cas(l)
    907         enddo
     924!     else
     925!       do l=1,llm
     926!          u(l) = u_mod_cas(l)
     927!       enddo
    908928      endif
    909929
     
    912932           v(l)=v(l)+timestep*(v_mod_cas(l)-v(l))/(nudge_v)
    913933        enddo
    914       else
    915         do l=1,llm
    916         vg(l) = v_mod_cas(l)
    917         enddo
     934!     else
     935!       do l=1,llm
     936!          v(l) = v_mod_cas(l)
     937!       enddo
    918938      endif
    919939
     
    922942           w(l)=w(l)+timestep*(w_mod_cas(l)-w(l))/(nudge_w)
    923943        enddo
    924       else
    925         do l=1,llm
    926         w(l) = w_mod_cas(l)
    927         enddo
     944 !    else
     945 !      do l=1,llm
     946        w(l) = w_mod_cas(l)
     947 !      enddo
    928948      endif
    929949
     
    941961
    942962      do l = 1, llm
    943        omega(l) = w_mod_cas(l)
     963! Modif w_mod_cas -> omega_mod_cas (MM+MPL 20170309)
     964       omega(l) = omega_mod_cas(l)
    944965       omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
    945966       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    946967
    947 !calcul advection
    948 !       if ((tend_u.eq.1).and.(tend_w.eq.0)) then
    949 !          d_u_adv(l)=du_mod_cas(l)
    950 !       else if ((tend_u.eq.1).and.(tend_w.eq.1)) then
    951 !          d_u_adv(l)=hu_mod_cas(l)-d_u_dyn_z(l)
     968!calcul advections
     969        if ((forc_u.eq.1).and.(forc_w.eq.0)) then
     970           d_u_adv(l)=du_mod_cas(l)
     971        else if ((forc_u.eq.1).and.(forc_w.eq.1)) then
     972           d_u_adv(l)=hu_mod_cas(l)-d_u_dyn_z(l)
     973        endif
     974
     975        if ((forc_v.eq.1).and.(forc_w.eq.0)) then
     976           d_v_adv(l)=dv_mod_cas(l)
     977        else if ((forc_v.eq.1).and.(forc_w.eq.1)) then
     978           d_v_adv(l)=hv_mod_cas(l)-d_v_dyn_z(l)
     979        endif
     980
     981! Puisque dth a ete converti en dt, on traite de la meme facon
     982! les flags tadv et thadv
     983        if ((tadv.eq.1.or.thadv.eq.1) .and. (forc_w.eq.0)) then
     984!          d_t_adv(l)=alpha*omega(l)/rcpd-dt_mod_cas(l)
     985           d_t_adv(l)=alpha*omega(l)/rcpd+dt_mod_cas(l)
     986        else if ((tadv.eq.1.or.thadv.eq.1) .and. (forc_w.eq.1)) then
     987!          d_t_adv(l)=alpha*omega(l)/rcpd-ht_mod_cas(l)-d_t_dyn_z(l)
     988           d_t_adv(l)=alpha*omega(l)/rcpd+ht_mod_cas(l)-d_t_dyn_z(l)
     989        endif
     990
     991!       if ((thadv.eq.1) .and. (forc_w.eq.0)) then
     992!          d_t_adv(l)=alpha*omega(l)/rcpd-dth_mod_cas(l)
     993!          d_t_adv(l)=alpha*omega(l)/rcpd+dth_mod_cas(l)
     994!       else if ((thadv.eq.1) .and. (forc_w.eq.1)) then
     995!          d_t_adv(l)=alpha*omega(l)/rcpd-hth_mod_cas(l)-d_t_dyn_z(l)
     996!          d_t_adv(l)=alpha*omega(l)/rcpd+hth_mod_cas(l)-d_t_dyn_z(l)
    952997!       endif
    953 !
    954 !       if ((tend_v.eq.1).and.(tend_w.eq.0)) then
    955 !          d_v_adv(l)=dv_mod_cas(l)
    956 !       else if ((tend_v.eq.1).and.(tend_w.eq.1)) then
    957 !          d_v_adv(l)=hv_mod_cas(l)-d_v_dyn_z(l)
    958 !       endif
    959 !
    960 !-----------------------------------------------------
    961         if (tadv.eq.1 .or. tadvh.eq.1) then
    962            d_t_adv(l)=alpha*omega(l)/rcpd-dt_mod_cas(l)
    963         else if (tadvv.eq.1) then
    964 ! ATTENTION d_t_dyn_z pas calcule (voir twpice)
    965            d_t_adv(l)=alpha*omega(l)/rcpd-ht_mod_cas(l)-d_t_dyn_z(l)
    966         endif
    967         print *,'interp_case d_t_dyn_z=',d_t_dyn_z(l),d_q_dyn_z(l)
    968 
    969 ! Verifier le signe !!
    970         if (thadv.eq.1 .or. thadvh.eq.1) then
    971            d_th_adv(l)=dth_mod_cas(l)
    972            print *,'dthadv=',d_th_adv(l)*86400.
    973         else if (thadvv.eq.1) then
    974            d_th_adv(l)=hth_mod_cas(l)-d_th_dyn_z(l)
    975         endif
    976  
    977 ! Verifier le signe !!
    978         if ((qadv.eq.1).and.(forc_w.eq.0)) then
     998
     999        if ((qadv.eq.1) .and. (forc_w.eq.0)) then
    9791000           d_q_adv(l,1)=dq_mod_cas(l)
    980         else if ((qadvh.eq.1).and.(forc_w.eq.1)) then
     1001!          d_q_adv(l,1)=-1*dq_mod_cas(l)
     1002        else if ((qadv.eq.1) .and. (forc_w.eq.1)) then
    9811003           d_q_adv(l,1)=hq_mod_cas(l)-d_q_dyn_z(l)
     1004!          d_q_adv(l,1)=-1*hq_mod_cas(l)-d_q_dyn_z(l)
    9821005        endif
    9831006         
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/dyn1d/1D_nudge_sandu_astex.h

    r2239 r2925  
    1818!     print *,'l dq1 relax dqadv',l,dq(l,1),relax_q(l,1),d_q_adv(l,1)
    1919!     print *,'l dq2 relax dqadv',l,dq(l,2),relax_q(l,2),d_q_adv(l,2)
    20 !     print *,'l dt  relax dtadv',l,dt_phys(l),relax_thl(l),d_th_adv(l)
     20!     print *,'l dt  relax dtadv',l,dt_phys(l),relax_thl(l),d_t_adv(l)
    2121      enddo
    2222        u(1:mxcalc)=u(1:mxcalc) + timestep*(                                &         
     
    3232     &      dq(1:mxcalc,2) - relax_q(1:mxcalc,2)+d_q_adv(1:mxcalc,2))
    3333        temp(1:mxcalc)=temp(1:mxcalc)+timestep*(                            &
    34      &      dt_phys(1:mxcalc)-relax_thl(1:mxcalc)+d_th_adv(1:mxcalc))
     34     &      dt_phys(1:mxcalc)-relax_thl(1:mxcalc)+d_t_adv(1:mxcalc))
     35
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/dyn1d/1D_read_forc_cases.h

    r2716 r2925  
    7676        dt_cooling(l)=thlpcar(kmax)-frac*(thlpcar(kmax)-thlpcar(kmax-1))
    7777        do k=2,kmax
     78          print *,'k l height(k) height(k-1) zlay(l) frac=',k,l,height(k),height(k-1),zlay(l),frac
    7879          frac = (height(k)-zlay(l))/(height(k)-height(k-1))
    7980          if(l==1) print*,'k, height, tttprof',k,height(k),tttprof(k)
     
    227228!?       omega2(l)=-rho(l)*omega(l)
    228229       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    229        d_th_adv(l) = alpha*omega(l)/rcpd-(ht_mod(l)+vt_mod(l))
     230       d_t_adv(l) = alpha*omega(l)/rcpd-(ht_mod(l)+vt_mod(l))
    230231       d_q_adv(l,1) = -(hq_mod(l)+vq_mod(l))
    231232       d_q_adv(l,2) = 0.0
     
    280281!on applique le forcage total au premier pas de temps
    281282!attention: signe different de toga
    282        d_th_adv(l) = alpha*omega(l)/rcpd+(ht_mod(l)+vt_mod(l))
     283       d_t_adv(l) = alpha*omega(l)/rcpd+(ht_mod(l)+vt_mod(l))
    283284       d_q_adv(l,1) = (hq_mod(l)+vq_mod(l))
    284285       d_q_adv(l,2) = 0.0
     
    341342!on applique le forcage total au premier pas de temps
    342343!attention: signe different de toga
    343        d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)
     344       d_t_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)
    344345!forcage en th
    345 !       d_th_adv(l) = ht_mod(l)
     346!       d_t_adv(l) = ht_mod(l)
    346347       d_q_adv(l,1) = hq_mod(l)
    347348       d_q_adv(l,2) = 0.0
     
    442443!on applique le forcage total au premier pas de temps
    443444!attention: signe different de toga
    444        d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)
     445       d_t_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)
    445446!forcage en th
    446 !       d_th_adv(l) = ht_mod(l)
     447!       d_t_adv(l) = ht_mod(l)
    447448       d_q_adv(l,1) = hq_mod(l)
    448449       d_q_adv(l,2) = 0.0
     
    557558!on applique le forcage total au premier pas de temps
    558559!attention: signe different de toga
    559 !      d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)
     560!      d_t_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)
    560561!forcage en th
    561        d_th_adv(l) = ht_mod(l)
     562       d_t_adv(l) = ht_mod(l)
    562563       d_q_adv(l,1) = hq_mod(l)
    563564       d_q_adv(l,2) = 0.0
     
    657658! Advective forcings are given in K or g/kg ... per HOUR
    658659!      IF(height(l).LT.1000) THEN
    659 !        d_th_adv(l) = (adv_theta_prof + rad_theta_prof)/3600.
     660!        d_t_adv(l) = (adv_theta_prof + rad_theta_prof)/3600.
    660661!        d_q_adv(l,1) = adv_qt_prof/1000./3600.
    661662!        d_q_adv(l,2) = 0.0
    662663!      ELSEIF (height(l).GE.1000.AND.height(l).LT.3000) THEN
    663 !        d_th_adv(l) = (adv_theta_prof + rad_theta_prof)*
     664!        d_t_adv(l) = (adv_theta_prof + rad_theta_prof)*
    664665!    :               (1-(height(l)-1000.)/2000.)
    665 !        d_th_adv(l) = d_th_adv(l)/3600.
     666!        d_t_adv(l) = d_t_adv(l)/3600.
    666667!        d_q_adv(l,1) = adv_qt_prof*(1-(height(l)-1000.)/2000.)
    667668!        d_q_adv(l,1) = d_q_adv(l,1)/1000./3600.
    668669!        d_q_adv(l,2) = 0.0
    669670!      ELSE
    670 !        d_th_adv(l) = 0.0
     671!        d_t_adv(l) = 0.0
    671672!        d_q_adv(l,1) = 0.0
    672673!        d_q_adv(l,2) = 0.0
     
    751752!?       omega2(l)=-rho(l)*omega(l)
    752753       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    753 !      d_th_adv(l) = alpha*omega(l)/rcpd+vt_mod(l)
     754!      d_t_adv(l) = alpha*omega(l)/rcpd+vt_mod(l)
    754755!      d_q_adv(l,1) = vq_mod(l)
    755        d_th_adv(l) = alpha*omega(l)/rcpd
     756       d_t_adv(l) = alpha*omega(l)/rcpd
    756757       d_q_adv(l,1) = 0.0
    757758       d_q_adv(l,2) = 0.0
     
    827828!      omega2(l)=-rho(l)*omega(l)
    828829       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    829 !      d_th_adv(l) = alpha*omega(l)/rcpd+vt_mod(l)
     830!      d_t_adv(l) = alpha*omega(l)/rcpd+vt_mod(l)
    830831!      d_q_adv(l,1) = vq_mod(l)
    831        d_th_adv(l) = alpha*omega(l)/rcpd
     832       d_t_adv(l) = alpha*omega(l)/rcpd
    832833       d_q_adv(l,1) = 0.0
    833834       d_q_adv(l,2) = 0.0
     
    890891!on applique le forcage total au premier pas de temps
    891892!attention: signe different de toga
    892        d_th_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l))
     893       d_t_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l))
    893894       d_q_adv(l,1) = (hq_mod_cas(l)+vq_mod_cas(l))
    894895       d_q_adv(l,2) = 0.0
    895896       d_u_adv(l) = (hu_mod_cas(l)+vu_mod_cas(l))
    896        d_u_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l))
     897! correction bug d_u -> d_v (MM+MPL 20170310)
     898       d_v_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l))
    897899      enddo     
    898900
     
    971973       q(l,2) = ql_mod_cas(l)
    972974       u(l) = u_mod_cas(l)
    973        ug(l)= u_mod_cas(l)
     975       ug(l)= ug_mod_cas(l)
    974976       v(l) = v_mod_cas(l)
    975        vg(l)= v_mod_cas(l)
    976        omega(l) = w_mod_cas(l)
     977       vg(l)= vg_mod_cas(l)
     978! Modif w_mod_cas -> omega_mod_cas (MM+MPL 20170309)
     979       omega(l) = omega_mod_cas(l)
    977980       omega2(l)=omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
    978981
     
    980983!on applique le forcage total au premier pas de temps
    981984!attention: signe different de toga
    982        d_th_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l))
     985       d_t_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l))
    983986       d_t_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l))
    984987!      d_q_adv(l,1) = (hq_mod_cas(l)+vq_mod_cas(l))
     
    987990!      d_u_adv(l) = (hu_mod_cas(l)+vu_mod_cas(l))
    988991       d_u_adv(l) = du_mod_cas(l)
    989 !      d_u_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l))
    990        d_u_adv(l) = dv_mod_cas(l)
     992!      d_v_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l))
     993! correction bug d_u -> d_v (MM+MPL 20170310)
     994       d_v_adv(l) = dv_mod_cas(l)
    991995      enddo     
    992996
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/dyn1d/compar1d.h

    r2716 r2925  
    4141
    4242      integer :: tadv, tadvv, tadvh, qadv, qadvv, qadvh, thadv, thadvv, thadvh, trad
    43       integer :: forc_omega, forc_w, forc_geo, forc_ustar
     43      integer :: forc_omega, forc_u, forc_v, forc_w, forc_geo, forc_ustar
    4444      real    :: nudging_u, nudging_v, nudging_w, nudging_t, nudging_q
    4545      common/com_par1d/                                                 &
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/dyn1d/lmdz1d.F90

    r2789 r2925  
    194194      real :: du_phys(llm),dv_phys(llm),dt_phys(llm)
    195195      real :: dt_dyn(llm)
    196       real :: dt_cooling(llm),d_t_adv(llm),d_th_adv(llm),d_t_nudge(llm)
     196      real :: dt_cooling(llm),d_t_adv(llm),d_t_nudge(llm)
    197197      real :: d_u_nudge(llm),d_v_nudge(llm)
    198198      real :: du_adv(llm),dv_adv(llm)
     
    206206      REAL, ALLOCATABLE, DIMENSION(:,:):: d_q_adv
    207207      REAL, ALLOCATABLE, DIMENSION(:,:):: d_q_nudge
    208 !      REAL, ALLOCATABLE, DIMENSION(:):: d_th_adv
     208!      REAL, ALLOCATABLE, DIMENSION(:):: d_t_adv
    209209
    210210!---------------------------------------------------------------------
     
    271271      dt_dyn(:)=0.
    272272      dt_cooling(:)=0.
    273       d_th_adv(:)=0.
     273      d_t_adv(:)=0.
    274274      d_t_nudge(:)=0.
    275275      d_u_nudge(:)=0.
     
    404404       heure_ini_cas=0.
    405405       pdt_cas=1800.         ! forcing frequency
     406      elseif (forcing_type .eq.105) THEN ! bomex starts 16-12-2004 0h
     407       forcing_case2 = .true.
     408       year_ini_cas=1969
     409       mth_ini_cas=6
     410       day_deb=24
     411       heure_ini_cas=0.
     412       pdt_cas=1800.         ! forcing frequency
    406413      elseif (forcing_type .eq.40) THEN
    407414       forcing_GCSSold = .true.
     
    573580      allocate(d_q_adv(llm,nqtot))
    574581      allocate(d_q_nudge(llm,nqtot))
    575 !      allocate(d_th_adv(llm))
     582!      allocate(d_t_adv(llm))
    576583
    577584      q(:,:) = 0.
     
    10571064         fcoriolis=0.0
    10581065         dt_cooling=0.0
    1059          d_th_adv=0.0
     1066         d_t_adv=0.0
    10601067         d_q_adv=0.0
    10611068       endif
     
    11431150        if (prt_level.ge.3) then
    11441151          print *,                                                          &
    1145      &    'physiq-> temp(1),dt_phys(1),d_th_adv(1),dt_cooling(1) ',         &
    1146      &              temp(1),dt_phys(1),d_th_adv(1),dt_cooling(1)
     1152     &    'physiq-> temp(1),dt_phys(1),d_t_adv(1),dt_cooling(1) ',         &
     1153     &              temp(1),dt_phys(1),d_t_adv(1),dt_cooling(1)
    11471154           print* ,'dv_phys=',dv_phys
    11481155           print* ,'dv_age=',dv_age
     
    11551162        temp(1:mxcalc)=temp(1:mxcalc)+timestep*(                            &
    11561163     &              dt_phys(1:mxcalc)                                       &
    1157      &             +d_th_adv(1:mxcalc)                                      &
     1164     &             +d_t_adv(1:mxcalc)                                      &
    11581165     &             +d_t_nudge(1:mxcalc)                                      &
    11591166     &             +dt_cooling(1:mxcalc))  ! Taux de chauffage ou refroid.
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/ener_conserv.F90

    r2895 r2925  
    2727USE phys_state_var_mod, ONLY : du_gwd_front,du_gwd_rando
    2828USE phys_output_var_mod, ONLY : bils_ec,bils_ech,bils_tke,bils_kinetic,bils_enthalp,bils_latent,bils_diss
     29USE add_phys_tend_mod, ONLY : fl_cor_ebil
     30
    2931
    3032IMPLICIT none
     
    6264   DO k = 1, klev
    6365   DO i = 1, klon
    64       ZRCPD = RCPD*(1.0+RVTMP2*(pqn(i,k)+pqln(i,k)+pqsn(i,k)))
    65       d_t_ec(i,k)=0.5/ZRCPD &
    66  &      *(puo(i,k)**2+pvo(i,k)**2-pun(i,k)**2-pvn(i,k)**2)
    67       ENDDO
    68       ENDDO
     66     IF (fl_cor_ebil .GT. 0) then
     67       ZRCPD = RCPD*(1.0+RVTMP2*(pqn(i,k)+pqln(i,k)+pqsn(i,k)))
     68     ELSE
     69       ZRCPD = RCPD*(1.0+RVTMP2*pqn(i,k))
     70     ENDIF
     71     d_t_ec(i,k)=0.5/ZRCPD &
     72 &     *(puo(i,k)**2+pvo(i,k)**2-pun(i,k)**2-pvn(i,k)**2)
     73   ENDDO
     74   ENDDO
    6975!-jld ec_conser
    7076
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/fisrtilp.F90

    r2885 r2925  
    103103  ! Options du programme:
    104104  !
    105   REAL seuil_neb ! un nuage existe vraiment au-dela
    106   PARAMETER (seuil_neb=0.001)
     105  REAL, SAVE :: seuil_neb=0.001 ! un nuage existe vraiment au-dela
     106  !$OMP THREADPRIVATE(seuil_neb)
     107
    107108
    108109  INTEGER ninter ! sous-intervals pour la precipitation
    109110  PARAMETER (ninter=5)
    110   LOGICAL evap_prec ! evaporation de la pluie
    111   PARAMETER (evap_prec=.TRUE.)
     111  INTEGER,SAVE :: iflag_evap_prec=1 ! evaporation de la pluie
     112  !$OMP THREADPRIVATE(iflag_evap_prec)
    112113  !
    113114  LOGICAL cpartiel ! condensation partielle
     
    124125  !
    125126  INTEGER i, k, n, kk
     127  INTEGER,save::itap=0
     128  !$OMP THREADPRIVATE(itap)
     129
    126130  REAL qsl, qsi
    127131  real zct      ,zcl
     
    154158  REAL t_glace_min_old
    155159  REAL zdz(klon),zrho(klon),ztot      , zrhol(klon)
    156   REAL zchau      ,zfroi      ,zfice(klon),zneb(klon)
     160  REAL zchau      ,zfroi      ,zfice(klon),zneb(klon),znebprecip(klon)
    157161  REAL zmelt, zpluie, zice
    158162  REAL dzfice(klon)
     
    212216!CR: pour iflag_ice_thermo=2, on active que la convection
    213217!  ice_thermo = iflag_ice_thermo .GE. 1
     218
     219itap=itap+1
     220znebprecip(:)=0.
     221
    214222   ice_thermo = (iflag_ice_thermo .EQ. 1).OR.(iflag_ice_thermo .GE. 3)
    215223  zdelq=0.0
     
    218226  IF (appel1er) THEN
    219227     CALL getin_p('iflag_oldbug_fisrtilp',iflag_oldbug_fisrtilp)
     228     CALL getin_p('iflag_evap_prec',iflag_evap_prec)
     229     CALL getin_p('seuil_neb',seuil_neb)
    220230     write(lunout,*)' iflag_oldbug_fisrtilp =',iflag_oldbug_fisrtilp
    221231     !
    222232     WRITE(lunout,*) 'fisrtilp, ninter:', ninter
    223      WRITE(lunout,*) 'fisrtilp, evap_prec:', evap_prec
     233     WRITE(lunout,*) 'fisrtilp, iflag_evap_prec:', iflag_evap_prec
    224234     WRITE(lunout,*) 'fisrtilp, cpartiel:', cpartiel
     235     
    225236     IF (ABS(dtime/REAL(ninter)-360.0).GT.0.001) THEN
    226237        WRITE(lunout,*) 'fisrtilp: Ce n est pas prevu, voir Z.X.Li', dtime
     
    387398     !   - zmqc: masse de precip qui doit etre thermalisee
    388399     !
    389      IF (evap_prec) THEN
     400     IF (iflag_evap_prec>=1) THEN
    390401        DO i = 1, klon
    391402!          S'il y a des precipitations
     
    480491!        S'il y a des precipitations
    481492         IF (zrfl(i)+zifl(i).GT.0.) THEN
     493
     494         IF (iflag_evap_prec==1) THEN
     495            znebprecip(i)=zneb(i)
     496         ELSE
     497            znebprecip(i)=MAX(zneb(i),znebprecip(i))
     498         ENDIF
    482499     
    483500        ! Evap max pour ne pas saturer la fraction sous le nuage
    484          zqev0 = MAX (0.0, (zqs(i)-zq(i))*zneb(i) )
     501         zqev0 = MAX (0.0, (zqs(i)-zq(i))*znebprecip(i) )
    485502
    486503         !JAM
     
    567584         zrfl(i) = zrfln(i)
    568585         zifl(i) = zifln(i)
     586!        print*,'REEVAP ',itap,k,znebprecip(1),zqev0,zqev,zqevi,zrfl(1)
    569587
    570588!CR ATTENTION: deplacement de la fonte de la glace
     
    595613           end if
    596614
     615           ELSE
     616              ! Si on n'a plus de pluies, on reinitialise a 0 la farcion
     617              ! sous nuageuse utilisee pour la pluie.
     618              znebprecip(i)=0.
    597619           ENDIF ! (zrfl(i)+zifl(i).GT.0.)
    598620        ENDDO
     
    603625     ! Fin evaporation de la precipitation
    604626     ! ----------------------------------------------------------------
    605      ENDIF ! (evap_prec)
     627     ENDIF ! (iflag_evap_prec>=1)
    606628     !
    607629     ! Calculer Qs et L/Cp*dQs/dT:
     
    702724                   qcloud,ctot,zpspsk,paprs,ztla,zthl, &
    703725                   ratqs,zqs,t)
    704               elseif (iflag_cloudth_vert==3) then
     726              elseif (iflag_cloudth_vert>=3) then
    705727               call cloudth_v3(klon,klev,k,ztv, &
    706728                   zq,zqta,fraca, &
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/oasis.F90

    r2872 r2925  
    8080
    8181  TYPE(FLD_CPL), DIMENSION(maxsend), SAVE, PUBLIC :: infosend   ! Information for sending coupling fields
     82!$OMP THREADPRIVATE(infosend)
    8283  TYPE(FLD_CPL), DIMENSION(maxrecv), SAVE, PUBLIC :: inforecv   ! Information for receiving coupling fields
     84!$OMP THREADPRIVATE(inforecv)
    8385 
    8486  LOGICAL,SAVE :: cpl_current
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/physiq_mod.F90

    r2924 r2925  
    245245
    246246    USE cmp_seri_mod
    247     USE add_phys_tend_mod, only : add_pbl_tend, add_phys_tend, prt_enerbil, &
     247    USE add_phys_tend_mod, only : add_pbl_tend, add_phys_tend, diag_phys_tend, prt_enerbil, &
    248248  &      fl_ebil, fl_cor_ebil
    249249
     
    16821682               pdtphys, &
    16831683               annee_ref, &
     1684               year_cur, &
    16841685               day_ref,  &
    16851686               day_ini, &
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/sisvat/sisvat.F

    r1990 r2925  
    70707070
    70717071      DO ig = 1,knonv
    7072         ztherm_i(ig)   = inertie_ice
     7072        ztherm_i(ig)   = inertie_lic
    70737073        IF (isnoSV(ig) > 0) ztherm_i(ig)   = inertie_sno
    70747074      ENDDO
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/soil.F90

    r2311 r2925  
    152152!-----------------------------------------------------------------------
    153153!   Calcul de l'inertie thermique a partir de la variable rnat.
    154 !   on initialise a inertie_ice meme au-dessus d'un point de mer au cas
     154!   on initialise a inertie_sic meme au-dessus d'un point de mer au cas
    155155!   ou le point de mer devienne point de glace au pas suivant
    156156!   on corrige si on a un point de terre avec ou sans glace
    157157!
    158 !-----------------------------------------------------------------------
     158!   iophys can be used to write the ztherm_i variable in a phys.nc file
     159!   and check the results; to do so, add "CALL iophys_ini" in physiq_mod
     160!   and add knindex to the list of inputs in all the calls to soil.F90
     161!   (and to soil.F90 itself !)
     162!-----------------------------------------------------------------------
     163
    159164  IF (indice == is_sic) THEN
    160165     DO ig = 1, knon
    161         ztherm_i(ig)   = inertie_ice
     166        ztherm_i(ig)   = inertie_sic
     167     ENDDO
     168     IF (iflag_sic == 0) THEN
     169       DO ig = 1, knon
     170         IF (snow(ig) > 0.0) ztherm_i(ig)   = inertie_sno
     171       ENDDO
     172!      Otherwise sea-ice keeps the same inertia, even when covered by snow
     173     ENDIF
     174!    CALL iophys_ecrit_index('ztherm_sic', 1, 'ztherm_sic', 'USI', &
     175!      knon, knindex, ztherm_i)
     176  ELSE IF (indice == is_lic) THEN
     177     DO ig = 1, knon
     178        ztherm_i(ig)   = inertie_lic
    162179        IF (snow(ig) > 0.0) ztherm_i(ig)   = inertie_sno
    163180     ENDDO
    164   ELSE IF (indice == is_lic) THEN
    165      DO ig = 1, knon
    166         ztherm_i(ig)   = inertie_ice
    167         IF (snow(ig) > 0.0) ztherm_i(ig)   = inertie_sno
    168      ENDDO
     181!    CALL iophys_ecrit_index('ztherm_lic', 1, 'ztherm_lic', 'USI', &
     182!      knon, knindex, ztherm_i)
    169183  ELSE IF (indice == is_ter) THEN
    170184     DO ig = 1, knon
     
    172186        IF (snow(ig) > 0.0) ztherm_i(ig)   = inertie_sno
    173187     ENDDO
     188!    CALL iophys_ecrit_index('ztherm_ter', 1, 'ztherm_ter', 'USI', &
     189!      knon, knindex, ztherm_i)
    174190  ELSE IF (indice == is_oce) THEN
    175191     DO ig = 1, knon
    176         ztherm_i(ig)   = inertie_ice
    177      ENDDO
     192!       This is just in case, but SST should be used by the model anyway
     193        ztherm_i(ig)   = inertie_sic
     194     ENDDO
     195!    CALL iophys_ecrit_index('ztherm_oce', 1, 'ztherm_oce', 'USI', &
     196!      knon, knindex, ztherm_i)
    178197  ELSE
    179198     WRITE(lunout,*) "valeur d indice non prevue", indice
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/surf_land_orchidee_mod.F90

    r2924 r2925  
    3232  PUBLIC  :: surf_land_orchidee
    3333
    34   LOGICAL, ALLOCATABLE, SAVE :: flag_omp(:)
    3534CONTAINS
    3635!
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/surf_land_orchidee_noopenmp_mod.F90

    r2572 r2925  
    144144    INTEGER                                   :: error
    145145    INTEGER, SAVE                             :: nb_fields_cpl ! number of fields for the climate-carbon coupling (between ATM and ORCHIDEE).
     146    !$OMP THREADPRIVATE(nb_fields_cpl)
    146147    REAL, SAVE, ALLOCATABLE, DIMENSION(:,:)   :: fields_cpl    ! Fluxes for the climate-carbon coupling
     148    !$OMP THREADPRIVATE(fields_cpl)
    147149    REAL, DIMENSION(klon)                     :: swdown_vrai
    148150    CHARACTER (len = 20)                      :: modname = 'surf_land_orchidee'
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/surf_land_orchidee_noz0h_mod.F90

    r2571 r2925  
    3333  PUBLIC  :: surf_land_orchidee
    3434
    35   LOGICAL, ALLOCATABLE, SAVE :: flag_omp(:)
    3635CONTAINS
    3736!
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/surface_data.F90

    r2538 r2925  
    99 
    1010  LOGICAL, SAVE          :: ok_veget      ! true for use of vegetation model ORCHIDEE
    11   CHARACTER(len=10), SAVE :: type_veget    ! orchidee/y/bucket/n/betaclim
    1211  !$OMP THREADPRIVATE(ok_veget)
    13   ! Martin
     12
     13  CHARACTER(len=10), SAVE :: type_veget   ! orchidee/y/bucket/n/betaclim
     14  !$OMP THREADPRIVATE(type_veget)
     15
    1416  LOGICAL, SAVE          :: ok_snow       ! true for coupling to snow model SISVAT
    1517  !$OMP THREADPRIVATE(ok_snow)
    16   ! Martin
    1718
    1819  CHARACTER(len=6), SAVE :: type_ocean    ! force/slab/couple
  • LMDZ5/branches/LMDZ_tree_FC/libf/phylmd/wake.F90

    r2761 r2925  
    178178  LOGICAL, SAVE                                         :: flag_wk_check_trgl
    179179  !$OMP THREADPRIVATE(flag_wk_check_trgl)
     180  INTEGER, SAVE                                         :: iflag_wk_check_trgl
     181  !$OMP THREADPRIVATE(iflag_wk_check_trgl)
    180182
    181183  REAL                                                  :: delta_t_min
     
    347349  CALL getin_p('flag_wk_check_trgl ', flag_wk_check_trgl)
    348350  WRITE(*,*) 'flag_wk_check_trgl=', flag_wk_check_trgl
     351  WRITE(*,*) 'flag_wk_check_trgl OBSOLETE. Utilisr iflag_wk_check_trgl plutot'
     352  iflag_wk_check_trgl=0 ; IF (flag_wk_check_trgl) iflag_wk_check_trgl=1
     353  CALL getin_p('iflag_wk_check_trgl ', iflag_wk_check_trgl)
     354  WRITE(*,*) 'iflag_wk_check_trgl=', iflag_wk_check_trgl
    349355
    350356  first=.false.
     
    17931799  ! Filter out bad wakes
    17941800
    1795   IF (flag_wk_check_trgl) THEN
     1801  IF (iflag_wk_check_trgl>=1) THEN
    17961802    ! Check triangular shape of dth profile
    17971803    DO i = 1, klon
     
    18051811          wape2(i) = -1.
    18061812          !! print *,'wake, rej 1'
    1807         ELSE IF (abs(2.*sum_dth(i)/(hw0(i)*dthmin(i)) - 1.) > 0.5) THEN
     1813        ELSE IF (iflag_wk_check_trgl==1.AND.abs(2.*sum_dth(i)/(hw0(i)*dthmin(i)) - 1.) > 0.5) THEN
    18081814          wape2(i) = -1.
    18091815          !! print *,'wake, rej 2'
Note: See TracChangeset for help on using the changeset viewer.