Changeset 142 in lmdz_wrf


Ignore:
Timestamp:
Jul 29, 2014, 12:19:12 PM (10 years ago)
Author:
lfita
Message:

Incorporing all the work done in LMDZ_WRFmeas about:

1.- Chekcing for NaNs?
2.- Instabilities in thermcell (negative SQRT & robustness IF)

Location:
trunk/WRFV3/lmdz
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • trunk/WRFV3/lmdz/calltherm.F90

    r1 r142  
    126126      endif
    127127
     128! L. Fita, LMD July 2014. Initializing tendencies and variables
     129      d_t_the = 0.
     130      d_q_the = 0.
     131      d_u_the = 0.
     132      d_v_the = 0.
     133      Alp = 0.
     134      Ale = 0.
     135
    128136! Incrementer le compteur de la physique
    129137     itap   = itap + 1
     
    132140!  ===================
    133141!         print*,'thermiques: WARNING on passe t au lieu de t_seri'
    134 
    135142
    136143! On prend comme valeur initiale des thermiques la valeur du pas
     
    237244
    238245!            print*,'THERM iflag_thermas_ed=',iflag_thermals_ed
     246
    239247            CALL thermcell_main(itap,klon,klev,zdt  &
    240248     &      ,pplay,paprs,pphi,debut  &
     
    270278! Ce serait bien de changer, mai en prenant le temps de vérifier que ca
    271279! fait bien ce qu'on croit.
    272 
    273280       flag_bidouille_stratocu=iflag_thermals<=12.or.iflag_thermals==14.or.iflag_thermals==16.or.iflag_thermals==18
    274281
     
    300307      ENDDO
    301308       fm_therm(:,klev+1)=0.
    302 
    303 
    304309
    305310!  accumulation de la tendance
  • trunk/WRFV3/lmdz/diagphy_mod.F90

    r31 r142  
    7979!C
    8080      integer pas
    81 ! L. Fita, LMD July 2014
    82       CHARACTER(LEN=50)                                  :: errmsg, fname
    83 
    84       fname = 'diagphy'
    85       errmsg = 'ERROR -- error -- ERROR -- error'
    8681
    8782      save pas
     
    8984!$OMP THREADPRIVATE(pas)
    9085!C
     86
     87! L. Fita, LMD July 2014
     88      CHARACTER(LEN=50)                                  :: errmsg, fname, varname
     89      LOGICAL                                            :: found
     90      REAL                                               :: largest
     91
     92      fname = 'diagphy'
     93      errmsg = 'ERROR -- error -- ERROR -- error'
     94      largest = 10.e4
     95
    9196      pas=pas+1
    9297      stops=0.
     
    153158
    154159! L. Fita, LMD July 2014. Checking for consistency
    155       IF (fs_bound .NE. fs_bound .OR. ABS(fs_bound) > 100000.) THEN
    156         PRINT *,TRIM(errmsg)
    157         PRINT *,'  ' + TRIM(fname) + ': Wrong fs_bound= ',fs_bound,' !!!'
     160      IF (fs_bound .NE. fs_bound .OR. ABS(fs_bound) > largest) THEN
     161        PRINT *,TRIM(errmsg)
     162        PRINT *,'  ' // TRIM(fname) // ': Wrong fs_bound= ',fs_bound,' !!!'
    158163        PRINT *,'    fs_bound: Heat flux at atm. boundaries'
    159164        PRINT *,'    fs_bound = stops-stopl - (ssols+ssoll)+ssens+sfront'
    160         PRINT *,'    stops-stopl= ',stops-stopl
     165        PRINT *,'    airetot= ',airetot
     166        IF (airetot .NE. airetot .OR. ABS(airetot) > largest*10.e6) THEN
     167          varname = 'airephy'
     168          CALL check_var(fname, varname, airephy, klon, largest*10.e6, .FALSE.)
     169        END IF
     170        PRINT *,'    stops= ',stops
     171        IF (stops .NE. stops .OR. ABS(stops) > largest) THEN
     172          varname = 'tops'
     173          CALL check_var(fname, varname, tops, klon, largest, .FALSE.)
     174        END IF
     175        PRINT *,'    stopl= ',stopl
     176        IF (stopl .NE. stopl .OR. ABS(stopl) > largest) THEN
     177          varname = 'topl'
     178          CALL check_var(fname, varname, topl, klon, largest, .FALSE.)
     179        END IF
    161180        PRINT *,'    ssols= ',ssols
     181        IF (ssols .NE. ssols .OR. ABS(ssols) > largest) THEN
     182          varname = 'sols'
     183          CALL check_var(fname, varname, sols, klon, largest, .FALSE.)
     184        END IF
    162185        PRINT *,'    ssoll= ',ssoll
     186        IF (ssoll .NE. ssoll .OR. ABS(ssoll) > largest) THEN
     187          varname = 'soll'
     188          CALL check_var(fname, varname, soll, klon, largest, .FALSE.)
     189        END IF
    163190        PRINT *,'    ssens= ',ssens
     191        IF (ssens .NE. ssens .OR. ABS(ssens) > largest) THEN
     192          varname = 'sens'
     193          CALL check_var(fname, varname, sens, klon, largest, .FALSE.)
     194        END IF
    164195        PRINT *,'    sfront= ',sfront
    165         STOP
    166       END IF
    167 
    168       IF (fq_bound .NE. fq_bound .OR. ABS(fq_bound) > 100000.) THEN
    169         PRINT *,TRIM(errmsg)
    170         PRINT *,'  ' + TRIM(fname) + ': Wrong fs_bound= ',fs_bound,' !!!'
     196        IF (sfront .NE. sfront .OR. ABS(sfront) > largest) THEN
     197          varname = 'evap'
     198          CALL check_var(fname, varname, evap, klon, largest, .FALSE.)
     199          varname = 'rain_fall'
     200          CALL check_var(fname, varname, rain_fall, klon, largest, .FALSE.)
     201          varname = 'snow_fall'
     202          CALL check_var(fname, varname, snow_fall, klon, largest, .FALSE.)
     203          varname = 'ts'
     204          CALL check_var(fname, varname, ts, klon, largest, .FALSE.)
     205        END IF
     206        STOP
     207      END IF
     208
     209      IF (d_etp_tot .NE. d_etp_tot .OR. ABS(d_etp_tot) > largest) THEN
     210        PRINT *,TRIM(errmsg)
     211        PRINT *,'  ' // TRIM(fname) // ': Wrong d_etp_tot= ',d_etp_tot,' !!!'
     212        PRINT *,'    d_etp_tot: Heat flux equivalent to atmospheric enthalpy' //     &
     213          ' change (W/m2)'
     214        PRINT *,'    d_etp_tot = input value !'
     215        STOP
     216      END IF
     217
     218      IF (fq_bound .NE. fq_bound .OR. ABS(fq_bound) > largest) THEN
     219        PRINT *,TRIM(errmsg)
     220        PRINT *,'  ' // TRIM(fname) // ': Wrong fq_bound= ',fs_bound,' !!!'
    171221        PRINT *,'    fq_bound: Watter flux at atm. boundaries'
    172222        PRINT *,'    fq_bound = evap_tot - rain_fall_tot -snow_fall_tot'
    173223        PRINT *,'    evap_tot= ',evap_tot
     224        PRINT *,'    airetot= ',airetot
     225        IF (airetot .NE. airetot .OR. ABS(airetot) > largest*10.e5) THEN
     226          varname = 'airephy'
     227          CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.)
     228        END IF
     229        IF (evap_tot .NE. evap_tot .OR. ABS(evap_tot) > largest) THEN
     230          varname = 'evap'
     231          CALL check_var(fname, varname, evap, klon, largest, .FALSE.)
     232        END IF
    174233        PRINT *,'    rain_fall_tot= ',rain_fall_tot
     234        IF (rain_fall_tot .NE. rain_fall_tot .OR. ABS(rain_fall_tot) > largest) THEN
     235          varname = 'rain_fall'
     236          CALL check_var(fname, varname, rain_fall, klon, largest, .FALSE.)
     237        END IF
    175238        PRINT *,'    snow_fall_tot= ',snow_fall_tot
     239        IF (snow_fall_tot .NE. snow_fall_tot .OR. ABS(snow_fall_tot) > largest) THEN
     240          varname = 'snow_fall'
     241          CALL check_var(fname, varname, snow_fall, klon, largest, .FALSE.)
     242        END IF
     243        STOP
     244      END IF
     245
     246      IF (d_qt_tot .NE. d_qt_tot .OR. ABS(d_qt_tot) > largest) THEN
     247        PRINT *,TRIM(errmsg)
     248        PRINT *,'  ' // TRIM(fname) // ': Wrong d_qt_tot= ',d_qt_tot,' !!!'
     249        PRINT *,'    d_qt_tot: Mass flux equivalent to atmospheric watter mass' //   &
     250          ' change (kg/m2/s)'
     251        PRINT *,'    d_qt_tot = input value !'
    176252        STOP
    177253      END IF
     
    184260
    185261      end SUBROUTINE diagphy
     262
     263!L. Fita, LMD 2004. Check variable
     264SUBROUTINE check_var(funcn, varn, var, sizev, bigvalue, stoprun)
     265!  Subroutine to check the consistency of a variable
     266!    * NaN value: by definition is variable /= variable
     267!    * bigvalue: threshold for the variable
     268
     269  IMPLICIT NONE
     270
     271#include "dimensions.h"
     272
     273  INTEGER, INTENT(IN)                                    :: sizev
     274  CHARACTER(LEN=50), INTENT(IN)                          :: funcn, varn
     275  REAL, DIMENSION(sizev), INTENT(IN)                     :: var
     276  REAL, INTENT(IN)                                       :: bigvalue
     277  LOGICAL, INTENT(IN)                                    :: stoprun
     278
     279! Local
     280  INTEGER                                                :: i, wrongi, xpt, ypt
     281  CHARACTER(LEN=50)                                      :: errmsg
     282  LOGICAL                                                :: found
     283  REAL, DIMENSION(sizev)                                 :: wrongvalues
     284  INTEGER, DIMENSION(sizev)                              :: wronggridpt
     285
     286!!!!!!! Variables
     287! funcn: at which functino of part of the program variable is checked
     288! varn: name of the variable
     289! var: variable to check
     290! sizev: size of the variable
     291! bigvalue: biggest attenaible value for the variable
     292! stoprun: Should the run stop if it founds a problem?
     293
     294  errmsg = 'ERROR -- error -- ERROR -- error'
     295
     296  found = .FALSE.
     297  wrongi = 0
     298  DO i=1,sizev
     299    IF (var(i) /= var(i) .OR. ABS(var(i)) > bigvalue ) THEN
     300      IF (wrongi == 0) found = .TRUE.
     301      wrongi = wrongi + 1
     302      wrongvalues(wrongi) = var(i)
     303      wronggridpt(wrongi) = i
     304    END IF
     305  END DO
     306
     307  IF (found) THEN
     308    PRINT *,TRIM(errmsg)
     309    PRINT *,"  at '" // TRIM(funcn) // "' variable '" //TRIM(varn)//                 &
     310      "' is wrong in Nvalues= ",wrongi,' at i (x, y) value___'
     311    DO i=1,wrongi
     312       ypt = INT(wronggridpt(i)/wiim) + 1
     313       xpt = wronggridpt(i) - (ypt-1)*wiim
     314      PRINT *,wronggridpt(i), '(',xpt,', ',ypt,')', wrongvalues(i)
     315    END DO
     316    IF (stoprun) THEN
     317      STOP
     318    END IF
     319  END IF
     320
     321  RETURN
     322
     323END SUBROUTINE check_var
     324
     325SUBROUTINE check_var3D(funcn, varn, var, sizev, zsize, bigvalue, stoprun)
     326!  Subroutine to check the consistency of a 3D LMDSZ - variable (klon, klev) !
     327!    * NaN value: by definition is variable /= variable
     328!    * bigvalue: threshold for the variable
     329
     330  IMPLICIT NONE
     331
     332#include "dimensions.h"
     333
     334  INTEGER, INTENT(IN)                                    :: sizev, zsize
     335  CHARACTER(LEN=50), INTENT(IN)                          :: funcn, varn
     336  REAL, DIMENSION(sizev,zsize), INTENT(IN)               :: var
     337  REAL, INTENT(IN)                                       :: bigvalue
     338  LOGICAL, INTENT(IN)                                    :: stoprun
     339
     340! Local
     341  INTEGER                                                :: i, k, wrongi, xpt, ypt
     342  CHARACTER(LEN=50)                                      :: errmsg
     343  LOGICAL                                                :: found
     344  REAL, DIMENSION(sizev*zsize)                           :: wrongvalues
     345  INTEGER, DIMENSION(sizev*zsize,2)                      :: wronggridpt
     346
     347!!!!!!! Variables
     348! funcn: at which functino of part of the program variable is checked
     349! varn: name of the variable
     350! var: variable to check
     351! sizev: size of the variable
     352! zsize: vertical size of the variable
     353! bigvalue: biggest attenaible value for the variable
     354! stoprun: Should the run stop if it founds a problem?
     355
     356  errmsg = 'ERROR -- error -- ERROR -- error'
     357
     358  found = .FALSE.
     359  wrongi = 0
     360  DO i=1,sizev
     361    DO k=1,zsize
     362      IF (var(i,k) /= var(i,k) .OR. ABS(var(i,k)) > bigvalue ) THEN
     363        IF (wrongi == 0) found = .TRUE.
     364        wrongi = wrongi + 1
     365        wrongvalues(wrongi) = var(i,k)
     366        wronggridpt(wrongi,1) = i
     367        wronggridpt(wrongi,2) = k
     368      END IF
     369    END DO
     370  END DO
     371
     372  IF (found) THEN
     373    PRINT *,TRIM(errmsg)
     374    PRINT *,"  at '" // TRIM(funcn) // "' variable '" //TRIM(varn)//                 &
     375      "' is wrong in Nvalues= ",wrongi,' at i (x,y) k value___'
     376    DO i=1,wrongi
     377       ypt = INT(wronggridpt(i,1)/wiim) + 1
     378       xpt = wronggridpt(i,1) - (ypt-1)*wiim
     379      PRINT *,wronggridpt(i,1), '(',xpt,', ',ypt,')', wronggridpt(i,2), wrongvalues(i)
     380    END DO
     381    IF (stoprun) THEN
     382      STOP
     383    END IF
     384  END IF
     385
     386  RETURN
     387
     388END SUBROUTINE check_var3D
     389
    186390
    187391!C======================================================================
     
    307511!$OMP THREADPRIVATE(h_vcol_pre, h_dair_pre, h_qw_pre, h_ql_pre)
    308512!$OMP THREADPRIVATE(h_qs_pre, qw_pre, ql_pre, qs_pre , ec_pre)
     513
     514! L. Fita, LMD July 2014
     515      CHARACTER(LEN=50)                                  :: errmsg, fname, varname
     516      LOGICAL                                            :: found
     517      REAL                                               :: largest
     518
     519      fname = 'diagetpq'
     520      errmsg = 'ERROR -- error -- ERROR -- error'
     521      largest = 10.e4
    309522!c======================================================================
    310523!C
     
    449662      ec_pre (idiag)    = ec_tot
    450663!C
     664      IF (d_h_vcol .NE. d_h_vcol .OR. ABS(d_h_vcol) > largest) THEN
     665        PRINT *,TRIM(errmsg)
     666        PRINT *,'  ' // TRIM(fname) // ': Wrong d_h_vcol= ',d_h_vcol,' !!!'
     667        PRINT *,'    d_h_vcol: Heat flux (W/m2) define as the Enthalpy change' //    &
     668          ' (J/m2) during one time step (dtime) for the whole atmosphere (air,' //   &
     669          ' watter vapour, liquid and solid)'
     670        PRINT *,'    d_h_vcol = (h_vcol_tot - h_vcol_pre(idiag2) )'
     671        PRINT *,'    h_vcol_tot= ',h_vcol_tot
     672        IF (h_vcol_tot .NE. h_vcol_tot .OR. ABS(h_vcol_tot) > largest) THEN
     673          PRINT *,'      h_vcol_tot = h_dair_tot+h_qw_tot+h_ql_tot+h_qs_tot'
     674          PRINT *,'    airetot= ',airetot
     675          IF (airetot .NE. airetot .OR. ABS(airetot) > largest*10.e5) THEN
     676            varname = 'airephy'
     677            CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.)
     678          END IF
     679          PRINT *,'    h_dair_tot= ',h_dair_tot
     680          IF (h_dair_tot .NE. h_dair_tot .OR. ABS(h_dair_tot) > largest*10.e6) THEN
     681            varname = 'zh_dair_col'
     682            CALL check_var(fname, varname, zh_dair_col, klon, largest*10.e6, .FALSE.)
     683            PRINT *,'      zh_dair_col = func(q,ql,qs,paprs,t)'
     684            varname = 'q'
     685            CALL check_var3D(fname, varname, q, klon, klev, largest, .FALSE.)
     686            varname = 'ql'
     687            CALL check_var3D(fname, varname, ql, klon, klev, largest, .FALSE.)
     688            varname = 'qs'
     689            CALL check_var3D(fname, varname, qs, klon, klev, largest, .FALSE.)
     690            varname = 'paprs'
     691            CALL check_var3D(fname, varname, paprs, klon, klev, largest*10.e6, .FALSE.)
     692            varname = 't'
     693            CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.)
     694          END IF
     695          PRINT *,'    h_qw_tot= ',h_qw_tot
     696          IF (h_qw_tot .NE. h_qw_tot .OR. ABS(h_qw_tot) > largest*10.e4) THEN
     697            varname = 'zh_qw_col'
     698            CALL check_var(fname, varname, zh_qw_col, klon, largest*10.e4, .FALSE.)
     699          END IF
     700          PRINT *,'    h_ql_tot= ',h_ql_tot
     701          IF (h_ql_tot .NE. h_ql_tot .OR. ABS(h_ql_tot) > largest*10.e2) THEN
     702            varname = 'zh_ql_col'
     703            CALL check_var(fname, varname, zh_ql_col, klon, largest*10.e2, .FALSE.)
     704          END IF
     705          PRINT *,'    h_qs_tot= ',h_qs_tot
     706          IF (h_qs_tot .NE. h_qs_tot .OR. ABS(h_qs_tot) > largest) THEN
     707            varname = 'zh_qs_col'
     708            CALL check_var(fname, varname, zh_qs_col, klon, largest, .FALSE.)
     709          END IF
     710        END IF
     711        STOP
     712      END IF
     713
     714      IF (qw_tot .NE. qw_tot .OR. ABS(qw_tot) > largest) THEN
     715        PRINT *,TRIM(errmsg)
     716        PRINT *,'  ' // TRIM(fname) // ': Wrong qw_tot= ',qw_tot,' !!!'
     717        PRINT *,'    qw_tot: total mass of watter vapour (kg/m2)'
     718        PRINT *,'    qw_tot = f(q,ql,qs,airephy,t,paprs)'
     719        varname = 'airephy'
     720        CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.)
     721        varname = 'q'
     722        CALL check_var3D(fname, varname, q, klon, klev, largest, .FALSE.)
     723        varname = 'ql'
     724        CALL check_var3D(fname, varname, ql, klon, klev, largest, .FALSE.)
     725        varname = 'qs'
     726        CALL check_var3D(fname, varname, qs, klon, klev, largest, .FALSE.)
     727        varname = 't'
     728        CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.)
     729        varname = 'paprs'
     730        CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.)
     731        STOP
     732      END IF
     733
     734      IF (ql_tot .NE. ql_tot .OR. ABS(ql_tot) > largest) THEN
     735        PRINT *,TRIM(errmsg)
     736        PRINT *,'  ' // TRIM(fname) // ': Wrong ql_tot= ',ql_tot,' !!!'
     737        PRINT *,'    ql_tot: total mass of liquid watter (kg/m2)'
     738        PRINT *,'    ql_tot = f(ql,airephy,t,paprs)'
     739        varname = 'airephy'
     740        CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.)
     741        varname = 'ql'
     742        CALL check_var3D(fname, varname, ql, klon, klev, largest, .FALSE.)
     743        varname = 't'
     744        CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.)
     745        varname = 'paprs'
     746        CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.)
     747        STOP
     748      END IF
     749
     750      IF (qs_tot .NE. qs_tot .OR. ABS(qs_tot) > largest) THEN
     751        PRINT *,TRIM(errmsg)
     752        PRINT *,'  ' // TRIM(fname) // ': Wrong qs_tot= ',qs_tot,' !!!'
     753        PRINT *,'    qs_tot: total mass of solid watter (kg/m2)'
     754        PRINT *,'    qs_tot = f(qs,airephy,t,paprs)'
     755        varname = 'airephy'
     756        CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.)
     757        varname = 'qs'
     758        CALL check_var3D(fname, varname, qs, klon, klev, largest, .FALSE.)
     759        varname = 't'
     760        CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.)
     761        varname = 'paprs'
     762        CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.)
     763        STOP
     764      END IF
     765
     766      IF (ec_tot .NE. ec_tot .OR. ABS(ec_tot) > largest*10.e4) THEN
     767        PRINT *,TRIM(errmsg)
     768        PRINT *,'  ' // TRIM(fname) // ': Wrong ec_tot= ',ec_tot,' !!!'
     769        PRINT *,'    ec_tot: total cinetic energy (kg/m2)'
     770        PRINT *,'    ec_tot = f(u,v,airephy,t,paprs)'
     771        varname = 'airephy'
     772        CALL check_var(fname, varname, airephy, klon, largest*10.e5, .FALSE.)
     773        varname = 'u'
     774        CALL check_var3D(fname, varname, u, klon, klev, largest, .FALSE.)
     775        varname = 'v'
     776        CALL check_var3D(fname, varname, v, klon, klev, largest, .FALSE.)
     777        varname = 't'
     778        CALL check_var3D(fname, varname, t, klon, klev, largest, .FALSE.)
     779        varname = 'paprs'
     780        CALL check_var3D(fname, varname, paprs, klon, klev, largest, .FALSE.)
     781        STOP
     782      END IF
     783
    451784      RETURN
     785
    452786      END SUBROUTINE diagetpq
    453787
  • trunk/WRFV3/lmdz/physiq.F90

    r6 r142  
    12991299      integer iostat
    13001300
    1301 ! Lluis
     1301! L. Fita, LMD. Point for checkings...
    13021302      INTEGER                                            :: llp
    13031303      llp = 644
     
    13071307!
    13081308      CALL phys_cal_update(jD_cur,jH_cur)
    1309       PRINT *,'  Lluis physiq: jD_cur: ',jD_cur, ' jH_cur: ',jH_cur, &
    1310         ' days_elapsed: ',days_elapsed
    13111309
    13121310!c======================================================================
     
    19621960       &      , fs_bound, fq_bound )
    19631961      END IF
     1962      PRINT *,'  Lluis 1 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    19641963
    19651964!c Diagnostiquer la tendance dynamique
     
    21222121!C
    21232122      END IF
    2124 
     2123      PRINT *,'  Lluis 2 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    21252124!c
    21262125!c=========================================================================
     
    21682167           zdtime=dtime*REAL(radpas) ! pas de temps du rayonnement (s)
    21692168           CALL zenang(zlongi,jH_cur,zdtime,rlat,rlon,rmu0,fract)
    2170            PRINT *,'  Lluis: radpas: ',radpas,' fract: ',fract(llp),  &
    2171              ' zdtime: ',zdtime,' jH_cur: ',jH_cur,' rmu0: ',rmu0(llp),&
    2172              ' dtime: ',dtime,' nbapp_rad: ',nbapp_rad
    2173            PRINT *,' Lluis @rlat@: ',rlat
    2174            PRINT *,' Lluis @rlon@: ',rlon
    2175            PRINT *,' Lluis @rmu0@: ',rmu0
    2176            PRINT *,' Lluis @fract@: ',fract
    21772169         ELSE
    21782170           CALL angle(zlongi, rlat, fract, rmu0)
    21792171         ENDIF
    21802172      ENDIF
    2181       PRINT *,'  Lluis zenang_an rlat: ',rlat(llp), ' rlon:',rlon(llp), &
    2182         ' rmu0: ',rmu0(llp),' fract: ',fract(llp),' jH_cur: ',jH_cur,  &
    2183         ' cycle_diurne: ',cycle_diurne
    21842173
    21852174      if (mydebug) then
     
    22752264
    22762265      ENDIF
     2266      PRINT *,'  Lluis 3 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    22772267!c =================================================================== c
    22782268!c   Calcul de Qsat
     
    26352625      END IF
    26362626!C
     2627      PRINT *,'  Lluis 4 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    26372628      IF (check) THEN
    26382629          za = qcheck(klon,klev,paprs,q_seri,ql_seri,airephy)
     
    27472738       &      , fs_bound, fq_bound )
    27482739      END IF
    2749 
     2740      PRINT *,'  Lluis 5 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    27502741!c      print*,'apres callwake iflag_cldcon=', iflag_cldcon
    27512742!c
     
    27712762!c      detr_therm(:,:)=0.
    27722763!c
     2764
    27732765      IF(prt_level>9)WRITE(lunout,*)                                                 &
    27742766       &    'AVANT LA CONVECTION SECHE , iflag_thermals='                            &
     
    29582950
    29592951      endif
     2952
    29602953!c
    29612954!c===================================================================
     
    29722965       &      , fs_bound, fq_bound )
    29732966      END IF
    2974 
     2967      PRINT *,'  Lluis 6 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    29752968
    29762969!c-------------------------------------------------------------------------
     
    30393032       &      , fs_bound, fq_bound )
    30403033      END IF
     3034      PRINT *,'  Lluis 7 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    30413035      if (mydebug) then
    30423036        call writefield_phy('u_seri',u_seri,llm)
     
    30453039        call writefield_phy('q_seri',q_seri,llm)
    30463040      endif
    3047 
    30483041!c
    30493042!c-------------------------------------------------------------------
     
    32713264       &      , fs_bound, fq_bound )
    32723265      END IF
     3266      PRINT *,'  Lluis 8 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    32733267!c
    32743268!c Calculer l'humidite relative pour diagnostique
     
    36433637       &      , fs_bound, fq_bound )
    36443638      END IF
     3639      PRINT *,'  Lluis 9 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    36453640!c
    36463641!c
     
    38173812       &      , fs_bound, fq_bound )
    38183813      END IF
     3814      PRINT *,'  Lluis 10 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    38193815!c
    38203816!c
     
    39533949!C     Donc la somme de ces 2 variations devrait etre nulle.
    39543950
     3951       PRINT *,'  LLuis sending: d_h_vcol: ',d_h_vcol,' d_qt ',d_qt,' d_ec ',d_ec
     3952
    39553953        call diagphy(airephy,ztit,ip_ebil_phy                                        &
    39563954       &      , topsw, toplw, solsw, sollw, sens                                     &
     
    39623960!C
    39633961      END IF
     3962      PRINT *,'  Lluis 11 d_h_vcol: ',d_h_vcol,' d_h_vcol_phy: ',d_h_vcol_phy
    39643963      PRINT *,'Lluis Reaching the SORTIES point'
    39653964
  • trunk/WRFV3/lmdz/thermcell_dq.F90

    r1 r142  
    3737      CHARACTER (LEN=20) :: modname='thermcell_dq'
    3838      CHARACTER (LEN=80) :: abort_message
    39 
    4039
    4140! Old explicite scheme
     
    115114      fqa(:,1)=0. ; fqa(:,nlay)=0.
    116115
    117 
    118116! Trace species evolution
    119117   if (impl==0) then
  • trunk/WRFV3/lmdz/thermcell_flux2.F90

    r1 r142  
    150150      enddo
    151151
    152 
    153 
    154152! Test provisoire : pour comprendre pourquoi on corrige plein de fois
    155153! le cas fm6, on commence par regarder une premiere fois avant les
     
    192190            endif
    193191         enddo
    194 
    195192
    196193!-------------------------------------------------------------------------
     
    461458         enddo
    462459      enddo
     460
    463461      if (labort_gcm) then
    464462                         ig=igout
  • trunk/WRFV3/lmdz/thermcell_height.F90

    r1 r142  
    6464      enddo
    6565
     66! L. Fita, LMD July 2014. Fixing the If when zw2 == 0
    6667      do l=1,nlay
    6768         do ig=1,ngrid
    6869            if (l.le.lmax(ig)) then
    69                 if (zw2(ig,l).lt.0.)then
    70                   print*,'pb2 zw2<0'
    71                 endif
     70! Previous version
     71!                if (zw2(ig,l).lt.0.)then
     72!                  print*,'pb2 zw2<0'
     73!                endif
     74!                zw2(ig,l)=sqrt(zw2(ig,l))
     75!                wmax(ig)=max(wmax(ig),zw2(ig,l))
     76!            else
     77!                 zw2(ig,l)=0.
     78!            endif
     79! New version
     80              if (zw2(ig,l).lt.0.)then
     81                print*,'pb2 zw2<0'
     82                zw2(ig,l)=0.
     83                wmax(ig)=max(wmax(ig),zw2(ig,l))
     84              else
    7285                zw2(ig,l)=sqrt(zw2(ig,l))
    7386                wmax(ig)=max(wmax(ig),zw2(ig,l))
     87              endif 
    7488            else
    7589                 zw2(ig,l)=0.
    7690            endif
     91
    7792          enddo
    7893      enddo
  • trunk/WRFV3/lmdz/thermcell_main.F90

    r1 r142  
    105105
    106106      INTEGER ig,k,l,ll,ierr
    107       real zsortie1d(klon)
    108       INTEGER lmax(klon),lmin(klon),lalim(klon)
    109       INTEGER lmix(klon)
    110       INTEGER lmix_bis(klon)
    111       real linter(klon)
    112       real zmix(klon)
    113       real zmax(klon),zw2(klon,klev+1),ztva(klon,klev),zw_est(klon,klev+1),ztva_est(klon,klev)
    114 !      real fraca(klon,klev)
    115 
    116       real zmax_sec(klon)
     107      real zsortie1d(ngrid)
     108      INTEGER lmax(ngrid),lmin(ngrid),lalim(ngrid)
     109      INTEGER lmix(ngrid)
     110      INTEGER lmix_bis(ngrid)
     111      real linter(ngrid)
     112      real zmix(ngrid)
     113      real zmax(ngrid),zw2(ngrid,nlay+1),ztva(ngrid,nlay),zw_est(ngrid,nlay+1),ztva_est(ngrid,nlay)
     114!      real fraca(ngrid,nlay)
     115
     116      real zmax_sec(ngrid)
    117117!on garde le zmax du pas de temps precedent
    118       real zmax0(klon)
     118      real zmax0(ngrid)
    119119!FH/IM     save zmax0
    120120
    121121      real lambda
    122122
    123       real zlev(klon,klev+1),zlay(klon,klev)
    124       real deltaz(klon,klev)
    125       REAL zh(klon,klev)
    126       real zthl(klon,klev),zdthladj(klon,klev)
    127       REAL ztv(klon,klev)
    128       real zu(klon,klev),zv(klon,klev),zo(klon,klev)
    129       real zl(klon,klev)
    130       real zsortie(klon,klev)
    131       real zva(klon,klev)
    132       real zua(klon,klev)
    133       real zoa(klon,klev)
    134 
    135       real zta(klon,klev)
    136       real zha(klon,klev)
    137       real fraca(klon,klev+1)
     123      real zlev(ngrid,nlay+1),zlay(ngrid,nlay)
     124      real deltaz(ngrid,nlay)
     125      REAL zh(ngrid,nlay)
     126      real zthl(ngrid,nlay),zdthladj(ngrid,nlay)
     127      REAL ztv(ngrid,nlay)
     128      real zu(ngrid,nlay),zv(ngrid,nlay),zo(ngrid,nlay)
     129      real zl(ngrid,nlay)
     130      real zsortie(ngrid,nlay)
     131      real zva(ngrid,nlay)
     132      real zua(ngrid,nlay)
     133      real zoa(ngrid,nlay)
     134
     135      real zta(ngrid,nlay)
     136      real zha(ngrid,nlay)
     137      real fraca(ngrid,nlay+1)
    138138      real zf,zf2
    139       real thetath2(klon,klev),wth2(klon,klev),wth3(klon,klev)
    140       real q2(klon,klev)
     139      real thetath2(ngrid,nlay),wth2(ngrid,nlay),wth3(ngrid,nlay)
     140      real q2(ngrid,nlay)
    141141! FH probleme de dimensionnement avec l'allocation dynamique
    142142!     common/comtherm/thetath2,wth2
    143       real wq(klon,klev)
    144       real wthl(klon,klev)
    145       real wthv(klon,klev)
     143      real wq(ngrid,nlay)
     144      real wthl(ngrid,nlay)
     145      real wthv(ngrid,nlay)
    146146   
    147       real ratqscth(klon,klev)
     147      real ratqscth(ngrid,nlay)
    148148      real var
    149149      real vardiff
    150       real ratqsdiff(klon,klev)
     150      real ratqsdiff(ngrid,nlay)
    151151
    152152      logical sorties
    153       real rho(klon,klev),rhobarz(klon,klev),masse(klon,klev)
    154       real zpspsk(klon,klev)
    155 
    156       real wmax(klon)
    157       real wmax_tmp(klon)
    158       real wmax_sec(klon)
    159       real fm0(klon,klev+1),entr0(klon,klev),detr0(klon,klev)
    160       real fm(klon,klev+1),entr(klon,klev),detr(klon,klev)
    161 
    162       real ztla(klon,klev),zqla(klon,klev),zqta(klon,klev)
     153      real rho(ngrid,nlay),rhobarz(ngrid,nlay),masse(ngrid,nlay)
     154      real zpspsk(ngrid,nlay)
     155
     156      real wmax(ngrid)
     157      real wmax_tmp(ngrid)
     158      real wmax_sec(ngrid)
     159      real fm0(ngrid,nlay+1),entr0(ngrid,nlay),detr0(ngrid,nlay)
     160      real fm(ngrid,nlay+1),entr(ngrid,nlay),detr(ngrid,nlay)
     161
     162      real ztla(ngrid,nlay),zqla(ngrid,nlay),zqta(ngrid,nlay)
    163163!niveau de condensation
    164       integer nivcon(klon)
    165       real zcon(klon)
     164      integer nivcon(ngrid)
     165      real zcon(ngrid)
    166166      REAL CHI
    167       real zcon2(klon)
    168       real pcon(klon)
    169       real zqsat(klon,klev)
    170       real zqsatth(klon,klev)
    171 
    172       real f_star(klon,klev+1),entr_star(klon,klev)
    173       real detr_star(klon,klev)
    174       real alim_star_tot(klon)
    175       real alim_star(klon,klev)
    176       real alim_star_clos(klon,klev)
    177       real f(klon), f0(klon)
     167      real zcon2(ngrid)
     168      real pcon(ngrid)
     169      real zqsat(ngrid,nlay)
     170      real zqsatth(ngrid,nlay)
     171
     172      real f_star(ngrid,nlay+1),entr_star(ngrid,nlay)
     173      real detr_star(ngrid,nlay)
     174      real alim_star_tot(ngrid)
     175      real alim_star(ngrid,nlay)
     176      real alim_star_clos(ngrid,nlay)
     177      real f(ngrid), f0(ngrid)
    178178!FH/IM     save f0
    179       real zlevinter(klon)
     179      real zlevinter(ngrid)
    180180      logical debut
    181181       real seuil
    182       real csc(klon,klev)
     182      real csc(ngrid,nlay)
    183183
    184184!!! nrlmd le 10/04/2012
    185185
    186186!------Entrées
    187       real pbl_tke(klon,klev+1,nbsrf)
    188       real pctsrf(klon,nbsrf)
    189       real omega(klon,klev)
    190       real airephy(klon)
     187      real pbl_tke(ngrid,nlay+1,nbsrf)
     188      real pctsrf(ngrid,nbsrf)
     189      real omega(ngrid,nlay)
     190      real airephy(ngrid)
    191191!------Sorties
    192       real zlcl(klon),fraca0(klon),w0(klon),w_conv(klon)
    193       real therm_tke_max0(klon),env_tke_max0(klon)
    194       real n2(klon),s2(klon)
    195       real ale_bl_stat(klon)
    196       real therm_tke_max(klon,klev),env_tke_max(klon,klev)
    197       real alp_bl_det(klon),alp_bl_fluct_m(klon),alp_bl_fluct_tke(klon),alp_bl_conv(klon),alp_bl_stat(klon)
     192      real zlcl(ngrid),fraca0(ngrid),w0(ngrid),w_conv(ngrid)
     193      real therm_tke_max0(ngrid),env_tke_max0(ngrid)
     194      real n2(ngrid),s2(ngrid)
     195      real ale_bl_stat(ngrid)
     196      real therm_tke_max(ngrid,nlay),env_tke_max(ngrid,nlay)
     197      real alp_bl_det(ngrid),alp_bl_fluct_m(ngrid),alp_bl_fluct_tke(ngrid),alp_bl_conv(ngrid),alp_bl_stat(ngrid)
    198198!------Local
    199199      integer nsrf
    200       real rhobarz0(klon)                    ! Densité au LCL
    201       logical ok_lcl(klon)                   ! Existence du LCL des thermiques
    202       integer klcl(klon)                     ! Niveau du LCL
    203       real interp(klon)                      ! Coef d'interpolation pour le LCL
     200      real rhobarz0(ngrid)                    ! Densité au LCL
     201      logical ok_lcl(ngrid)                   ! Existence du LCL des thermiques
     202      integer klcl(ngrid)                     ! Niveau du LCL
     203      real interp(ngrid)                      ! Coef d'interpolation pour le LCL
    204204!--Triggering
    205205      real Su                                ! Surface unité: celle d'un updraft élémentaire
     
    215215      real zmax_moy_coef
    216216      parameter(zmax_moy_coef=0.33)
    217       real depth(klon)                       ! Epaisseur moyenne du cumulus
    218       real w_max(klon)                       ! Vitesse max statistique
    219       real s_max(klon)
     217      real depth(ngrid)                       ! Epaisseur moyenne du cumulus
     218      real w_max(ngrid)                       ! Vitesse max statistique
     219      real s_max(ngrid)
    220220!--Closure
    221       real pbl_tke_max(klon,klev)            ! Profil de TKE moyenne
    222       real pbl_tke_max0(klon)                ! TKE moyenne au LCL
    223       real w_ls(klon,klev)                   ! Vitesse verticale grande échelle (m/s)
     221      real pbl_tke_max(ngrid,nlay)            ! Profil de TKE moyenne
     222      real pbl_tke_max0(ngrid)                ! TKE moyenne au LCL
     223      real w_ls(ngrid,nlay)                   ! Vitesse verticale grande échelle (m/s)
    224224      real coef_m                            ! On considère un rendement pour alp_bl_fluct_m
    225225      parameter(coef_m=1.)
     
    231231!
    232232      !nouvelles variables pour la convection
    233       real Ale_bl(klon)
    234       real Alp_bl(klon)
    235       real alp_int(klon),dp_int(klon),zdp
    236       real ale_int(klon)
    237       integer n_int(klon)
    238       real fm_tot(klon)
    239       real wght_th(klon,klev)
    240       integer lalim_conv(klon)
     233      real Ale_bl(ngrid)
     234      real Alp_bl(ngrid)
     235      real alp_int(ngrid),dp_int(ngrid),zdp
     236      real ale_int(ngrid)
     237      integer n_int(ngrid)
     238      real fm_tot(ngrid)
     239      real wght_th(ngrid,nlay)
     240      integer lalim_conv(ngrid)
    241241!v1d     logical therm
    242242!v1d     save therm
     
    250250      EXTERNAL SCOPY
    251251!
     252
     253! L. Fita, LMD July 2014. Initializing variables.
     254!   Some not initializated according to values: iflag_trig_bl, iflag_clos_bl
     255     zdthladj = 0.
     256     pbl_tke_max0 = 0.
     257     fraca0 = 0.
     258     w_conv = 0.
     259     w0 = 0.
     260     therm_tke_max0 = 0.
     261     env_tke_max0 = 0.
     262     alp_bl_det = 0.
     263     alp_bl_fluct_m = 0.
     264     alp_bl_fluct_tke = 0.
     265     alp_bl_conv = 0.
     266     alp_bl_stat = 0.
     267     interp = 0.
     268     klcl = 0
    252269
    253270!-----------------------------------------------------------------------
     
    306323
    307324       sorties=.true.
    308       IF(ngrid.NE.klon) THEN
     325      IF(ngrid.NE.ngrid) THEN
    309326         PRINT*
    310327         PRINT*,'STOP dans convadj'
    311328         PRINT*,'ngrid    =',ngrid
    312          PRINT*,'klon  =',klon
     329         PRINT*,'ngrid  =',ngrid
    313330      ENDIF
    314331!
    315332!     write(lunout,*)'WARNING thermcell_main f0=max(f0,1.e-2)'
    316      do ig=1,klon
     333     do ig=1,ngrid
    317334         f0(ig)=max(f0(ig),1.e-2)
    318335         zmax0(ig)=max(zmax0(ig),40.)
     
    331348      CALL thermcell_env(ngrid,nlay,po,pt,pu,pv,pplay,  &
    332349     &           pplev,zo,zh,zl,ztv,zthl,zu,zv,zpspsk,zqsat,lev_out)
    333        
     350
    334351      if (prt_level.ge.1) print*,'thermcell_main apres thermcell_env'
    335352
     
    360377      enddo
    361378         zlev(:,1)=0.
    362          zlev(:,nlay+1)=(2.*pphi(:,klev)-pphi(:,klev-1))/RG
     379         zlev(:,nlay+1)=(2.*pphi(:,nlay)-pphi(:,nlay-1))/RG
    363380      do l=1,nlay
    364381         zlay(:,l)=pphi(:,l)/RG
     
    516533!     print*,"ZMAX ",lalim,lmin,linter,lmix,lmax,zmax,zmax0,zmix,wmax
    517534
    518 
    519 
    520535      call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lalim ')
    521536      call test_ltherm(ngrid,nlay,pplev,pplay,lmin ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmin  ')
     
    547562      endif
    548563
    549 
    550 
    551 
    552564! Choix de la fonction d'alimentation utilisee pour la fermeture.
    553565! Apparemment sans importance
     
    584596!deduction des flux
    585597!-------------------------------------------------------------------------------
    586 
    587598      CALL thermcell_flux2(ngrid,nlay,ptimestep,masse, &
    588599     &       lalim,lmax,alim_star,  &
     
    620631      call thermcell_dq(ngrid,nlay,dqimpl,ptimestep,fm0,entr0,masse,  &
    621632     &                   po,pdoadj,zoa,lev_out)
    622 
    623633!------------------------------------------------------------------
    624634! Calcul de la fraction de l'ascendance
    625635!------------------------------------------------------------------
    626       do ig=1,klon
     636      do ig=1,ngrid
    627637         fraca(ig,1)=0.
    628638         fraca(ig,nlay+1)=0.
    629639      enddo
    630640      do l=2,nlay
    631          do ig=1,klon
     641         do ig=1,ngrid
    632642            if (zw2(ig,l).gt.1.e-10) then
    633643            fraca(ig,l)=fm(ig,l)/(rhobarz(ig,l)*zw2(ig,l))
     
    753763         enddo
    754764      enddo
     765
    755766!calcul des flux: q, thetal et thetav
    756767      do l=1,nlay
     
    768779    do ig=1,ngrid
    769780      ok_lcl(ig)=.false.
    770       if ( (pcon(ig) .gt. pplay(ig,klev-1)) .and. (pcon(ig) .lt. pplay(ig,1)) ) ok_lcl(ig)=.true.
     781      if ( (pcon(ig) .gt. pplay(ig,nlay-1)) .and. (pcon(ig) .lt. pplay(ig,1)) ) ok_lcl(ig)=.true.
    771782    enddo
    772783
     
    775786      do ig=1,ngrid
    776787        if (ok_lcl(ig)) then
    777           if ((pplay(ig,l) .ge. pcon(ig)) .and. (pplay(ig,l+1) .le. pcon(ig))) then
    778           klcl(ig)=l
    779           interp(ig)=(pcon(ig)-pplay(ig,klcl(ig)))/(pplay(ig,klcl(ig)+1)-pplay(ig,klcl(ig)))
     788!          if ((pplay(ig,l) .ge. pcon(ig)) .and. (pplay(ig,l+1) .le. pcon(ig))) then
     789! L. Fita, LMD July 2014. Adding avoiding divisons by zero...
     790          if ((pplay(ig,l).ge.pcon(ig)) .and. (pplay(ig,l+1).le.pcon(ig)) .and.      &
     791            (klcl(ig).gt.0).and.(klcl(ig)+1.le.nlay-1) .and. (klcl(ig)+1.gt.0) ) then
     792            klcl(ig)=l
     793            interp(ig)=(pcon(ig)-pplay(ig,klcl(ig)))/(pplay(ig,klcl(ig)+1)-pplay(ig,klcl(ig)))
    780794          endif
    781795        endif
     
    794808    do ig =1,ngrid
    795809     zmax(ig)=pphi(ig,lmax(ig))/rg
    796      if (ok_lcl(ig)) then
    797       rhobarz0(ig)=rhobarz(ig,klcl(ig))+(rhobarz(ig,klcl(ig)+1) &
     810!     if (ok_lcl(ig)) then
     811! L. Fita, LMD July 2014. Adding avoiding divisons by zero...
     812     if ( ok_lcl(ig) .and. (klcl(ig).gt.0) .and. (klcl(ig)+1.le.nlay-1) .and.        &
     813       (klcl(ig)+1.gt.0) ) then
     814       rhobarz0(ig)=rhobarz(ig,klcl(ig))+(rhobarz(ig,klcl(ig)+1) &
    798815 &               -rhobarz(ig,klcl(ig)))*interp(ig)
    799       zlcl(ig)=(pplev(ig,1)-pcon(ig))/(rhobarz0(ig)*RG)
    800       zlcl(ig)=min(zlcl(ig),zmax(ig))   ! Si zlcl > zmax alors on pose zlcl = zmax
     816       zlcl(ig)=(pplev(ig,1)-pcon(ig))/(rhobarz0(ig)*RG)
     817       zlcl(ig)=min(zlcl(ig),zmax(ig))   ! Si zlcl > zmax alors on pose zlcl = zmax
    801818     else
    802819      rhobarz0(ig)=0.
     
    868885     fraca0(ig)=0.
    869886     w0(ig)=0.
     887! L. Fita, LMD July 2014. Adding zero value for stability issues
     888     w_conv(ig) = 0.
    870889!!jyg le 27/04/2012
    871890!!     zlcl(ig)=0.
     
    933952  ENDIF ! iflag_trig_bl
    934953!  print *,'ENDIF  iflag_trig_bl'    !!jyg
    935 
    936954!------------Closure------------------
    937955
     
    9911009      lalim_conv(:)=lalim(:)
    9921010
    993       do k=1,klev
     1011      do k=1,nlay
    9941012         do ig=1,ngrid
    9951013            if (k<=lalim_conv(ig)) fm_tot(ig)=fm_tot(ig)+fm(ig,k)
     
    9991017! assez bizarre car, si on est dans la couche d'alim et que alim_star et
    10001018! plus petit que 1.e-10, on prend wght_th=1.
    1001       do k=1,klev
     1019      do k=1,nlay
    10021020         do ig=1,ngrid
    10031021            if (k<=lalim_conv(ig).and.alim_star(ig,k)>1.e-10) then
     
    10591077      ratqsdiff(:,:)=0.
    10601078
    1061       do l=1,klev
     1079      do l=1,nlay
    10621080         do ig=1,ngrid
    10631081            if (l<=lalim(ig)) then
     
    10691087      if (prt_level.ge.1) print*,'14f OK convect8'
    10701088
    1071       do l=1,klev
     1089      do l=1,nlay
    10721090         do ig=1,ngrid
    10731091            if (l<=lalim(ig)) then
     
    11051123!-----------------------------------------------------------------------------
    11061124
    1107       subroutine test_ltherm(klon,klev,pplev,pplay,long,seuil,ztv,po,ztva,zqla,f_star,zw2,comment)
     1125      subroutine test_ltherm(ngrid,nlay,pplev,pplay,long,seuil,ztv,po,ztva,zqla,f_star,zw2,comment)
    11081126      IMPLICIT NONE
    11091127#include "iniprint.h"
    11101128
    1111       integer i, k, klon,klev
    1112       real pplev(klon,klev+1),pplay(klon,klev)
    1113       real ztv(klon,klev)
    1114       real po(klon,klev)
    1115       real ztva(klon,klev)
    1116       real zqla(klon,klev)
    1117       real f_star(klon,klev)
    1118       real zw2(klon,klev)
    1119       integer long(klon)
     1129      integer i, k, ngrid,nlay
     1130      real pplev(ngrid,nlay+1),pplay(ngrid,nlay)
     1131      real ztv(ngrid,nlay)
     1132      real po(ngrid,nlay)
     1133      real ztva(ngrid,nlay)
     1134      real zqla(ngrid,nlay)
     1135      real f_star(ngrid,nlay)
     1136      real zw2(ngrid,nlay)
     1137      integer long(ngrid)
    11201138      real seuil
    11211139      character*21 comment
     
    11271145
    11281146!  test sur la hauteur des thermiques ...
    1129          do i=1,klon
     1147         do i=1,ngrid
    11301148!IMtemp           if (pplay(i,long(i)).lt.seuil*pplev(i,1)) then
    11311149           if (prt_level.ge.10) then
    11321150               print*,'WARNING ',comment,' au point ',i,' K= ',long(i)
    11331151               print*,'  K  P(MB)  THV(K)     Qenv(g/kg)THVA        QLA(g/kg)   F*        W2'
    1134                do k=1,klev
     1152               do k=1,nlay
    11351153                  write(6,'(i3,7f10.3)') k,pplay(i,k),ztv(i,k),1000*po(i,k),ztva(i,k),1000*zqla(i,k),f_star(i,k),zw2(i,k)
    11361154               enddo
  • trunk/WRFV3/lmdz/thermcell_plume.F90

    r1 r142  
    9090      REAL zw2fact,expa
    9191      Zsat=.false.
     92
    9293! Initialisation
    9394      RLvCp = RLVTT/RCPD
     
    167168      enddo
    168169      alim_star_tot(:)=1.
    169 
    170170
    171171!------------------------------------------------------------------------------
     
    192192enddo
    193193!
    194 
    195194!==============================================================================
    196195!boucle de calcul de la vitesse verticale dans le thermique
     
    539538      REAL fact_gamma,fact_epsilon,fact_gamma2,fact_epsilon2
    540539      REAL c2(ngrid,klev)
     540
    541541      Zsat=.false.
    542542! Initialisation
     
    603603         do ig=1,ngrid
    604604            if (ztv(ig,l)> ztv(ig,l+1) .and. ztv(ig,1)>=ztv(ig,l) ) then
    605                alim_star(ig,l)=MAX((ztv(ig,l)-ztv(ig,l+1)),0.)  &
     605! L. Fita, LMD July 2014. Adding IF for zlev == 0.
     606               IF ( zlev(ig,l+1) < 0.) THEN
     607                 alim_star(ig,l) = 0.
     608               ELSE
     609                 alim_star(ig,l)=MAX((ztv(ig,l)-ztv(ig,l+1)),0.)  &
    606610     &                       *sqrt(zlev(ig,l+1))
    607                lalim(ig)=l+1
     611                 lalim(ig)=l+1
     612               END IF
    608613               alim_star_tot(ig)=alim_star_tot(ig)+alim_star(ig,l)
    609614            endif
     
    618623      enddo
    619624      alim_star_tot(:)=1.
    620 
    621 
    622 
    623625!------------------------------------------------------------------------------
    624626! Calcul dans la premiere couche
     
    686688
    687689!------------------------------------------------
    688 !AJAM:nouveau calcul de w² 
     690!AJAM:nouveau calcul de w\B2 
    689691!------------------------------------------------
    690692              zdz=zlev(ig,l+1)-zlev(ig,l)
     
    741743      endif
    742744   enddo
    743 
    744745
    745746!----------------------------------------------------------------------------
Note: See TracChangeset for help on using the changeset viewer.