Ignore:
Timestamp:
Apr 3, 2016, 12:09:34 AM (8 years ago)
Author:
Laurent Fairhead
Message:

Merged trunk changes r2457:2487 into testing branch

Location:
LMDZ5/branches/testing
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/branches/testing

  • LMDZ5/branches/testing/libf/phylmd/cv3_routines.F90

    r2435 r2488  
    77SUBROUTINE cv3_param(nd, k_upper, delt)
    88
     9  USE ioipsl_getin_p_mod, ONLY : getin_p
    910  use mod_phys_lmdz_para
    1011  IMPLICIT NONE
    … …  
    3940  INTEGER, INTENT(IN)              :: k_upper
    4041  REAL, INTENT(IN)                 :: delt ! timestep (seconds)
    41 
    4242
    4343! Local variables
    … …  
    6565
    6666  IF (first) THEN
    67 
    6867! -- "microphysical" parameters:
    69     sigdz = 0.01
    70     spfac = 0.15
    71     pbcrit = 150.0
    72     ptcrit = 500.0
    7368! IM beg: ajout fis. reglage ep
    74     flag_epkeorig = 1
    75     elcrit = 0.0003
    76     tlcrit = -55.0
     69! CR+JYG: shedding coefficient (used when iflag_mix_adiab=1)
    7770! IM lu dans physiq.def via conf_phys.F90     epmax  = 0.993
    7871
    7972    omtrain = 45.0 ! used also for snow (no disctinction rain/snow)
    80 
    8173! -- misc:
    82 
    8374    dtovsh = -0.2 ! dT for overshoot
    84     dpbase = -40. ! definition cloud base (400m above LCL)
    8575! cc      dttrig = 5.   ! (loose) condition for triggering
    8676    dttrig = 10. ! (loose) condition for triggering
    87     flag_wb = 1
    88     wbmax = 6. ! (m/s) adiab updraught speed at LFC (used in cv3p1_closure)
    89 
    90 ! -- rate of approach to quasi-equilibrium:
    91 
    9277    dtcrit = -2.0
    93     tau = 8000.
    94 
    9578! -- end of convection
    96 
    97     tau_stop = 15000.
    98     ok_convstop = .False.
    99 
    100     ok_intermittent = .False.
    101 
    10279! -- interface cloud parameterization:
    103 
    10480    delta = 0.01 ! cld
    105 
    10681! -- interface with boundary-layer (gust factor): (sb)
    107 
    10882    betad = 10.0 ! original value (from convect 4.3)
    10983
    110    !$OMP MASTER
    111     OPEN (99, FILE='conv_param.data', STATUS='old', FORM='formatted', ERR=9999)
    112     READ (99, *, END=9998) dpbase
    113     READ (99, *, END=9998) pbcrit
    114     READ (99, *, END=9998) ptcrit
    115     READ (99, *, END=9998) sigdz
    116     READ (99, *, END=9998) spfac
    117     READ (99, *, END=9998) tau
    118     READ (99, *, END=9998) flag_wb
    119     READ (99, *, END=9998) wbmax
    120     READ (99, *, END=9998) ok_convstop
    121     READ (99, *, END=9998) tau_stop
    122     READ (99, *, END=9998) ok_intermittent
    123 9998 CONTINUE
    124     CLOSE (99)
    125 9999 CONTINUE
     84! Var interm pour le getin
     85     dpbase=-40.
     86     CALL getin_p('dpbase',dpbase)
     87     pbcrit=150.0
     88     CALL getin_p('pbcrit',pbcrit)
     89     ptcrit=500.0
     90     CALL getin_p('ptcrit',ptcrit)
     91     sigdz=0.01
     92     CALL getin_p('sigdz',sigdz)
     93     spfac=0.15
     94     CALL getin_p('spfac',spfac)
     95     tau=8000.
     96     CALL getin_p('tau',tau)
     97     flag_wb=1
     98     CALL getin_p('flag_wb',flag_wb)
     99     wbmax=6.
     100     CALL getin_p('wbmax',wbmax)
     101     ok_convstop=.False.
     102     CALL getin_p('ok_convstop',ok_convstop)
     103     tau_stop=15000.
     104     CALL getin_p('tau_stop',tau_stop)
     105     ok_intermittent=.False.
     106     CALL getin_p('ok_intermittent',ok_intermittent)
     107     coef_peel=0.25
     108     CALL getin_p('coef_peel',coef_peel)
     109
     110     flag_epKEorig=1
     111     CALL getin_p('flag_epKEorig',flag_epKEorig)
     112     elcrit=0.0003
     113     CALL getin_p('elcrit',elcrit)
     114     tlcrit=-55.0
     115     CALL getin_p('tlcrit',tlcrit)
     116
    126117    WRITE (*, *) 'dpbase=', dpbase
    127118    WRITE (*, *) 'pbcrit=', pbcrit
    … …  
    130121    WRITE (*, *) 'spfac=', spfac
    131122    WRITE (*, *) 'tau=', tau
    132     WRITE (*, *) 'flag_wb =', flag_wb
    133     WRITE (*, *) 'wbmax =', wbmax
    134     WRITE (*, *) 'ok_convstop =', ok_convstop
    135     WRITE (*, *) 'tau_stop =', tau_stop
    136     WRITE (*, *) 'ok_intermittent =', ok_intermittent
    137 
    138 ! IM Lecture du fichier ep_param.data
    139     OPEN (79, FILE='ep_param.data', STATUS='old', FORM='formatted', ERR=7999)
    140     READ (79, *, END=7998) flag_epkeorig
    141     READ (79, *, END=7998) elcrit
    142     READ (79, *, END=7998) tlcrit
    143 7998 CONTINUE
    144     CLOSE (79)
    145 7999 CONTINUE
    146     WRITE (*, *) 'flag_epKEorig', flag_epkeorig
     123    WRITE (*, *) 'flag_wb=', flag_wb
     124    WRITE (*, *) 'wbmax=', wbmax
     125    WRITE (*, *) 'ok_convstop=', ok_convstop
     126    WRITE (*, *) 'tau_stop=', tau_stop
     127    WRITE (*, *) 'ok_intermittent=', ok_intermittent
     128    WRITE (*, *) 'coef_peel=', coef_peel
     129
     130    WRITE (*, *) 'flag_epKEorig=', flag_epKEorig
    147131    WRITE (*, *) 'elcrit=', elcrit
    148132    WRITE (*, *) 'tlcrit=', tlcrit
    149 ! IM end: ajout fis. reglage ep
    150   !$OMP END MASTER
    151 
    152    CALL bcast(dpbase)
    153    CALL bcast(pbcrit)
    154    CALL bcast(ptcrit)
    155    CALL bcast(sigdz)
    156    CALL bcast(spfac)
    157    CALL bcast(tau)
    158    CALL bcast(flag_wb)
    159    CALL bcast(wbmax)
    160    CALL bcast(ok_convstop)
    161    CALL bcast(tau_stop)
    162    CALL bcast(ok_intermittent)
    163 
    164    CALL bcast(flag_epkeorig)
    165    CALL bcast(elcrit)
    166    CALL bcast(tlcrit)
    167 
    168133    first = .FALSE.
    169 
    170134  END IF ! (first)
    171135
    … …  
    41784142                          ft, fq, fu, fv, ftra, &
    41794143                          Ma, upwd, dnwd, dnwd0, qcondc, wd, cape, &
     4144                          epmax_diag, & ! epmax_cape
    41804145                          iflag1, &
    41814146                          precip1, sig1, w01, &
    41824147                          ft1, fq1, fu1, fv1, ftra1, &
    4183                           Ma1, upwd1, dnwd1, dnwd01, qcondc1, wd1, cape1)
     4148                          Ma1, upwd1, dnwd1, dnwd01, qcondc1, wd1, cape1, &
     4149                          epmax_diag1) ! epmax_cape
    41844150  IMPLICIT NONE
    41854151
    … …  
    41984164  REAL qcondc(nloc, nd)
    41994165  REAL wd(nloc), cape(nloc)
     4166  REAL epmax_diag(nloc)
    42004167
    42014168!outputs:
    … …  
    42094176  REAL qcondc1(nloc, nd)
    42104177  REAL wd1(nloc), cape1(nloc)
     4178  REAL epmax_diag1(len) ! epmax_cape
    42114179
    42124180!local variables:
    … …  
    42184186    wd1(idcum(i)) = wd(i)
    42194187    cape1(idcum(i)) = cape(i)
     4188    epmax_diag1(idcum(i))=epmax_diag(i) ! epmax_cape
    42204189  END DO
    42214190
    … …  
    42524221  RETURN
    42534222END SUBROUTINE cv3_uncompress
     4223
     4224
     4225        subroutine cv3_epmax_fn_cape(nloc,ncum,nd &
     4226                 , ep,hp,icb,inb,clw,nk,t,h,hnk,lv,lf,frac &
     4227                 , pbase, p, ph, tv, buoy, sig, w0,iflag &
     4228                 , epmax_diag)
     4229        implicit none
     4230
     4231        ! On fait varier epmax en fn de la cape
     4232        ! Il faut donc recalculer ep, et hp qui a déjà été calculé et
     4233        ! qui en dépend
     4234        ! Toutes les autres variables fn de ep sont calculées plus bas.
     4235
     4236  include "cvthermo.h"
     4237  include "cv3param.h" 
     4238  include "conema3.h"
     4239  include "cvflag.h"
     4240
     4241! inputs:
     4242      INTEGER, INTENT (IN)                               :: ncum, nd, nloc
     4243      INTEGER, DIMENSION (nloc), INTENT (IN)             :: icb, inb, nk
     4244      REAL, DIMENSION (nloc), INTENT (IN)                :: hnk,pbase
     4245      REAL, DIMENSION (nloc, nd), INTENT (IN)            :: t, lv, lf, tv, h
     4246      REAL, DIMENSION (nloc, nd), INTENT (IN)            :: clw, buoy,frac
     4247      REAL, DIMENSION (nloc, nd), INTENT (IN)            :: sig,w0
     4248      INTEGER, DIMENSION (nloc), INTENT (IN)             :: iflag(nloc)
     4249      REAL, DIMENSION (nloc, nd), INTENT (IN)            :: p
     4250      REAL, DIMENSION (nloc, nd+1), INTENT (IN)          :: ph
     4251! inouts:
     4252      REAL, DIMENSION (nloc, nd), INTENT (INOUT)         :: ep,hp 
     4253! outputs
     4254      REAL, DIMENSION (nloc), INTENT (OUT)           :: epmax_diag
     4255
     4256! local
     4257      integer i,k   
     4258!      real hp_bak(nloc,nd)
     4259!      real ep_bak(nloc,nd)
     4260      real m_loc(nloc,nd)
     4261      real sig_loc(nloc,nd)
     4262      real w0_loc(nloc,nd)
     4263      integer iflag_loc(nloc)
     4264      real cape(nloc)
     4265       
     4266        if (coef_epmax_cape.gt.1e-12) then
     4267
     4268        ! il faut calculer la cape: on fait un calcule simple car tant qu'on ne
     4269        ! connait pas ep, on ne connait pas les mélanges, ddfts etc... qui sont
     4270        ! necessaires au calcul de la cape dans la nouvelle physique
     4271       
     4272!        write(*,*) 'cv3_routines check 4303'
     4273        do i=1,ncum
     4274        do k=1,nd
     4275          sig_loc(i,k)=sig(i,k)
     4276          w0_loc(i,k)=w0(i,k)
     4277          iflag_loc(i)=iflag(i)
     4278!          ep_bak(i,k)=ep(i,k)
     4279        enddo ! do k=1,nd
     4280        enddo !do i=1,ncum
     4281
     4282!        write(*,*) 'cv3_routines check 4311'
     4283!        write(*,*) 'nl=',nl
     4284        CALL cv3_closure(nloc, ncum, nd, icb, inb, & ! na->nd
     4285          pbase, p, ph, tv, buoy, &
     4286          sig_loc, w0_loc, cape, m_loc,iflag_loc)
     4287
     4288!        write(*,*) 'cv3_routines check 4316'
     4289!        write(*,*) 'ep(1,:)=',ep(1,:)
     4290        do i=1,ncum
     4291           epmax_diag(i)=epmax-coef_epmax_cape*sqrt(cape(i))
     4292           epmax_diag(i)=amax1(epmax_diag(i),0.0)
     4293!           write(*,*) 'i,icb,inb,cape,epmax_diag=', &
     4294!                i,icb(i),inb(i),cape(i),epmax_diag(i)
     4295           do k=1,nl
     4296                ep(i,k)=ep(i,k)/epmax*epmax_diag(i)
     4297                ep(i,k)=amax1(ep(i,k),0.0)
     4298                ep(i,k)=amin1(ep(i,k),epmax_diag(i))
     4299           enddo
     4300        enddo
     4301 !       write(*,*) 'ep(1,:)=',ep(1,:)
     4302
     4303      !write(*,*) 'cv3_routines check 4326'
     4304! On recalcule hp:
     4305!      do k=1,nl
     4306!        do i=1,ncum
     4307!         hp_bak(i,k)=hp(i,k)
     4308!       enddo
     4309!      enddo
     4310      do k=1,nl
     4311        do i=1,ncum
     4312          hp(i,k)=h(i,k)
     4313        enddo
     4314      enddo
     4315
     4316  IF (cvflag_ice) THEN
     4317
     4318      do k=minorig+1,nl
     4319       do i=1,ncum
     4320        if((k.ge.icb(i)).and.(k.le.inb(i)))then
     4321          hp(i, k) = hnk(i) + (lv(i,k)+(cpd-cpv)*t(i,k)+frac(i,k)*lf(i,k))* &
     4322                              ep(i, k)*clw(i, k)
     4323        endif
     4324       enddo
     4325      enddo !do k=minorig+1,n
     4326  ELSE !IF (cvflag_ice) THEN
     4327
     4328      DO k = minorig + 1, nl
     4329       DO i = 1, ncum
     4330        IF ((k>=icb(i)) .AND. (k<=inb(i))) THEN
     4331          hp(i,k)=hnk(i)+(lv(i,k)+(cpd-cpv)*t(i,k))*ep(i,k)*clw(i,k)
     4332        endif
     4333       enddo
     4334      enddo !do k=minorig+1,n
     4335
     4336  ENDIF !IF (cvflag_ice) THEN     
     4337      !write(*,*) 'cv3_routines check 4345'
     4338!      do i=1,ncum 
     4339!       do k=1,nl
     4340!        if ((abs(hp_bak(i,k)-hp(i,k))/hp_bak(i,k).gt.1e-1).or. &
     4341!            ((abs(hp_bak(i,k)-hp(i,k))/hp_bak(i,k).gt.1e-4).and. &
     4342!            (ep(i,k)-ep_bak(i,k).lt.1e-4))) then
     4343!           write(*,*) 'i,k=',i,k
     4344!           write(*,*) 'coef_epmax_cape=',coef_epmax_cape
     4345!           write(*,*) 'epmax_diag(i)=',epmax_diag(i)
     4346!           write(*,*) 'ep(i,k)=',ep(i,k)
     4347!           write(*,*) 'ep_bak(i,k)=',ep_bak(i,k)
     4348!           write(*,*) 'hp(i,k)=',hp(i,k)
     4349!           write(*,*) 'hp_bak(i,k)=',hp_bak(i,k)
     4350!           write(*,*) 'h(i,k)=',h(i,k)
     4351!           write(*,*) 'nk(i)=',nk(i)
     4352!           write(*,*) 'h(i,nk(i))=',h(i,nk(i))
     4353!           write(*,*) 'lv(i,k)=',lv(i,k)
     4354!           write(*,*) 't(i,k)=',t(i,k)
     4355!           write(*,*) 'clw(i,k)=',clw(i,k)
     4356!           write(*,*) 'cpd,cpv=',cpd,cpv
     4357!           stop
     4358!        endif
     4359!       enddo !do k=1,nl
     4360!      enddo !do i=1,ncum 
     4361      endif !if (coef_epmax_cape.gt.1e-12) then
     4362      !write(*,*) 'cv3_routines check 4367'
     4363
     4364      return
     4365      end subroutine cv3_epmax_fn_cape
     4366
     4367
     4368
Note: See TracChangeset for help on using the changeset viewer.