Ignore:
Timestamp:
Jun 6, 2016, 4:04:57 PM (8 years ago)
Author:
Laurent Fairhead
Message:

Merged trunk changes r2487:2541 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

    r2488 r2542  
    105105     ok_intermittent=.False.
    106106     CALL getin_p('ok_intermittent',ok_intermittent)
     107     ok_optim_yield=.False.
     108     CALL getin_p('ok_optim_yield',ok_optim_yield)
    107109     coef_peel=0.25
    108110     CALL getin_p('coef_peel',coef_peel)
    … …  
    126128    WRITE (*, *) 'tau_stop=', tau_stop
    127129    WRITE (*, *) 'ok_intermittent=', ok_intermittent
     130    WRITE (*, *) 'ok_optim_yield =', ok_optim_yield
    128131    WRITE (*, *) 'coef_peel=', coef_peel
    129132
    … …  
    580583
    581584  DO i = 1, len                                           !convect3
    582     icb1(i) = max(icb(i), 2)                              !convect3
    583     icb1(i) = min(icb(i), nl)                             !convect3
     585    icb1(i) = min(max(icb(i), 2), nl)
    584586! if icb is below LCL, start loop at ICB+1:
    585587! (icbs est le premier niveau au-dessus du LCL)
    … …  
    26052607!
    26062608
    2607           d6 = bfac*wdtrain(il) - 100.*sigd(il)*bfac*(ph(il,i)-ph(il,i+1))*evap(il, i)
    2608           e6 = bfac*wdtrain(il)
    2609           f6 = -100.*sigd(il)*bfac*(ph(il,i)-ph(il,i+1))*evap(il, i)
    2610 
     2609!jyg<
     2610          d6 = prec(il,i)-prec(il,i+1)
     2611
     2612!!          d6 = bfac*wdtrain(il) - 100.*sigd(il)*bfac*(ph(il,i)-ph(il,i+1))*evap(il, i)
     2613!!          e6 = bfac*wdtrain(il)
     2614!!          f6 = -100.*sigd(il)*bfac*(ph(il,i)-ph(il,i+1))*evap(il, i)
     2615!>jyg
    26112616!CR:tmax_fonte_cv: T for which ice is totally melted (used to be 275.15)
    26122617          thaw = (t(il,i)-273.15)/(tmax_fonte_cv-273.15)
    26132618          thaw = min(max(thaw,0.0), 1.0)
     2619!jyg<
    26142620          water(il, i) = water(il, i+1) + (1-fraci(il,i))*d6
    2615           water(il, i) = max(water(il,i), 0.)
    2616           ice(il, i) = ice(il, i+1) + fraci(il, i)*d6
    2617           ice(il, i) = max(ice(il,i), 0.)
     2621          ice(il, i)   = ice(il, i+1)   + fraci(il, i)*d6
     2622          water(il, i) = min(prec(il,i), max(water(il,i), 0.))
     2623          ice(il, i)   = min(prec(il,i), max(ice(il,i),   0.))
     2624
     2625!!          water(il, i) = water(il, i+1) + (1-fraci(il,i))*d6
     2626!!          water(il, i) = max(water(il,i), 0.)
     2627!!          ice(il, i) = ice(il, i+1) + fraci(il, i)*d6
     2628!!          ice(il, i) = max(ice(il,i), 0.)
     2629!>jyg
    26182630          fondue(il, i) = ice(il, i)*thaw
    26192631          water(il, i) = water(il, i) + fondue(il, i)
    … …  
    29672979!
    29682980!local variables:
    2969       INTEGER i, k, il, n, j, num1
    2970       REAL rat, delti
    2971       REAL ax, bx, cx, dx, ex
    2972       REAL cpinv, rdcp, dpinv
    2973       REAL awat(nloc)
    2974       REAL lvcp(nloc, na), lfcp(nloc, na)                  ! , mke(nloc, na) ! unused . jyg
    2975       REAL am(nloc), work(nloc), ad(nloc), amp1(nloc)
     2981      INTEGER                                            :: i, k, il, n, j, num1
     2982      REAL                                               :: rat, delti
     2983      REAL                                               :: ax, bx, cx, dx, ex
     2984      REAL                                               :: cpinv, rdcp, dpinv
     2985      REAL, DIMENSION (nloc)                             ::  awat
     2986      REAL, DIMENSION (nloc, nd)                         :: lvcp, lfcp              ! , mke ! unused . jyg
     2987      REAL, DIMENSION (nloc)                             :: am, work, ad, amp1
    29762988!!      real up1(nloc), dn1(nloc)
    2977       REAL up1(nloc, nd, nd), dn1(nloc, nd, nd)
    2978       REAL asum(nloc), bsum(nloc), csum(nloc), dsum(nloc)
    2979       REAL esum(nloc), fsum(nloc), gsum(nloc), hsum(nloc)
    2980       REAL th_wake(nloc, nd)
    2981       REAL alpha_qpos(nloc), alpha_qpos1(nloc)
    2982       REAL qcond(nloc, nd), nqcond(nloc, nd), wa(nloc, nd) ! cld
    2983       REAL siga(nloc, nd), sax(nloc, nd), mac(nloc, nd) ! cld
    2984       REAL sument(nloc), sigment(nloc,nd), qtment(nloc,nd) ! cld
    2985       REAL qnk(nloc)
     2989      REAL, DIMENSION (nloc, nd, nd)                     :: up1, dn1
     2990!jyg<
     2991      REAL, DIMENSION (nloc, nd)                         :: up_to, up_from
     2992      REAL, DIMENSION (nloc, nd)                         :: dn_to, dn_from
     2993!>jyg
     2994      REAL, DIMENSION (nloc)                             :: asum, bsum, csum, dsum
     2995      REAL, DIMENSION (nloc)                             :: esum, fsum, gsum, hsum
     2996      REAL, DIMENSION (nloc, nd)                         :: th_wake
     2997      REAL, DIMENSION (nloc)                             :: alpha_qpos, alpha_qpos1
     2998      REAL, DIMENSION (nloc, nd)                         :: qcond, nqcond, wa           ! cld
     2999      REAL, DIMENSION (nloc, nd)                         :: siga, sax, mac              ! cld
     3000      REAL, DIMENSION (nloc)                             :: sument
     3001      REAL, DIMENSION (nloc, nd)                         :: sigment, qtment             ! cld
     3002      REAL, DIMENSION (nloc)                             :: qnk
    29863003      REAL sumdq !jyg
    29873004!
    … …  
    32233240! print*,'cv3_yield apres ft'
    32243241
     3242!jyg<
     3243!-----------------------------------------------------------
     3244           IF (ok_optim_yield) THEN                       !|
     3245!-----------------------------------------------------------
     3246!
     3247!***                                                      ***
     3248!***    Compute convective mass fluxes upwd and dnwd      ***
     3249
     3250upwd(:,:) = 0.
     3251up_to(:,:) = 0.
     3252up_from(:,:) = 0.
     3253dnwd(:,:) = 0.
     3254dn_to(:,:) = 0.
     3255dn_from(:,:) = 0.
     3256!
     3257! =================================================
     3258!              upward fluxes                      |
     3259! ------------------------------------------------
     3260DO i = 2, nl
     3261  DO il = 1, ncum
     3262    IF (i<=inb(il)) THEN
     3263      up_to(il,i) = m(il,i)
     3264    ENDIF
     3265  ENDDO
     3266  DO j = 1, i-1
     3267    DO il = 1, ncum
     3268      IF (i<=inb(il)) THEN
     3269        up_to(il,i) = up_to(il,i) + ment(il,j,i)
     3270      ENDIF
     3271    ENDDO
     3272  ENDDO
     3273ENDDO
     3274!
     3275DO i = 1, nl
     3276  DO il = 1, ncum
     3277    IF (i<=inb(il)) THEN
     3278      up_from(il,i) = cbmf(il)*wghti(il,i)
     3279    ENDIF
     3280  ENDDO
     3281ENDDO
     3282!!DO i = 2, nl
     3283!!  DO j = i+1, nl          !! Permuter les boucles i et j
     3284DO j = 3, nl
     3285  DO i = 2, j-1
     3286    DO il = 1, ncum
     3287      IF (j<=inb(il)) THEN
     3288        up_from(il,i) = up_from(il,i) + ment(il,i,j)
     3289      ENDIF
     3290    ENDDO
     3291  ENDDO
     3292ENDDO
     3293!
     3294! The difference between upwd(il,i) and upwd(il,i-1) is due to updrafts ending in layer
     3295!(i-1) (theses drafts cross interface (i-1) but not interface(i)) and to updrafts starting
     3296!from layer (i-1) (theses drafts cross interface (i) but not interface(i-1)):
     3297!
     3298DO i = 2, nlp
     3299  DO il = 1, ncum
     3300    upwd(il,i) = max(0., upwd(il,i-1) - up_to(il,i-1) + up_from(il,i-1))
     3301  ENDDO
     3302ENDDO
     3303!
     3304! =================================================
     3305!              downward fluxes                    |
     3306! ------------------------------------------------
     3307DO i = 1, nl
     3308  DO j = i+1, nl
     3309    DO il = 1, ncum
     3310      IF (j<=inb(il)) THEN
     3311        dn_to(il,i) = dn_to(il,i) + ment(il,j,i)
     3312      ENDIF
     3313    ENDDO
     3314  ENDDO
     3315ENDDO
     3316!
     3317!!DO i = 2, nl
     3318!!  DO j = 1, i-1          !! Permuter les boucles i et j
     3319DO j = 1, nl
     3320  DO i = j+1, nl
     3321    DO il = 1, ncum
     3322      IF (i<=inb(il)) THEN
     3323        dn_from(il,i) = dn_from(il,i) + ment(il,i,j)
     3324      ENDIF
     3325    ENDDO
     3326  ENDDO
     3327ENDDO
     3328!
     3329! The difference between dnwd(il,i) and dnwd(il,i+1) is due to downdrafts ending in layer
     3330!(i) (theses drafts cross interface (i+1) but not interface(i)) and to downdrafts
     3331!starting from layer (i) (theses drafts cross interface (i) but not interface(i+1)):
     3332!
     3333DO i = nl-1, 1, -1
     3334  DO il = 1, ncum
     3335    dnwd(il,i) = max(0., dnwd(il,i+1) - dn_to(il,i) + dn_from(il,i))
     3336  ENDDO
     3337ENDDO
     3338! =================================================
     3339!
     3340!-----------------------------------------------------------
     3341        ENDIF !(ok_optim_yield)                           !|
     3342!-----------------------------------------------------------
     3343!>jyg
     3344
    32253345! ***  calculate tendencies of potential temperature and mixing ratio  ***
    32263346! ***               at levels above the lowest level                   ***
    … …  
    32303350
    32313351
    3232   DO i = 2, nl + 1 ! newvecto: mettre nl au lieu nl+1?
     3352!jyg<
     3353!!  DO i = 2, nl + 1 ! newvecto: mettre nl au lieu nl+1?
     3354  DO i = 2, nl
     3355!>jyg
    32333356
    32343357    num1 = 0
    … …  
    32383361    IF (num1<=0) GO TO 500
    32393362
     3363!
    32403364!jyg<
    3241 !!    CALL zilch(amp1, ncum)
    3242 !!    CALL zilch(ad, ncum)
     3365!-----------------------------------------------------------
     3366           IF (ok_optim_yield) THEN                       !|
     3367!-----------------------------------------------------------
     3368DO il = 1, ncum
     3369   amp1(il) = upwd(il,i+1)
     3370   ad(il) = dnwd(il,i)
     3371ENDDO
     3372!-----------------------------------------------------------
     3373        ELSE !(ok_optim_yield)                            !|
     3374!-----------------------------------------------------------
     3375!>jyg
    32433376    DO il = 1,ncum
    32443377      amp1(il) = 0.
    32453378      ad(il) = 0.
    32463379    ENDDO
    3247 !>jyg
    32483380
    32493381    DO k = 1, nl + 1
    … …  
    32623394    END DO
    32633395
    3264     DO k = 1, i
    3265       DO j = i + 1, nl + 1
     3396    DO j = i + 1, nl + 1         
     3397       DO k = 1, i
     3398          !yor! reverted j and k loops
     3399          DO il = 1, ncum
     3400!yor!        IF (i<=inb(il) .AND. j<=(inb(il)+1)) THEN ! the second condition implies the first !
     3401             IF (j<=(inb(il)+1)) THEN 
     3402                amp1(il) = amp1(il) + ment(il, k, j)
     3403             END IF
     3404          END DO
     3405       END DO
     3406    END DO
     3407
     3408    DO k = 1, i - 1
     3409!jyg<
     3410!!      DO j = i, nl + 1 ! newvecto: nl au lieu nl+1?
     3411      DO j = i, nl
     3412!>jyg
    32663413        DO il = 1, ncum
    3267           IF (i<=inb(il) .AND. j<=(inb(il)+1)) THEN
    3268             amp1(il) = amp1(il) + ment(il, k, j)
    3269           END IF
    3270         END DO
    3271       END DO
    3272     END DO
    3273 
    3274     DO k = 1, i - 1
    3275       DO j = i, nl + 1 ! newvecto: nl au lieu nl+1?
    3276         DO il = 1, ncum
    3277           IF (i<=inb(il) .AND. j<=inb(il)) THEN
     3414!yor!        IF (i<=inb(il) .AND. j<=inb(il)) THEN ! the second condition implies the 1st !
     3415             IF (j<=inb(il)) THEN   
    32783416            ad(il) = ad(il) + ment(il, j, k)
    32793417          END IF
    … …  
    32813419      END DO
    32823420    END DO
     3421!
     3422!-----------------------------------------------------------
     3423        ENDIF !(ok_optim_yield)                           !|
     3424!-----------------------------------------------------------
     3425!
     3426!!   print *,'yield, i, amp1, ad', i, amp1(1), ad(1)
    32833427
    32843428    DO il = 1, ncum
    … …  
    34283572!AC!      enddo
    34293573
    3430     DO k = i, nl + 1
     3574!jyg<
     3575!!    DO k = i, nl + 1
     3576    DO k = i, nl
     3577!>jyg
    34313578
    34323579      IF (iflag_mix/=0) THEN
    … …  
    37313878    END DO
    37323879  END DO
     3880!jyg<
     3881!-----------------------------------------------------------
     3882           IF (ok_optim_yield) THEN                       !|
     3883!-----------------------------------------------------------
    37333884  DO i = 1, nl
    3734     DO j = 1, nl
     3885    DO il = 1, ncum
     3886      IF (iflag(il)<=1) THEN
     3887        upwd(il, i) = upwd(il, i)/alpha_qpos(il)
     3888        dnwd(il, i) = dnwd(il, i)/alpha_qpos(il)
     3889      END IF
     3890    END DO
     3891  END DO
     3892!-----------------------------------------------------------
     3893        ENDIF !(ok_optim_yield)                           !|
     3894!-----------------------------------------------------------
     3895!>jyg
     3896  DO j = 1, nl !yor! inverted i and j loops
     3897     DO i = 1, nl
    37353898      DO il = 1, ncum
    37363899        IF (iflag(il)<=1) THEN
    … …  
    37663929  DO i = 1, nl
    37673930    DO il = 1, ncum
    3768       upwd(il, i) = 0.0
    3769       dnwd(il, i) = 0.0
    3770     END DO
    3771   END DO
    3772 
    3773   DO i = 1, nl
    3774     DO il = 1, ncum
    37753931      dnwd0(il, i) = -mp(il, i)
    37763932    END DO
    … …  
    37853941
    37863942
     3943!jyg<
     3944!-----------------------------------------------------------
     3945           IF (.NOT.ok_optim_yield) THEN                  !|
     3946!-----------------------------------------------------------
    37873947  DO i = 1, nl
    37883948    DO il = 1, ncum
    3789       IF (i>=icb(il) .AND. i<=inb(il)) THEN
    3790         upwd(il, i) = 0.0
    3791         dnwd(il, i) = 0.0
    3792       END IF
    3793     END DO
    3794   END DO
     3949      upwd(il, i) = 0.0
     3950      dnwd(il, i) = 0.0
     3951    END DO
     3952  END DO
     3953
     3954!!  DO i = 1, nl                                           ! useless; jyg
     3955!!    DO il = 1, ncum                                      ! useless; jyg
     3956!!      IF (i>=icb(il) .AND. i<=inb(il)) THEN              ! useless; jyg
     3957!!        upwd(il, i) = 0.0                                ! useless; jyg
     3958!!        dnwd(il, i) = 0.0                                ! useless; jyg
     3959!!      END IF                                             ! useless; jyg
     3960!!    END DO                                               ! useless; jyg
     3961!!  END DO                                                 ! useless; jyg
    37953962
    37963963  DO i = 1, nl
    … …  
    38033970  END DO
    38043971
     3972!yor! commented original
     3973!  DO i = 1, nl
     3974!    DO k = i, nl
     3975!      DO n = 1, i - 1
     3976!        DO il = 1, ncum
     3977!          IF (i>=icb(il) .AND. i<=inb(il) .AND. k<=inb(il)) THEN
     3978!            up1(il, k, i) = up1(il, k, i) + ment(il, n, k)
     3979!            dn1(il, k, i) = dn1(il, k, i) - ment(il, k, n)
     3980!          END IF
     3981!        END DO
     3982!      END DO
     3983!    END DO
     3984!  END DO
     3985!yor! replaced with
    38053986  DO i = 1, nl
    38063987    DO k = i, nl
    38073988      DO n = 1, i - 1
    38083989        DO il = 1, ncum
    3809           IF (i>=icb(il) .AND. i<=inb(il) .AND. k<=inb(il)) THEN
    3810             up1(il, k, i) = up1(il, k, i) + ment(il, n, k)
    3811             dn1(il, k, i) = dn1(il, k, i) - ment(il, k, n)
     3990          IF (i>=icb(il) .AND. k<=inb(il)) THEN ! yor ! as i always <= k
     3991             up1(il, k, i) = up1(il, k, i) + ment(il, n, k)
    38123992          END IF
    38133993        END DO
    … …  
    38153995    END DO
    38163996  END DO
     3997  DO i = 1, nl
     3998    DO n = 1, i - 1
     3999      DO k = i, nl
     4000        DO il = 1, ncum
     4001          IF (i>=icb(il) .AND. k<=inb(il)) THEN ! yor !  i always <= k
     4002             dn1(il, k, i) = dn1(il, k, i) - ment(il, k, n)
     4003          END IF
     4004        END DO
     4005      END DO
     4006    END DO
     4007  END DO
     4008!yor! end replace
    38174009
    38184010  DO i = 1, nl
    … …  
    38654057!!!!
    38664058!!!!      ENDDO
     4059!-----------------------------------------------------------
     4060        ENDIF !(.NOT.ok_optim_yield)                      !|
     4061!-----------------------------------------------------------
     4062!>jyg
    38674063
    38684064! ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
Note: See TracChangeset for help on using the changeset viewer.