spml/ae_module レファレンスマニュアル

spml/ae_module モジュールは 1 次元周期境界条件の下での 流体運動をスペクトル法により数値計算するための Fortran90 関数を提供する. 2 次元データの 1 次元に関して同時にスペクトル計算を実行するための 関数も提供しており, 2, 3 次元領域での計算のベースも提供する. このモジュールは内部で ISPACK/FTPACK の Fortran77 サブルーチンを呼んでいる. スペクトルデータおよび格子点データの格納方法については ISPACK/FTPACK のマニュアルを参照されたい.

目次


サブルーチン・関数一覧(機能別)

初期化 機能
ae_initial スペクトル変換の格子点数, 波数, 領域の大きさの設定
座標変数 機能
g_X 格子点座標(X)を格納した 1 次元配列.
g_X_Weigtht 重み座標を格納した 1 次元配列.
基本変換 機能
g_e, ag_ae スペクトルデータから格子データへの変換
e_g, ae_ag 格子データからスペクトルデータへの変換
微分 機能
e_Dx_e, ae_Dx_ae スペクトルデータに X 微分を作用させる
積分・平均 機能
a_Int_ag, a_Avr_ag 1 次元格子点データの並んだ 2 次元配列の積分および平均.
Int_g, Avr_g 1 次元格子点データの積分および平均.


関数・変数の名前と型について

命名法

各データの種類の説明


サブルーチンの説明

subroutine ae_initial(i,k,xmin,xmax)

  1. 機能 : スペクトル変換の格子点数, 波数, 領域の大きさを設定する.
  2. 引数の説明
        integer,intent(in) :: i              ! 格子点の設定(X)
        integer,intent(in) :: k              ! 切断波数の設定(X)
        real(8),intent(in) :: xmin, xmax     ! X 座標の範囲
        
  3. 備考
    他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで初期設定を しなければならない.

変数の説明

g_X

  1. 説明 : 格子点座標(X)を格納した 1 次元配列.
  2. 変数の型
          real(8), dimension(0:im-1) :: g_X
        
  3. 備考

g_X_Weigtht

  1. 説明 : 重み座標を格納した 1 次元配列.
  2. 変数の型
          real(8), dimension(0:im-1) :: g_X_Weigtht
        
  3. 備考
    g_X_Weight には X 方向の格子点の間隔が格納してある.


各関数の説明

function g_e(e), ag_ae(ae)

  1. 機能 : スペクトルデータから格子データへ変換する.
  2. 変数の型
          real(8), dimension(0:im-1)             :: g_e
          real(8), dimension(-km:km), intent(in) :: e
    
          real(8), dimension(size(ae,1),0:im-1)     :: ag_ae
          real(8), dimension(:,-km:km), intent(in)  :: ae
        
  3. 備考

function e_g(g), ae_ag(ag)

  1. 機能 : 格子データからスペクトルデータへ変換する.
  2. 変数の型
          real(8), dimension(-km:km)              :: e_g
          real(8), dimension(0:im-1), intent(in)  :: g
    
          real(8), dimension(size(ag,1),-km:km)     :: ae_ag
          real(8), dimension(:,0:im-1), intent(in)  :: ag
        
  3. 備考

function e_Dx_e(e), function ae_Dx_ae

  1. 機能 : 入力スペクトルデータに X 微分を作用する.
  2. 変数の型
          real(8), dimension(-km:km)                 :: e_Dx_e
          real(8), dimension(-km:km), intent(in)     :: e
    
          real(8), dimension(size(ae,1),-km:km)      :: ae_dx_ae
          real(8), dimension(:,-km:km), intent(in)   :: ae
        
  3. 備考
    スペクトルデータの X 微分とは, 対応する格子点データに X 微分を作用させたデータのスペクトル変換のことである.

function a_Int_ag(xy), function a_Avr_ag(xy)

  1. 機能 : 1 次元格子点データが並んだ 2 次元配列の積分および平均.
  2. 変数の型
          real(8), dimension(:,0:im-1), intent(in) :: ag
          real(8), dimension(size(ag,1))           :: a_Int_ag, a_Avr_ag
        
  3. 備考

function Int_g(x), function Avr_g(x)

  1. 機能 : 1 次元格子点データの積分および平均.
  2. 変数の型
          real(8), dimension(0:im-1), intent(in)   :: g
          real(8)                                  :: Int_g, Avr_g
        
  3. 備考


地球流体電脳倶楽部 SPMODEL プロジェクト
spmodel@gfd-dennou.org

2002/08/18 作成 (竹広真一)
2002/08/20 更新 (竹広真一)