* PACKAGE P2CCLL    !" 物理過程(2) 対流雲層
*
***********************************************************************
      SUBROUTINE CCLLEV     !" 対流雲層
     O         ( GDTC  , GDQC  , KBTM  , KTOP  ,
     I           GDT   , GDQ   , GDP   , GDZ   ,
     I           IREUSE,
     W           ADIAF                           )
*
*"       雲底下端 = LCL
*"       エントレインメント
*
*   [PARAM] 
#ifdef SYS_IBMS
      INCLUDE   (ZCDIM)                      !" 格子点数, 波数
      INCLUDE   (ZCCOM)                      !" 標準物理定数
#else
#include        "zcdim.F"                    !" 格子点数, 波数
#include        "zccom.F"                    !" 標準物理定数
#endif
*
*   [OUTPUT] 
      REAL       GDTC  ( IDIM*JDIM, KMAX )   !" 雲の温度
      REAL       GDQC  ( IDIM*JDIM, KMAX )   !" 雲の比湿
      INTEGER    KBTM  ( IDIM*JDIM )         !" 雲底レベル
      INTEGER    KTOP  ( IDIM*JDIM )         !" 雲頂レベル
*
*   [INPUT] 
      REAL       GDT   ( IDIM*JDIM, KMAX )   !" 温度  Ｔ
      REAL       GDQ   ( IDIM*JDIM, KMAX )   !" 比湿  ｑ
      REAL       GDP   ( IDIM*JDIM, KMAX )   !" 気圧  Ｐ
      REAL       GDZ   ( IDIM*JDIM, KMAX )   !" 高度  Ｚ 
*
      INTEGER    IREUSE                      !" 前値再利用フラグ
*
*   [WORK] 
      REAL       ADIAF ( IDIM*JDIM, KMAX )   !" 乾燥断熱降温率
*
*   [INTERNAL PARM] 
      INTEGER    KBMIN
      INTEGER    KBMAX
      DATA       KBMIN  / 1    /
      DATA       KBMAX  / KMAX /
      REAL       EPST                        !" 不安定の許容誤差
      DATA       EPST  / 0.0   /
      REAL       ENTR                        !" エントレイン率
      DATA       ENTR  / 0.0   /
      NAMELIST  /NMCCLL/  KBMIN, KBMAX, EPST, ENTR
*
      LOGICAL    OFIRST, OSET
      DATA       OFIRST, OSET  / .TRUE., .FALSE. /
      SAVE
*
*   [INTERNAL WORK] 
      INTEGER    KBTMZ ( IDIM*JDIM )         !" 雲底レベル
      INTEGER    KTOPZ ( IDIM*JDIM )         !" 雲頂レベル
*     
      INTEGER    IJ, K, KUNDEF
      INTEGER    IFPAR, JFPAR
      REAL       AKAPPA, ELF
      REAL       RENTR, QSAT, DQSAT, QEXE, DT, DQ
*
#ifdef SYS_IBMS
      INCLUDE   (ZQSAT)                      !" 飽和水蒸気圧
#else
#include        "zqsat.F"                    !" 飽和水蒸気圧
#endif
*
      IF ( OFIRST ) THEN
         WRITE ( 6,* ) ' CONVECTIVE CLOUD LAYER P2.3 DATE=92/02/07'
         OFIRST = .FALSE.
*
         CALL   REWNML ( IFPAR , JFPAR )
         READ   ( IFPAR, NMCCLL, END=90 )
   90    WRITE  ( JFPAR, NMCCLL )
*
         KUNDEF = -1
      ENDIF
*
*"         < 0. IREUSE = 1: 前値再利用 (shallow convection) >
*
      IF ( ( IREUSE .EQ. 1 ) .AND. OSET ) THEN
         DO 100 IJ = 1, IDIM*JDIM
            KBTM ( IJ ) = KBTMZ ( IJ )
            KTOP ( IJ ) = KTOPZ ( IJ )
  100    CONTINUE
         OSET = .FALSE.
      ENDIF
*
      IF ( IREUSE .NE. 1 ) OSET = .TRUE.
*
*
*"         < 1. ファクター計算 >
*
      AKAPPA = RAIR / CP
      DO 1100 K = 1, KMAX-1
         DO 1100 IJ = 1, IDIM*JDIM
            ADIAF( IJ,K ) = ( GDP( IJ,K+1 )/GDP( IJ,K ) ) **AKAPPA
 1100 CONTINUE
*
      ELF    = EL / CP
*
*"         < 2. 初期設定 >
*
      DO 2100 IJ = 1, IDIM*JDIM
         KBTM  ( IJ ) = KUNDEF
         KTOP  ( IJ ) = KUNDEF
 2100 CONTINUE
*
      DO 2200 IJ = 1, IDIM*JDIM
         GDTC( IJ,1 ) = GDT( IJ,1 )
         GDQC( IJ,1 ) = GDQ( IJ,1 )
 2200 CONTINUE
*
*"         < 3. 雲底・雲頂 >
*
      DO 3100 K = 2, KMAX
         DO 3110 IJ = 1, IDIM*JDIM
*
            IF ( KBTM( IJ ) .EQ. KUNDEF ) THEN
*
*"               < 3.1 持ち上げ凝結高度 >
*
               GDTC( IJ,K ) = GDTC( IJ,K-1 ) * ADIAF( IJ,K-1 )
               GDQC( IJ,K ) = GDQC( IJ,K-1 )
               QSAT         = FQSAT( GDTC( IJ,K ),GDP( IJ,K ) )
*
               IF (  ( GDQC( IJ,K ) .GE. QSAT ).AND.
     &               ( K .GE. KBMIN )          .AND.
     &               ( K .LE. KBMAX )                ) THEN
                  KBTM( IJ )   = K
                  GDTC( IJ,K ) = GDT( IJ,K )
                  GDQC( IJ,K ) = FQSAT( GDTC( IJ,K ), GDP( IJ,K ) )
               ENDIF
*
            ELSE IF ( KTOP( IJ ) .EQ. KUNDEF )  THEN
*
*"               < 3.2 湿潤上昇 >
*
               RENTR        = ENTR * ( GDZ( IJ,K ) - GDZ( IJ,K-1 ) )
               RENTR        = MIN( RENTR, 1. )
*
               GDTC( IJ,K ) = ( 1.-RENTR ) 
     &                         * GDTC( IJ,K-1 ) * ADIAF( IJ,K-1 )
     &                       + RENTR 
     &                         * ( GDT( IJ,K-1 ) + GDT( IJ,K ) )/2.  
               GDQC( IJ,K ) = ( 1.  -RENTR ) * GDQC( IJ,K-1 )
     &                       + RENTR 
     &                         * ( GDQ( IJ,K-1 ) + GDQ( IJ,K ) )/2.  
*
               QSAT         = FQSAT( GDTC( IJ,K ),GDP( IJ,K ) )
*
               IF ( GDQC( IJ,K ) .GT. QSAT ) THEN
                  DQSAT        = FDQSAT( GDTC( IJ,K ), QSAT  )
                  QEXE         = GDQC ( IJ,K ) - QSAT
                  DT           = ELF * QEXE / ( 1. + ELF * DQSAT )
                  DQ           = DQSAT * DT
*
                  GDTC( IJ,K ) = GDTC( IJ,K ) + DT
                  GDQC( IJ,K ) = QSAT         + DQ
               ENDIF
*
*"               < 3.3 安定度チェック >
*
               IF ( GDTC( IJ,K ) .LT. GDT( IJ,K )-EPST ) THEN
*
                  IF ( K .EQ. KBTM( IJ )+1  ) THEN
                     KBTM( IJ   ) = K
                     GDTC( IJ,K ) = GDT( IJ,K )
                     GDQC( IJ,K ) = FQSAT( GDTC( IJ,K ), GDP( IJ,K ) )
                  ELSE
                     KTOP( IJ ) = K - 1
                  ENDIF
*
               ELSE IF ( K .EQ. KMAX ) THEN
                  KTOP( IJ ) = K
               ENDIF
*
            ENDIF
*
 3110    CONTINUE
 3100 CONTINUE
*
      DO 3300 IJ = 1, IDIM*JDIM
         IF ( KTOP( IJ ) .EQ. KUNDEF ) THEN
            KBTM( IJ ) = KUNDEF
         ENDIF
*
         KBTMZ ( IJ ) = KBTM ( IJ )
         KTOPZ ( IJ ) = KTOP ( IJ )
 3300 CONTINUE
*
      RETURN
      END
