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

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

目次


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

初期化 機能
at_Initial チェビシェフ変換の格子点数, 波数, 領域の大きさの設定
座標変数 機能
g_X 格子点座標(X)を格納した 1 次元配列.
g_X_Weigtht 重み座標を格納した 1 次元配列.
基本変換 機能
g_t, ag_at チェビシェフデータから格子データへの変換
t_g, at_ag 格子データからチェビシェフデータへの変換
微分 機能
t_Dx_t, at_Dx_at チェビシェフデータに X 微分を作用させる
境界値問題 機能
at_Boundaries_DD, at_Boundaries_DN, at_Boundaries_ND, at_Boundaries_NN ディリクレ, ノイマン境界条件の適用
at_BoundariesTau_DD, at_BoundariesTau_DN, at_BoundariesTau_ND, at_BoundariesTau_NN ディリクレ, ノイマン境界条件の適用(タウ法)
at_BoundariesGrid_DD, at_BoundariesGrid_DN, at_BoundariesGrid_ND, at_BoundariesGrid_NN ディリクレ, ノイマン境界条件の適用(選点法)
積分・平均 機能
a_Int_ag, a_Avr_ag 1 次元格子点データの並んだ 2 次元配列の積分および平均.
Int_g, Avr_g 1 次元格子点データの積分および平均.


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

命名法

各データの種類の説明


サブルーチンの説明

subroutine at_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. 備考
    他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで初期設定を しなければならない.

subroutine at_Boundaries_DD(at_data,values), subroutine at_Boundaries_DD(t_data,values)
subroutine at_Boundaries_DN(at_data,values), subroutine at_Boundaries_DN(t_data,values)
subroutine at_Boundaries_ND(at_data,values), subroutine at_Boundaries_ND(t_data,values)
subroutine at_Boundaries_NN(at_data,values), subroutine at_Boundaries_NN(t_data,values)

  1. 備考 : これらのサブルーチンは at_BoundariesTau_?? に同じ. 以前の版との整合性のために残しており, 将来廃止予定である.

subroutine at_BoundariesTau_DD(at_data,values), subroutine at_BoundariesTau_DD(t_data,values)
subroutine at_BoundariesTau_DN(at_data,values), subroutine at_BoundariesTau_DN(t_data,values)
subroutine at_BoundariesTau_ND(at_data,values), subroutine at_BoundariesTau_ND(t_data,values)
subroutine at_BoundariesTau_NN(at_data,values), subroutine at_BoundariesTau_NN(t_data,values)

  1. 機能 : ディリクレ・ノイマン型境界条件の適用.
    at_Boundaries_DD : 両境界での値を与える.
    at_Boundaries_DN : i=0 で値, i=im で勾配の値を与える.
    at_Boundaries_DN : i=0 で勾配の値, i=im で値を与える.
    at_Boundaries_NN : 両境界で勾配の値を与える.
  2. 引数の説明
    1 次元チェビシェフデータの並んだ 2 次元配列の場合
         real(8), dimension(:,0:km),intent(inout)       :: at_data ! 境界条件を適用するデータ(2 次元)
         real(8), dimension(size(at_data,1),2), intent(in), optional :: values  ! 境界値
        
    1 次元チェビシェフデータの場合
        real(8), dimension(0:km),intent(inout)       :: t_data    ! 境界条件を適用するデータ(1 次元)
        real(8), dimension(2), intent(in), optional  :: values          ! 境界値
        
  3. 備考
    チェビシェフ空間において境界条件を満たすべく高次の係数を定める方法をとっている(タウ法). 境界値を省略した場合には値が 0 に設定される.

subroutine at_BoundariesGrid_DD(at_data,values), subroutine at_BoundariesGrid_DD(t_data,values)
subroutine at_BoundariesGrid_DN(at_data,values), subroutine at_BoundariesGrid_DN(t_data,values)
subroutine at_BoundariesGrid_ND(at_data,values), subroutine at_BoundariesGrid_ND(t_data,values)
subroutine at_BoundariesGrid_NN(at_data,values), subroutine at_BoundariesGrid_NN(t_data,values)

  1. 機能 : ディリクレ・ノイマン型境界条件の適用.
    at_Boundaries_DD : 両境界での値を与える.
    at_Boundaries_DN : i=0 で値, i=im で勾配の値を与える.
    at_Boundaries_DN : i=0 で勾配の値, i=im で値を与える.
    at_Boundaries_NN : 両境界で勾配の値を与える.
  2. 引数の説明
    1 次元チェビシェフデータの並んだ 2 次元配列の場合
         real(8), dimension(:,0:km),intent(inout)       :: at_data ! 境界条件を適用するデータ(2 次元)
         real(8), dimension(size(at_data,1),2), intent(in), optional :: values  ! 境界値
        
    1 次元チェビシェフデータの場合
        real(8), dimension(0:km),intent(inout)       :: t_data    ! 境界条件を適用するデータ(1 次元)
        real(8), dimension(2), intent(in), optional  :: values          ! 境界値
        
  3. 備考
    鉛直実格子点空間において内部領域の値と境界条件を満たすように条件を課している(選点法). このルーチンを用いるためには at_Initial にて設定するチェビシェフ切断波数(km)と鉛直格子点数(im)を等しくしておく必要がある. 境界値を省略した場合には値が 0 に設定される.


変数の説明

g_X

  1. 説明 : 格子点座標(X)を格納した 1 次元配列.
  2. 変数の型
          real(8), dimension(0:im) :: g_X
        
  3. 備考
    km 次のチェビシェフ多項式の零点から定まる格子点.

g_X_Weigtht

  1. 説明 : 重み座標を格納した 1 次元配列.
  2. 変数の型
          real(8), dimension(0:im) :: g_X_Weigtht
        
  3. 備考
    g_X_Weight には各格子点における積分のための重みが格納してある.


各関数の説明

function g_t(t), ag_at(at)

  1. 機能 : チェビシェフデータから格子データへ変換する.
  2. 変数の型
          real(8), dimension(0:im)               :: g_t
          real(8), dimension(0:km), intent(in)   :: t
    
          real(8), dimension(size(at,1),0:im)     :: ag_at
          real(8), dimension(:,0:km), intent(in)  :: at
        
  3. 備考

function t_g(g), at_ag(ag)

  1. 機能 : 格子データからチェビシェフデータへ変換する.
  2. 変数の型
          real(8), dimension(0:km)              :: t_g
          real(8), dimension(0:im), intent(in)  :: g
    
          real(8), dimension(size(ag,1),0:km)     :: at_ag
          real(8), dimension(:,0:im), intent(in)  :: ag
        
  3. 備考

function t_Dx_t(t), function at_Dx_at

  1. 機能 : 入力チェビシェフデータに X 微分を作用する.
  2. 変数の型
          real(8), dimension(0:km)                 :: t_Dx_t
          real(8), dimension(0:km), intent(in)     :: t
    
          real(8), dimension(size(at,1),0:km)      :: at_dx_at
          real(8), dimension(:,0:km), intent(in)   :: at
        
  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/11/19 更新 (竹広真一)