*" PACKAGE AEQADM 
*"                  ver 1.00                   94/04/28 ishiwatari
*"                  SRC:=$GC/model/agcm5/src/admin/aeqadm.F
***********************************************************************
      PROGRAM GCM5                           
*
*"    湿潤多層全球モデル
*"    Arakawa & Suarez 鉛直差分 σ座標
*"    平衡モデル
*"    標準物理過程, 簡略放射
*
*   [PARAM]






* PACKAGE  ZXT10L32 !" 格子点数，波数 X-Z
*
      INTEGER     IMAX           !" 経度格子数
      INTEGER     JMAX           !" 緯度格子数
      INTEGER     IDIM           !" 経度配列サイズ
      INTEGER     JDIM           !" 緯度配列サイズ
      INTEGER     MMAX           !" 最大東西波数
      INTEGER     LMAX           !" 最大南北波数
      INTEGER     NMAX           !" 最大全波数
      INTEGER     MINT           !" 東西波数間隔
      INTEGER     KMAX           !" 鉛直層数
      INTEGER     KDIM           !" 鉛直配列サイズ
      INTEGER     NMDIM          !" 水平波数成分数
      INTEGER     JMXHF          !" JMAX/2+1
      INTEGER     NWORK          !" ワーク領域大きさ
*
      PARAMETER ( IMAX  = 32,     JMAX =  1,
     &            IDIM  = IMAX+1, JDIM = JMAX,
     &            MMAX  = 10,     LMAX =  0,    NMAX = MMAX, MINT = 1,
     &            KMAX  = 32,     KDIM = KMAX,
     &            NMDIM = (MMAX/MINT+1)*(2*(NMAX+1)-MMAX)
     &                  - (NMAX-LMAX)/MINT*(NMAX-LMAX+1),
     &            JMXHF = JMAX/2+1,
     &            NWORK = IDIM*JDIM*KMAX*12                       )

* PACKAGE ZHDIM  !" 文字列文字数 
*
*"  [HIS]  90/08/31(numaguti)
*
      INTEGER     NCC
      INTEGER     NDC
      INTEGER     NFILN
*
      PARAMETER ( NCC   = 16, NDC   = 64 ,
     &            NFILN = 38              )

* PACKAGE ZHDIM  !" 標準物理定数(common)  
*
*"  [HIS]  90/08/31(numaguti)
*
      ReaL*8       ER             !" 地球半径
      ReaL*8       GRAV                        !" 重力加速度
      ReaL*8       CP                          !" 大気定圧比熱
      ReaL*8       RAIR                        !" 大気気体定数
      ReaL*8       EL                          !" 水の凝結の潜熱
      ReaL*8       CPVAP                       !" 水蒸気定圧比熱
      ReaL*8       RVAP                        !" 水蒸気気体定数
      ReaL*8       DH2O                        !" 水の密度
      ReaL*8       EPSV                        !" 水と大気の分子量比ε
      ReaL*8       ES0                         !" ０℃飽和蒸気圧：Ｐａ
      ReaL*8       STB                         !" Stefan-Boltzman
      ReaL*8       FKARM                       !" Karman 定数
      ReaL*8       EPSVT                       !" 1/ε-1
*
      COMMON   /COMCON/
     &          ER    , GRAV  ,
     &          CP    , RAIR  ,
     &          EL    , CPVAP , RVAP  , DH2O  , EPSV  , ES0   ,
     &          STB   , FKARM , EPSVT

*
*   [VAR]
      ReaL*8       GAU   ( IDIM, JDIM, KMAX )  
      ReaL*8       GAV   ( IDIM, JDIM, KMAX )  
      ReaL*8       GAT   ( IDIM, JDIM, KMAX )  
      ReaL*8       GAPS  ( IDIM, JDIM       )  
      ReaL*8       GAQ   ( IDIM, JDIM, KMAX )  
      ReaL*8       GAVOR ( IDIM, JDIM, KMAX )  
      ReaL*8       GADIV ( IDIM, JDIM, KMAX )  
*"       : 格子点データ(t) <DYNMCS>  格子点データ(t+Δt)
*
      ReaL*8       GBU   ( IDIM, JDIM, KMAX )  
      ReaL*8       GBV   ( IDIM, JDIM, KMAX )  
      ReaL*8       GBT   ( IDIM, JDIM, KMAX )  
      ReaL*8       GBPS  ( IDIM, JDIM       )  
      ReaL*8       GBQ   ( IDIM, JDIM, KMAX )  
      ReaL*8       GBVOR ( IDIM, JDIM, KMAX )  
      ReaL*8       GBDIV ( IDIM, JDIM, KMAX )  
      ReaL*8       GBTG  ( IDIM, JDIM )        
*"       : 格子点データ(t-Δt) <DYNMCS>  格子点データ(t)
*
      ReaL*8       GTUA  ( IDIM, JDIM, KMAX )  
      ReaL*8       GTVA  ( IDIM, JDIM, KMAX )  
      ReaL*8       GTH   ( IDIM, JDIM, KMAX )  
      ReaL*8       GTR   ( IDIM, JDIM, KMAX )  
*"       : <PHYSCS> 格子点変化項(P) <DYNMCS> 格子点変化項(P+D)
*
      INTEGER    ISTEP                 
      INTEGER    ITA                   
      INTEGER    ITB                   
      INTEGER    IDATEA( 6 )           
      INTEGER    IDATEB( 6 )           
      ReaL*8       DELT                  
      ReaL*8       DELTP                 
      LOGICAL    OADVNC                
      LOGICAL    ORSTR                 
*"       : 時刻等
*
*   [ONCE]
      INTEGER    ITSTRT                
      INTEGER    ITEND                 
      INTEGER    IORSTR                
*"       : 実験管理パラメータ
*
      INTEGER    IDELT                 
      INTEGER    INISTP                
      ReaL*8       TFIL                  
*"       : 時間ステップパラメータ
*
*"    COMMON /COMCON/ ( include ZCCOM )
*     ReaL*8       CP                    
*     ReaL*8       EL                    
*     ReaL*8       ER                    
*     ReaL*8       GRAV                  
*     ReaL*8       RAIR                  
*     ReaL*8       RVAP                  
*     ReaL*8       CPVAP                 
*     ReaL*8       DH2O                  
*     ReaL*8       EPSV                  
*     ReaL*8       ES0                   
*     ReaL*8       STB                   
*     ReaL*8       FKARM                 
*"    COMMON /COMCON/ end
*"       : 物理パラメータ
*
      ReaL*8       ALAT  ( JDIM )        
      ReaL*8       DLAT  ( JDIM )        
      ReaL*8       ALON  ( IDIM )        
      ReaL*8       DLON  ( IDIM )        
      ReaL*8       SIG   ( KMAX )        
      ReaL*8       SIGM  (KMAX+1)        
      ReaL*8       DSIG  ( KMAX )        
      ReaL*8       DSIGM (KMAX+1)        
*"       : 座標値
*
*   [WORK]
      COMMON    /COMWRK/ WORK
      ReaL*8       WORK  ( NWORK )       
*
      integer k
*
*" << SETPUP : 初期設定 >>
*
      CALL CLCSTR ( 'SETUP' )
      CALL YPREP                       
*
      CALL SETPAR                      
     O         ( ITSTRT, ITEND , IORSTR,
     O           IDELT , INISTP, TFIL    )
*
      CALL SETCOR                      
     O         ( ALON  , DLON  ,
     O           ALAT  , DLAT  ,
     O           SIG   , DSIG  ,
     O           SIGM  , DSIGM  )
*
      CALL RDSTRT                      
     O         ( GAU   , GAV   , GAT   , GAPS  , GAQ   ,
     O           GAVOR , GADIV ,
     O           GBU   , GBV   , GBT   , GBPS  , GBQ   ,
     O           GBVOR , GBDIV ,
     O           ITA   , ITB   , IDATEA, IDATEB, ISTEP ,
     C           ITSTRT, ALAT  , DLAT                    )
*
      CALL ADMXMN 
     I         ( GAU   , GAV   , GAT   , GAPS  , GAQ   ,
     I           GADIV , GAVOR , ITA   , 'start GA'      )
      CALL ADMXMN 
     I         ( GBU   , GBV   , GBT   , GBPS  , GBQ   ,
     I           GBDIV , GBVOR , ITB   , 'start GB'      )
*
      CALL RDSTRG                      
     O         ( GBTG  ,
     I           GBT   ,
     C           ITB     )
*
      CALL SETTIM                      
     I         ( ITA   , ITB   , IDATEA, IDATEB, ISTEP  )
      CALL AHSTRG                      
*
      CALL CLCEND ( 'SETUP' )
*
*  << LOOP : <************* メイン・ループ >>
*
      ITA=40
      ITB=20
      ITEND=20
      INISTP = 0
*
 5000 CONTINUE
*





         ISTEP = ISTEP + 1
*
*     << START : ステップはじめの処理 >>
*
         CALL CLCSTR ( 'TIMCON' )
*
         CALL TIMSTP                   
     O         ( DELT  , OADVNC, ORSTR ,
     I           ITA   , ITB   , IDATEA, IDATEB, ISTEP  ,
     C           IDELT , INISTP, ITEND , IORSTR          )
         write (6,*) 'MAIN after timstp:oadvnc,itend,ita,itb,inistep'
         write (6,*) oadvnc,itend,ita,itb,inistp
*
         CALL SAVEGB                   
     I         ( GBU   , GBV   , GBT   , GBPS  , GBQ   ,
     I           GBVOR , GBDIV                           )
         CALL CLCEND ( 'TIMCON' )
*
*"    << INTEGR : 時間積分 >>
* 
c         DELTP = 2 * DELT 
*
         write (6,*) 'in aequib:before physcs'
         write (6,*) 'in aequib:before physcs'
         write (6,*) 'in aequib:before physcs'
         write (6,*) 'in aequib:before physcs'
         write (6,*) 'in aequib:before physcs'
         CALL PHYSCS                   
     O         ( GTUA  , GTVA  , GTH   , GTR   , 
     M           GBTG  ,
     I           GBU   , GBV   , GBT   , GBPS  , GBQ   ,
     I           GBVOR , GBDIV , 
     I           ITB   , IDATEB, DELTP , DELT , 
     C           ALON  , ALAT  ,
     C           SIG   , SIGM  , DSIG  , DSIGM          )
*
*
*"    << OUTPUT : データ出力 >>
c         IF ( OADVNC ) THEN
            CALL CLCSTR ( 'HISTIN' )
            CALL HISTIN ( GBT   , 'T'   ) 
            CALL HISTIN ( GBPS  , 'PS'  )
            CALL HISTIN ( GBQ   , 'Q'   )
            CALL HISTIN ( GBTG  , 'TG'  )
            CALL CLCEND ( 'HISTIN' )
*
            CALL CLCSTR ( 'HISTOU' )
            CALL HISTRP
            CALL HISTOU                
            CALL CLCEND ( 'HISTOU' )
c         ENDIF
*
*"    << NEXT : 次のステップへ >>
*
         CALL CLCSTR ( 'MISC' )
*
         CALL    ADVSTP                
     O         ( ITA   , ITB   , IDATEA, IDATEB,
     I           OADVNC                          )
         write (6,*) 'MAIN after advstp: ita,itb'
         write (6,*) ita,itb
*
c         IF ( ORSTR ) THEN
            CALL WRRSTR                
     I         ( GAU   , GAV   , GAT   , GAPS  , GAQ   ,
     I           GAVOR , GADIV ,
     I           GBU   , GBV   , GBT   , GBPS  , GBQ   ,
     I           GBVOR , GBDIV ,
     I           ITA   , ITB   , IDATEA, IDATEB, ISTEP  )
*
            CALL WRRSTG   
     I         ( GBTG  ,
     I           ITB   , IDATEB, ISTEP  )
c         ENDIF
*
         CALL CLCEND ( 'MISC' )
*
*" << LEND : メイン・ループ終わり *************> >>
*
c      IF ( .NOT. ( ITB .GE. ITEND )  ) GOTO  5000
*
      CALL CLCOUT                      
      CALL YFINE                       
*
      STOP
      END
***********************************************************************
      SUBROUTINE AHSTRG                
*
*   [INTERN PARAM]
      INTEGER    ISTYPL
      DATA       ISTYPL / 1 /
      ReaL*8       VMISS
*
      CALL GZDBGT ( 'MISS', VMISS )   
*
      CALL     HISTRG
     I       ( 'T     ', 'temperature         ' ,'K     ', 'ALEV',
     I           VMISS , VMISS , 2.5   , 10.   , ISTYPL  ,
     I         'CON'   , '   ' , 0     , 0     , 'X'     ,'(F12.3)' )
*
      CALL     HISTRG
     I       ( 'PS    ', 'surface pressure    ' ,'mb    ', 'ASFC',
     I           VMISS , VMISS , 2.5   , 10.   , ISTYPL  ,
     I         'CON'   , '   ' , 0     , 0     , 'X'     ,'(F12.3)' )
*
      CALL     HISTRG
     I       ( 'Q     ', 'specific humidity   ' ,'g/g   ', 'ALEV',
     I           VMISS , VMISS , 1.E-3 , 5.E-3 , ISTYPL  ,
     I         'CON'   , '   ' , 0     , 0     , 'X'     ,'(3PF12.3)' )
*
      CALL     HISTRG
     I       ( 'TG    ', 'ground temperature  ' ,'K     ', 'ASFC',
     I           VMISS , VMISS , 2.5   , 10.   , ISTYPL  ,
     I         'CON'   , '   ' , 0     , 0     , 'X'     ,'(F12.3)' )
*
      RETURN
      END
*************************************************************************
      SUBROUTINE ADMXMN        
     I         ( GDU   , GDV   , GDT   , GDPS  , GDQ   ,
     I           GDDIV , GDVOR , IT    , HLABEL          )
*
*   [PARAM]




* PACKAGE  ZXT10L32 !" 格子点数，波数 X-Z
*
      INTEGER     IMAX           !" 経度格子数
      INTEGER     JMAX           !" 緯度格子数
      INTEGER     IDIM           !" 経度配列サイズ
      INTEGER     JDIM           !" 緯度配列サイズ
      INTEGER     MMAX           !" 最大東西波数
      INTEGER     LMAX           !" 最大南北波数
      INTEGER     NMAX           !" 最大全波数
      INTEGER     MINT           !" 東西波数間隔
      INTEGER     KMAX           !" 鉛直層数
      INTEGER     KDIM           !" 鉛直配列サイズ
      INTEGER     NMDIM          !" 水平波数成分数
      INTEGER     JMXHF          !" JMAX/2+1
      INTEGER     NWORK          !" ワーク領域大きさ
*
      PARAMETER ( IMAX  = 32,     JMAX =  1,
     &            IDIM  = IMAX+1, JDIM = JMAX,
     &            MMAX  = 10,     LMAX =  0,    NMAX = MMAX, MINT = 1,
     &            KMAX  = 32,     KDIM = KMAX,
     &            NMDIM = (MMAX/MINT+1)*(2*(NMAX+1)-MMAX)
     &                  - (NMAX-LMAX)/MINT*(NMAX-LMAX+1),
     &            JMXHF = JMAX/2+1,
     &            NWORK = IDIM*JDIM*KMAX*12                       )

*
*   [INPUT] 
      ReaL*8       GDU   ( IDIM, JDIM, KMAX )  
      ReaL*8       GDV   ( IDIM, JDIM, KMAX )  
      ReaL*8       GDT   ( IDIM, JDIM, KMAX )  
      ReaL*8       GDPS  ( IDIM, JDIM       )  
      ReaL*8       GDQ   ( IDIM, JDIM, KMAX )  
      ReaL*8       GDVOR ( IDIM, JDIM, KMAX )  
      ReaL*8       GDDIV ( IDIM, JDIM, KMAX )  
      INTEGER    IT
      CHARACTER  HLABEL *(*)
*
      WRITE ( 6,* ) HLABEL, ' IT= ', IT
      CALL MAXMIN( GDU,   IDIM, JDIM, KMAX, IDIM, JDIM, 'U'   )
      CALL MAXMIN( GDV,   IDIM, JDIM, KMAX, IDIM, JDIM, 'V'   )
      CALL MAXMIN( GDT,   IDIM, JDIM, KMAX, IDIM, JDIM, 'T'   )
      CALL MAXMIN( GDPS,  IDIM, JDIM, 1,    IDIM, JDIM, 'PS'  )
      CALL MAXMIN( GDQ,   IDIM, JDIM, KMAX, IDIM, JDIM, 'Q'   )
      CALL MAXMIN( GDVOR, IDIM, JDIM, KMAX, IDIM, JDIM, 'VOR' )
      CALL MAXMIN( GDDIV, IDIM, JDIM, KMAX, IDIM, JDIM, 'DIV' )
*
      RETURN
      END
