Ignore:
Timestamp:
Apr 15, 2025, 11:56:45 AM (2 months ago)
Author:
aborella
Message:

Merge with trunk testing r5597. We have convergence in prod and debug in NPv7.0.1c

Location:
LMDZ6/branches/contrails
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • LMDZ6/branches/contrails

  • LMDZ6/branches/contrails/libf/phylmd/cv3_routines.f90

    r5536 r5618  
    1212  USE conema3_mod_h
    1313  USE lmdz_cv_ini, ONLY : alpha,alpha1,beta,betad,coef_peel,cv_flag_feed,delta,dpbase,dtcrit,dtovsh,dttrig,ejectice,ejectliq,elcrit,flag_epkeorig,flag_wb,minorig,nl,nlm,nlp,noconv_stop,noff,omtrain,pbcrit,ptcrit,sigdz,spfac,t_top_max,tau,tau_stop,tlcrit,wbmax
     14  USE lmdz_cv_ini, ONLY : keep_bug_indices_cv3_tracer, keep_bug_q_nocons_cv
    1415
    1516
     
    139140     keepbug_ice_frac = .TRUE.
    140141     CALL getin_p('keepbug_ice_frac', keepbug_ice_frac)
     142     keep_bug_indices_cv3_tracer = .FALSE.
     143     CALL getin_p('keep_bug_indices_cv3_tracer', keep_bug_indices_cv3_tracer)
     144     keep_bug_q_nocons_cv = .TRUE.
     145     CALL getin_p('keep_bug_q_nocons_cv', keep_bug_q_nocons_cv)
     146
    141147
    142148    WRITE (*, *) 't_top_max=', t_top_max
     
    164170    WRITE (*, *) 'adiab_ascent_mass_flux_depends_on_ejectliq =', adiab_ascent_mass_flux_depends_on_ejectliq
    165171    WRITE (*, *) 'keepbug_ice_frac =', keepbug_ice_frac
     172    WRITE (*, *) 'keep_bug_indices_cv3_tracer =', keep_bug_indices_cv3_tracer
     173    WRITE (*, *) 'keep_bug_q_nocons_cv =', keep_bug_q_nocons_cv
    166174
    167175    first = .FALSE.
     
    26992707                     wdtrainA, wdtrainS, wdtrainM)                                      ! RomP
    27002708  USE lmdz_cv_ini, ONLY : cpd,ginv,grav,nl,nlp,sigdz
     2709  USE lmdz_cv_ini, ONLY : keep_bug_q_nocons_cv
    27012710  USE cvflag_mod_h
    27022711  USE print_control_mod, ONLY: prt_level, lunout
     
    29012910
    29022911
    2903     DO il = 1, ncum
    2904       IF (i<=inb(il) .AND. lwork(il)) THEN
    2905         wdtrain(il) = grav*ep(il, i)*m(il, i)*clw(il, i)
    2906         wdtrainS(il, i) = wdtrain(il)/grav                                            !   Ps   jyg
    2907 !!        wdtrainA(il, i) = wdtrain(il)/grav                                          !   Ps   RomP
    2908       END IF
    2909     END DO
    2910 
    2911     IF (i>1) THEN
    2912       DO j = 1, i - 1
     2912  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2913  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2914  IF (keep_bug_q_nocons_cv) THEN
     2915  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    29132916        DO il = 1, ncum
    29142917          IF (i<=inb(il) .AND. lwork(il)) THEN
    2915             awat = elij(il, j, i) - (1.-ep(il,i))*clw(il, i)
    2916             awat = max(awat, 0.0)
    2917             wdtrain(il) = wdtrain(il) + grav*awat*ment(il, j, i)
    2918             wdtrainM(il, i) = wdtrain(il)/grav - wdtrainS(il, i)    !   Pm  jyg
    2919 !!            wdtrainM(il, i) = wdtrain(il)/grav - wdtrainA(il, i)  !   Pm  RomP
     2918            wdtrain(il) = grav*ep(il, i)*m(il, i)*clw(il, i)
     2919            wdtrainS(il, i) = wdtrain(il)/grav                                            !   Ps   jyg
    29202920          END IF
    29212921        END DO
    2922       END DO
    2923     END IF
    2924 
    2925     IF (cvflag_prec_eject) THEN
    2926 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2927       IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
    2928 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2929 !!! Warning : this option leads to water conservation violation
    2930 !!!           Expert only
    2931 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2932           IF ( i > 1) THEN
     2922   
     2923        IF (i>1) THEN
     2924          DO j = 1, i - 1
    29332925            DO il = 1, ncum
    29342926              IF (i<=inb(il) .AND. lwork(il)) THEN
    2935                 wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))/(1. - qta(il, i-1))    !   Pa   jygprl
    2936                 wdtrain(il) = wdtrain(il) + grav*wdtrainA(il,i)
     2927                awat = elij(il, j, i) - (1.-ep(il,i))*clw(il, i)
     2928                awat = max(awat, 0.0)
     2929                wdtrain(il) = wdtrain(il) + grav*awat*ment(il, j, i)
     2930                wdtrainM(il, i) = wdtrain(il)/grav - wdtrainS(il, i)    !   Pm  jyg
    29372931              END IF
    29382932            END DO
    2939           ENDIF  ! ( i > 1)
    2940 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2941       ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
    2942 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2943           IF ( i > 1) THEN
     2933          END DO
     2934        END IF
     2935   
     2936        IF (cvflag_prec_eject) THEN
     2937    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2938          IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
     2939    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2940    !!! Warning : this option leads to water conservation violation
     2941    !!!           Expert only
     2942    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2943              IF ( i > 1) THEN
     2944                DO il = 1, ncum
     2945                  IF (i<=inb(il) .AND. lwork(il)) THEN
     2946                    wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))/(1. - qta(il, i-1))    !   Pa   jygprl
     2947                    wdtrain(il) = wdtrain(il) + grav*wdtrainA(il,i)
     2948                  END IF
     2949                END DO
     2950              ENDIF  ! ( i > 1)
     2951    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2952          ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
     2953    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2954              IF ( i > 1) THEN
     2955                DO il = 1, ncum
     2956                  IF (i<=inb(il) .AND. lwork(il)) THEN
     2957                    wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))                        !   Pa   jygprl
     2958                    wdtrain(il) = wdtrain(il) + grav*wdtrainA(il,i)
     2959                  END IF
     2960                END DO
     2961              ENDIF  ! ( i > 1)
     2962   
     2963          ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
     2964    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2965        ENDIF  ! (cvflag_prec_eject)
     2966   
     2967  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2968  ELSE ! (keep_bug_q_nocons_cv)
     2969  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2970        DO il = 1, ncum
     2971          IF (i<=inb(il) .AND. lwork(il)) THEN
     2972            wdtrainS(il, i) = ep(il, i)*m(il, i)*clw(il, i)                               ! jyg
     2973          END IF
     2974        END DO
     2975   
     2976        IF (i>1) THEN
     2977          DO j = 1, i - 1
    29442978            DO il = 1, ncum
    29452979              IF (i<=inb(il) .AND. lwork(il)) THEN
    2946                 wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))                        !   Pa   jygprl
    2947                 wdtrain(il) = wdtrain(il) + grav*wdtrainA(il,i)
     2980                awat = elij(il, j, i) - (1.-ep(il,i))*clw(il, i)
     2981                awat = max(awat, 0.0)
     2982                wdtrainM(il, i) = wdtrainM(il, i) + awat*ment(il, j, i)                   ! jyg
    29482983              END IF
    29492984            END DO
    2950           ENDIF  ! ( i > 1)
    2951 
    2952       ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
    2953 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    2954     ENDIF  ! (cvflag_prec_eject)
    2955 
     2985          END DO
     2986        END IF
     2987   
     2988        IF (cvflag_prec_eject) THEN
     2989    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2990          IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
     2991    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2992    !!! Warning : this option leads to water conservation violation
     2993    !!!           Expert only
     2994    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2995              IF ( i > 1) THEN
     2996                DO il = 1, ncum
     2997                  IF (i<=inb(il) .AND. lwork(il)) THEN
     2998                    wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))/(1. - qta(il, i-1))    !   Pa   jygprl
     2999                  END IF
     3000                END DO
     3001              ENDIF  ! ( i > 1)
     3002    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3003          ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
     3004    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3005              IF ( i > 1) THEN
     3006                DO il = 1, ncum
     3007                  IF (i<=inb(il) .AND. lwork(il)) THEN
     3008                    wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))                        !   Pa   jygprl
     3009                  END IF
     3010                END DO
     3011              ENDIF  ! ( i > 1)
     3012   
     3013          ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
     3014    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3015        ENDIF  ! (cvflag_prec_eject)
     3016   
     3017        IF ( i > 1) THEN
     3018          DO il = 1, ncum
     3019            IF (i<=inb(il) .AND. lwork(il)) THEN
     3020              wdtrain(il) = grav*(wdtrainS(il,i) + wdtrainM(il,i) + wdtrainA(il,i))
     3021            END IF
     3022          END DO
     3023        ENDIF  ! ( i > 1)
     3024   
     3025  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3026  ENDIF ! (keep_bug_q_nocons_cv)
     3027  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3028  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    29563029
    29573030! ***    find rain water and evaporation using provisional   ***
     
    31353208          ice(il, i) = ice(il, i) - fondue(il, i)
    31363209
    3137           IF (water(il,i)+ice(il,i)<1.E-30) THEN
    3138             faci(il, i) = 0.
    3139           ELSE
    3140             faci(il, i) = ice(il, i)/(water(il,i)+ice(il,i))
    3141           END IF
     3210!!          IF (water(il,i)+ice(il,i)<1.E-30) THEN
     3211!!            faci(il, i) = 0.
     3212!!          ELSE
     3213!!            faci(il, i) = ice(il, i)/(water(il,i)+ice(il,i))
     3214!!          END IF
     3215
     3216            faci(il,i) = ice(il, i)/max((water(il,i)+ice(il,i)), smallestreal)
    31423217
    31433218!           water(il,i)=water(il,i+1)+(1.-fraci(il,i))*e6+(1.-faci(il,i))*f6
     
    34193494! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    34203495
    3421 
    34223496  RETURN
    34233497
     
    34453519    USE cvflag_mod_h
    34463520   USE lmdz_cv_ini, ONLY : grav,minorig,nl,nlp,rowl,rrd,nl,ci,cl,cpd,cpv
     3521   USE lmdz_cv_ini, ONLY : keep_bug_q_nocons_cv
    34473522  IMPLICIT NONE
    34483523
     
    35273602      REAL, DIMENSION (nloc, nd)                         :: sigment, qtment             ! cld
    35283603      REAL, DIMENSION (nloc, nd, nd)                     :: qdet
    3529       REAL sumdq !jyg
     3604!!      REAL sumdq !jyg
    35303605!
    35313606! -------------------------------------------------------------
     3607
    35323608
    35333609! initialization:
     
    40034079! ***                      through each level                          ***
    40044080
    4005 
    40064081!jyg<
    40074082!!  DO i = 2, nl + 1 ! newvecto: mettre nl au lieu nl+1?
     
    40204095           IF (ok_optim_yield) THEN                       !|
    40214096!-----------------------------------------------------------
    4022 DO il = 1, ncum
    4023    amp1(il) = upwd(il,i+1)
    4024    ad(il) = dnwd(il,i)
    4025 ENDDO
     4097    IF (keep_bug_q_nocons_cv) THEN    !!jyg20250215
     4098      DO il = 1, ncum
     4099         amp1(il) = upwd(il,i+1)
     4100         ad(il) = dnwd(il,i)
     4101      ENDDO
     4102    ELSE  ! (keep_bug_q_nocons_cv)
     4103      DO il = 1, ncum
     4104         amp1(il) = upwd(il,i+1)
     4105         ad(il) = - dnwd(il,i)
     4106      ENDDO
     4107    ENDIF  ! (keep_bug_q_nocons_cv)
    40264108!-----------------------------------------------------------
    40274109        ELSE !(ok_optim_yield)                            !|
     
    43564438500 END DO
    43574439
    4358 !JYG<
    4359 !Conservation de l'eau
    4360 !   sumdq = 0.
    4361 !   DO k = 1, nl
    4362 !     sumdq = sumdq + fr(1, k)*100.*(ph(1,k)-ph(1,k+1))/grav
    4363 !   END DO
    4364 !   PRINT *, 'cv3_yield, apres 500, sum(dq), precip, somme ', sumdq, Vprecip(1, 1), sumdq + vprecip(1, 1)
    4365 !JYG>
     4440!!!JYG<
     4441!!!Conservation de l'eau
     4442!!   sumdq = 0.
     4443!!   DO k = 1, nl
     4444!!     sumdq = sumdq + fr(1, k)*100.*(ph(1,k)-ph(1,k+1))/grav
     4445!!   END DO
     4446!!   PRINT *, 'cv3_yield, apres 500, sum(dq), precip, somme ', sumdq, Vprecip(1, 1), sumdq + vprecip(1, 1)
     4447!!!JYG>
    43664448! ***   move the detrainment at level inb down to level inb-1   ***
    43674449! ***        in such a way as to preserve the vertically        ***
     
    43984480  END DO
    43994481
    4400 !JYG<
    4401 !Conservation de l'eau
    4402 !   sumdq = 0.
    4403 !   DO k = 1, nl
    4404 !     sumdq = sumdq + fr(1, k)*100.*(ph(1,k)-ph(1,k+1))/grav
    4405 !   END DO
    4406 !   PRINT *, 'cv3_yield, apres 503, sum(dq), precip, somme ', sumdq, Vprecip(1, 1), sumdq + vprecip(1, 1)
    4407 !JYG>
     4482!!!JYG<
     4483!!!Conservation de l'eau
     4484!!   sumdq = 0.
     4485!!   DO k = 1, nl
     4486!!     sumdq = sumdq + fr(1, k)*100.*(ph(1,k)-ph(1,k+1))/grav
     4487!!   END DO
     4488!!   PRINT *, 'cv3_yield, apres 503, sum(dq), precip, somme ', sumdq, Vprecip(1, 1), sumdq + vprecip(1, 1)
     4489!!!JYG>
    44084490
    44094491!AC!      do j=1,ntra
     
    49365018                      ep, Vprecip, elij, clw, epmlmMm, eplaMm, &
    49375019                      icb, inb)
    4938    USE lmdz_cv_ini, ONLY : nl
     5020  USE lmdz_cv_ini, ONLY : nl,keep_bug_indices_cv3_tracer
    49395021  USE cvflag_mod_h
     5022  USE ioipsl_getin_p_mod, ONLY : getin_p
    49405023  IMPLICIT NONE
    49415024
    49425025
    49435026!inputs:
     5027!------
    49445028  INTEGER, INTENT (IN)                               :: ncum, nd, na, nloc, len
    49455029  INTEGER, DIMENSION (len), INTENT (IN)              :: icb, inb
     
    49495033  REAL, DIMENSION (len, nd+1), INTENT (IN)           :: Vprecip
    49505034!ouputs:
     5035!------
    49515036  REAL, DIMENSION (len, na, na), INTENT (OUT)        :: phi, phi2, epmlmMm
    49525037  REAL, DIMENSION (len, na), INTENT (OUT)            :: da, d1a, dam, eplaMm
    49535038!
     5039!local variables:
     5040!---------------
    49545041! variables pour tracer dans precip de l'AA et des mel
    4955 !local variables:
    49565042  INTEGER i, j, k
    49575043  REAL epm(nloc, na, na)
    4958 
     5044!
    49595045! variables d'Emanuel : du second indice au troisieme
    49605046! --->    tab(i,k,j) -> de l origine k a l arrivee j
     
    49625048! variables personnelles : du troisieme au second indice
    49635049! --->    tab(i,j,k) -> de k a j
    4964 ! phi, phi2
    4965 
    4966 ! initialisations
     5050! phi, phi2, epm, epmlmMm
     5051
    49675052
    49685053  da(:, :) = 0.
     
    50225107        d1a(i, j) = d1a(i, j) + ment(i, k, j)*ep(i, k)*(1.-sigij(i,k,j))
    50235108        IF (k<=j) THEN
    5024           dam(i, j) = dam(i, j) + ment(i, k, j)*epm(i, k, j)*(1.-ep(i,k))*(1.-sigij(i,k,j))
    50255109          phi2(i, j, k) = phi(i, j, k)*epm(i, j, k)
    50265110        END IF
     
    50285112    END DO
    50295113  END DO
     5114
     5115  IF (keep_bug_indices_cv3_tracer) THEN
     5116    DO j = 1, nl
     5117      DO k = 1, nl
     5118        DO i = 1, ncum
     5119          IF (k<=j) THEN
     5120            dam(i, j) = dam(i, j) + ment(i, k, j)*epm(i, k, j)*(1.-ep(i,k))*(1.-sigij(i,k,j))
     5121          END IF ! (k<=j)
     5122        END DO
     5123      END DO
     5124    END DO
     5125  ELSE  ! (keep_bug_indices_cv3_tracer)
     5126    DO j = 1, nl
     5127      DO k = 1, nl
     5128        DO i = 1, ncum
     5129          IF (k<=j) THEN
     5130            dam(i, j) = dam(i, j) + ment(i, k, j)*epm(i, j, k)*(1.-ep(i,k))*(1.-sigij(i,k,j))
     5131          END IF ! (k<=j)
     5132        END DO
     5133      END DO
     5134    END DO
     5135  ENDIF ! (keep_bug_indices_cv3_tracer)
    50305136
    50315137  RETURN
Note: See TracChangeset for help on using the changeset viewer.