** PACKAGE P2ADBT  !"  物理過程(2) 平衡モデル用 湿潤断熱線
*
*"  [HIS]  94/05/02(ishiwatari)
**********************************************************************
      SUBROUTINE ADIABT       !" 湿潤断熱構造
     M         ( GDT   , GDQ   ,
     I           GDP   , GDPM  , GDZ   ,
     I           GDTG  , GDQG  )
*
*   [PARAM] 
#ifdef SYS_IBMS
      INCLUDE   (ZCDIM)                      !" 格子点数，波数
      INCLUDE   (ZCCOM)                      !" 標準物理定数
#else
#include        "zcdim.F"                    !" 格子点数，波数
#include        "zccom.F"                    !" 標準物理定数
#endif
*
*   [MODIFY] 
      REAL       GDT   ( IDIM*JDIM, KMAX )   !" 温度  Ｔ
      REAL       GDQ   ( IDIM*JDIM, KMAX )   !" 比湿  ｑ
      REAL       GQPR  ( IDIM*JDIM       )   !" 降水量
      REAL       GQCNV ( IDIM*JDIM, KMAX )   !" 凝結調節加熱率
      REAL       GSCNV ( IDIM*JDIM, KMAX )   !" 凝結調節比湿変化
*
*   [INPUT] 
      REAL       GDP    ( IDIM*JDIM, KMAX )  !" 気圧  Ｐ
      REAL       GDPM   ( IDIM*JDIM, KMAX+1 )!" 気圧（半整数）
      REAL       GDZ    ( IDIM*JDIM, KMAX )  !" 高度
      REAL       GDTG   ( IDIM*JDIM )  !" 地表面温度
      REAL       GDQG   ( IDIM*JDIM )  !" 地表面 q 
C      REAL       DELTA                       !" 2Δt
*
*   [INTERNAL PARM] 
      INTEGER    ITRMAX
      PARAMETER (ITRMAX=10)
      REAL       EPST ( ITRMAX )             !" 不安定の許容誤差
      REAL       RCRIT                       !" 臨界相対湿度
c      DATA       EPST / 0.01  , 0.02  , 0.02  , 0.05  , 0.05  ,
c     &                  0.10  , 0.10  , 0.20  , 0.20  , 0.40    /
c      DATA       RCRIT/ 0.990 /
*
      NAMELIST  /NMMADJ/ EPST, RCRIT
*
      LOGICAL    OFIRST
      DATA       OFIRST / .TRUE. /
      SAVE
*
*   [INTERNAL WORK] 
      COMMON    /COMWRK/
     &           ADJDA , ADJDB ,
     &           GDTB  , GDQB  ,
     &           DELP  , DPDP  , DPPMK ,
     &           GDQS  
      REAL       ADJDA ( IDIM*JDIM )         !" 今回調節されたか否か
      REAL       ADJDB ( IDIM*JDIM )         !" 前回調節されたか否か
      REAL       GDTB  ( IDIM*JDIM, KMAX )   !" 調節前の温度
      REAL       GDQB  ( IDIM*JDIM, KMAX )   !" 調節前の比湿
      REAL       DELP  ( IDIM*JDIM, KMAX )   !" ΔＰ
      REAL       DPDP  ( IDIM*JDIM, KMAX )   !" ファクター
      REAL       DPPMK ( IDIM*JDIM, KMAX )   !" ファクター
      REAL       GDQS  ( IDIM*JDIM, KMAX )   !" 飽和比湿
      REAL       GDPM34( IDIM*JDIM       )   !" P(K=3/4)
*
      INTEGER    IJ, K, ITR, I, J
      INTEGER    IFPAR, JFPAR
      REAL       ELF, AKAPPA
      REAL       ST, QEXE, QSAT, DQSAT, DQSATB, DT, STEXE
      REAL       DHDT, DHDTB, DTB, ADJPTS
*      
*   [INTERNAL FUNC]
#ifdef SYS_IBMS
      INCLUDE   (ZQSAT)                      !" 飽和水蒸気圧
#else
#include        "zqsat.F"                    !" 飽和水蒸気圧
#endif
*
*   [ONCE] 
      IF ( OFIRST ) THEN
         WRITE ( 6,* ) ' MOIST ADIABAT DATE=94/05/02 ishiwatari'
         OFIRST = .FALSE.
*
         CALL   REWNML ( IFPAR , JFPAR )
         READ   ( IFPAR, NMMADJ, END=190 )
  190    WRITE  ( JFPAR, NMMADJ )
*
      ENDIF
*
*"         < 2. ファクター計算 >
*
      ELF    = EL/CP
      AKAPPA = RAIR / CP
      DO 2000 IJ=1,IDIM*JDIM
         GDPM34( IJ ) = (   (   GDPM( IJ,1 ) ** (1. + AKAPPA )
     &                        - GDP ( IJ,1 ) ** (1. + AKAPPA ) )
     &                    / ( (   GDPM( IJ,1 ) - GDP( IJ,1 ) )
     &                      * (1. + AKAPPA) )
     &                  )
     &                   **( 1./AKAPPA )
 2000 CONTINUE 
*"             これインチキ. 本当は補間公式にすべし!?
*
c$$$      write(6,*) 'in p2adbt:k,gdp,gdpm'
c$$$      do 110 ij=1,idim*jdim
c$$$         do 100 k=1,5
c$$$            write(6,*) ij,k,gdp(ij,k),gdpm(ij,k)
c$$$ 100     continue 
c$$$            write(6,*) gdpm34( ij )
c$$$ 110  continue 
*
      DO 2100 K = 1, KMAX
         DO 2100 IJ = 1, IDIM*JDIM
            DELP  ( IJ,K ) = GDPM( IJ,K   ) - GDPM( IJ,K+1 )
 2100 CONTINUE
*
      DO 2300 IJ = 1, IDIM*JDIM
         DO 2200 K = 2, KMAX
            DPDP  ( IJ,K ) = DELP( IJ,K ) / DELP( IJ,K-1 )
*"              DPDP は後使ってない 
            DPPMK ( IJ,K ) = AKAPPA
     &                       * ( GDP( IJ,K-1 ) - GDP( IJ,K   ) )
     &                       /  GDPM( IJ,K   ) / 2.  
 2200    CONTINUE
         DPPMK ( IJ,1 ) = AKAPPA
     &                       * ( GDPM( IJ,1 ) - GDP( IJ,1   ) )
     &                       /  GDPM34( IJ ) / 2.  
 2300 CONTINUE
*
*"         < 3. T, q 構造決定>
*
*
*"         < 3.1 イテレーション(Newton 法) >
*
      DO 3100 IJ = 1, IDIM*JDIM
         DO 3110 K = 1, KMAX
c         DO 3110 K = 2, KMAX
            DO 3120 ITR = 1, ITRMAX
*
               IF(ITR .EQ. 1) THEN
                  IF( K .EQ. 1) THEN
                     GDT( IJ,1 ) = GDTG( IJ )
                  ELSE
                     GDT(IJ,K) = GDT(IJ,K-1)
                  END IF
               END IF
*
c$$$               if(k .eq. 1) then
c$$$                  write(6,*) 'in p2adbt:before setting qsat,dqsat'
c$$$                  write(6,*) 'gdt(ij,1), gdpm(ij,1),fqsat(Tg,Ps)'
c$$$                  write(6,*) gdt(ij,1),gdpm(ij,1),
c$$$     &                       fqsat(gdtg(ij),gdpm(ij,1))
c$$$               end if
               QSAT   = FQSAT ( GDT( IJ,K ), GDP( IJ,K ) )
               DQSAT  = FDQSAT( GDT( IJ,K   ), QSAT  )
*
c$$$                    write(6,*) 'IN ADJUSTMENT:bunshi,bunbo'
c$$$                    write(6,*) ( gdt(ij,k-1) - gdt(ij,k) + 
c$$$     &                           elf*(gdq(ij,k-1)-qsat)
c$$$     &               - dppmk(ij,K)*(GDT(IJ,K-1)+GDT(IJ,K)) ), 
c$$$     &              (1 + elf*dqsat + dppmk(IJ,K))
c$$$                    write(6,*) 'IN ADJUSTMENT:bunshi,bunbo'
c$$$                    write(6,*) GDT(IJ,K-1), GDT(IJ,K),
c$$$     &                           gdq(ij,k-1),qsat,
c$$$     &                         ELF*(GDQ(IJ,K-1)-QSAT),
c$$$     &               - DPPMK(IJ,K)*(GDT(IJ,K-1)+GDT(IJ,K)),
c$$$     &               elf*dqsat, dppmk(IJ,K)
*
               IF (K .EQ. 1) THEN
                  DT = (   GDTG( IJ ) - GDT( IJ,1 ) 
     &                   + ELF*( GDQG( IJ )-QSAT )
     &               - DPPMK( IJ,1 )*(GDTG( IJ )+GDT( IJ,1 )) )
     &              /(1 + ELF*DQSAT + DPPMK(IJ,1))
c$$$                  write(6,*) 'in p2adbt(k=1):ij,gdtg,gdt,gdqg,qsat'
c$$$                  write(6,*) 'ddpmk,dqsat'
c$$$                  write(6,*) 'del T, L*dq, dp*(Tg+T)'
c$$$                  write(6,*) 'bunbo'
c$$$                  write(6,*) ij,gdtg(ij),gdt(ij,1),gdqg(ij),qsat
c$$$                  write(6,*) dppmk(ij,1),dqsat
c$$$                  write(6,*) gdtg(ij)-gdt(ij,1),
c$$$     &                       elf*(gdqg(ij)-qsat),
c$$$     &                       -dppmk(ij,1)*(gdtg(ij)+gdt(ij,1))
c$$$                  write(6,*) 1+elf*dqsat+dppmk(ij,1)
               ELSE
                  DT = (   GDT( IJ,K-1 ) - GDT(IJ,K) 
     &                   + ELF*(GDQ( IJ,K-1 )-QSAT)
     &                   - DPPMK(IJ,K)*(GDT( IJ,K-1 )+GDT(IJ,K)) )
     &              /(1 + ELF*DQSAT + DPPMK(IJ,K))
               END IF
               GDT  ( IJ,K ) = GDT( IJ,K ) + DT
               GDQ  ( IJ,K ) = FQSAT ( GDT( IJ,K ), GDP( IJ,K ) )
*
c$$$               write(6,*) 'at the end of p2adbt:ij,k,dt,gdt,gdq'
c$$$               write(6,*) ij,k,dt,gdt( ij,k ), gdq( ij,k )
*
 3120       CONTINUE
 3110    CONTINUE
 3100 CONTINUE
*
      RETURN
      END
