* PACKAGE DADMN    !" 力学メインパート
*
*"  [HIS] 92/01/09 (momoko)
*"        95/08/14 (takepiro) 
*"        95/08/29 (takepiro) psi=const に対応
***********************************************************************
      SUBROUTINE DYNMCS
     M            (  PSI2  , ZETA2 , T2  , 
     M               PSI1  , ZETA1 , T1  , 
     I               HINTG , 
     C               DX    , DZ    , DT  )
*
* [PARAM]
#include        "zcdim.F"            !" 格子点数, 波数
#include        "zhdim.F"            !" 文字列文字数
*
* [MODIFY]
*" 格子点データ(t) <DYNMCS>  格子点データ(t+Δt)
      REAL PSI2  ( 0:NX+1 , 0:NZ+1 ) !" 流線関数(1STEP前の値)
      REAL ZETA2 ( 0:NX+1 , 0:NZ+1 ) !" 渦度(1STEP前の値)
      REAL T2    ( 0:NX+1 , 0:NZ+1 ) !" 温度(1STEP前の値)
*
*" 格子点データ(t-Δt) <DYNMCS>  格子点データ(t)
      REAL PSI1  ( 0:NX+1 , 0:NZ+1 ) !" 流線関数(2STEP前の値)
      REAL ZETA1 ( 0:NX+1 , 0:NZ+1 ) !" 渦度(2STEP前の値)
      REAL T1    ( 0:NX+1 , 0:NZ+1 ) !" 温度(2STEP前の値)
*
* [INPUT]
      CHARACTER  HINTG * (NCC)       !" 時間積分スイッチ
*
      REAL DX                        !" X 軸上の格子点間隔
      REAL DZ                        !" Z 軸上の格子点間隔
      REAL DT                        !" 時間間隔
*
* [INTERNAL WORK]
      REAL PSI3  ( 0:NX+1 , 0:NZ+1 ) !" 流線関数(現STEPの値)
      REAL ZETA3 ( 0:NX+1 , 0:NZ+1 ) !" 渦度(現STEPの値)
      REAL T3    ( 0:NX+1 , 0:NZ+1 ) !" 温度(現STEPの値)
*
      REAL TNDZ  ( 0:NX+1 , 0:NZ+1 ) !" 力学作業領域
      REAL TNDT  ( 0:NX+1 , 0:NZ+1 ) !" 力学作業領域
*
      REAL AJAC  ( 0:NX+1 , 0:NZ+1 ) !" 力学作業領域
      REAL ALAP  ( 0:NX+1 , 0:NZ+1 ) !" 力学作業領域
*
* [INTERNAL SAVE]
      REAL       RAYL           !" レイリー数
      REAL       PRND           !" プランドル数
      REAL       QINT           !" 内部加熱率
*
      LOGICAL OFIRST
      DATA    OFIRST  / .TRUE. /
*    
      SAVE  OFIRST , 
     &      RAYL   , PRND  , QINT
*
*"        < 0. パラメター設定 >
*
      IF ( OFIRST ) THEN
         OFIRST = .FALSE.
         CALL DPARM
     O        ( RAYL  , PRND  , QINT  )
*
         CALL DBNSET            !" 境界条件の設定
     I        ( PSI2 , 
     I          ZETA2, T2   ,
     I          DX   , DZ   )
      ENDIF
*
*"        < 1. 時間積分 >
*
      IF ( HINTG(1:1) .EQ. 'E' ) THEN
         CALL DEULER
     O        ( PSI3  , ZETA3 , T3  ,
     I          PSI2  , ZETA2 , T2  , 
     I          RAYL  , PRND  , QINT,
     C          DX    , DZ    , DT  ,
     W          AJAC  , ALAP          )
      ELSE IF ( HINTG(1:1) .EQ. 'H' ) THEN
         CALL DHEUN
     O        ( PSI3  , ZETA3 , T3  , 
     I          PSI2  , ZETA2 , T2  , 
     I          RAYL  , PRND  , QINT, 
     C          DX    , DZ    , DT  ,
     W          TNDZ  , TNDT  ,
     W          AJAC  , ALAP   ) 
      ELSE IF ( HINTG(1:1) .EQ. 'L' ) THEN 
         CALL DLPFRG
     O         ( PSI3  , ZETA3 , T3   ,
     I           PSI2  , ZETA2 , T2   , 
     I           PSI1  , ZETA1 , T1   , 
     I           RAYL  , PRND  , QINT ,
     C           DX    , DZ    , DT  ,
     W           AJAC  , ALAP   ) 
      ELSE IF ( HINTG(1:1) .EQ. 'R' ) THEN
         CALL DRUNGE
     O        ( PSI3  , ZETA3 , T3  ,
     I          PSI2  , ZETA2 , T2  , 
     I          RAYL  , PRND  , QINT,
     C          DX    , DZ    , DT  ,
     W          PSI1  , ZETA1 , T1  , 
     W          TNDZ  , TNDT  ,
     W          AJAC  , ALAP          )
      ELSE
         CALL MSGDMP( 'E', 'DYNMCS', 
     &                'SCHEME OF TIME INTEGRATION NOT SUPPORTED' )
      ENDIF
*
*"        < 2. PSI,T,ZETAのデータをずらす >
*
      CALL COPY( PSI1 , PSI2 , (NX+2)*(NZ+2) )
      CALL COPY( T1   , T2   , (NX+2)*(NZ+2) )
      CALL COPY( ZETA1, ZETA2, (NX+2)*(NZ+2) )
*
      CALL COPY( PSI2 , PSI3 , (NX+2)*(NZ+2) )
      CALL COPY( T2   , T3   , (NX+2)*(NZ+2) )
      CALL COPY( ZETA2, ZETA3, (NX+2)*(NZ+2) )
*
      RETURN
      END
***********************************************************************
      SUBROUTINE  DPARM                        !" 力学パラメター設定
     O             ( RAYL ,  PRND  , QINT )
*
*   [OUTPUT]
      REAL       RAYL           !" レイリー数
      REAL       PRND           !" プランドル数
      REAL       QINT           !" 内部加熱率
*
*   [NAMELIST]
      REAL       RAYLI, SRAYL   !" レイリー数
      REAL       PRNDL, SPRND   !" プランドル数
      REAL       QINTR, SQINT   !" 内部加熱率
*
*   [INTERNAL WORK]
      INTEGER    IFPAR, JFPAR
      LOGICAL    OFIRST
*
* [NAMELIST DEFALUT]
      NAMELIST    /NMDYN/  RAYLI, PRNDL, QINTR
      DATA  RAYLI , PRNDL , QINTR 
     &     /  0.0 ,   1.0 ,   0.0 /
*      
      DATA OFIRST  /.TRUE./
      SAVE SRAYL, SPRND, SQINT
*
*"         < 1. NAMELIST の読み込み >
*
      IF ( OFIRST )THEN
         OFIRST = .FALSE.
*
         CALL   REWNML ( IFPAR , JFPAR )
         READ   ( IFPAR, NMDYN , END=1190 )
 1190    WRITE  ( JFPAR, NMDYN  )
*
         SRAYL   = RAYLI
         SPRND   = PRNDL
         SQINT   = QINTR
      ENDIF
*
*"         < 2. パラメターの出力 >
*
      RAYL   = SRAYL
      PRND   = SPRND
      QINT   = SQINT
*
      RETURN
      END
