Class ae_module
In: libsrc/ae_module/ae_module.f90

ae_module

Copyright :2001-2013 Shin-ichi Takehiro, Youhei SASAKI, SPMODEL Development Group
License :MIT/X11. Please see COPYRIGHT

概要

spml/ae_module モジュールは 1 次元周期境界条件の下での流体運動を実 フーリエ変換によるスペクトル法で数値計算するための Fortran90 関数 を提供する.

2 次元データの 1 次元に関して同時にスペクトル計算を実行するための 関数も提供しており, 2, 3 次元領域での計算のベースも提供する.

このモジュールは内部で ISPACK/FTPACK の Fortran77 サブルーチンを呼 んでいる. スペクトルデータの格納方法については ISPACK/FTPACK と異なっ ているので以下のコメントに注意されたい.

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

命名法

  • 関数名の先頭 (e_, g_, ae_, ag_) は, 返す値の形を示している.
    e_ :スペクトルデータ,
    g_ :1 次元格子点データ,
    ae_ :1 次元スペクトルデータが複数並んだ 2 次元データ,
    ag_ :1 次元格子点データが複数並んだ 2 次元データ.
  • 関数名の間の文字列(Dx)は, その関数の作用を表している.
  • 関数名の最後 (_e, _ae, _g, _ag) は, 入力変数の形スペクトルデータおよび 格子点データであることを示している.
    _e :スペクトルデータ
    _g :1 次元格子点データ
    _ae :1 次元スペクトルデータが複数並んだ 2 次元データ
    _ag :1 次元格子点データが複数並んだ 2 次元データ

各データの種類の説明

  • 以下, DP = kind(1.0D) である
  • g : 1 次元格子点データ.
    • 変数の種類と次元は real(DP), dimension(0:im-1).
    • im は X 座標の格子点数であり, サブルーチン ae_Initial にて あらかじめ設定しておく.
  • e : スペクトルデータ.
    • 変数の種類と次元は real(8), dimension(-km:km).
    • km は X 方向の最大波数であり, サブルーチン ae_Initial にて あらかじめ設定しておく.
    • スペクトルデータの格納のされ方は 0:km までが cos(kx) の係数, -km:-1 が sin(kx) の係数となっている.
  • ag : 1 次元(X)格子点データの並んだ 2 次元データ.
    • 変数の種類と次元は real(8), dimension(:,0:im-1). 第 2 次元が X 方向を表す.
  • ae : 1 次元スペクトルデータの並んだ 2 次元データ.
    • 変数の種類と次元は real(8), dimension(:,-km:km). 第 2 次元がスペクトルを表す.
  • g_ で始まる関数が返す値は 1 次元格子点データに同じ.
  • e_ で始まる関数が返す値はスペクトルデータに同じ.
  • ag_ で始まる関数が返す値は 1 次元格子点データの並んだ 2 次元データに同じ.
  • ae_ で始まる関数が返す値は 1 次元スペクトルデータの並んだ 2 次元データに同じ.
  • スペクトルデータに対する微分等の作用とは, 対応する格子点データに 微分などを作用させたデータをスペクトル変換したもののことである.

変数・手続き群の要約

初期化

ae_Initial :スペクトル変換の格子点数, 波数, 領域の大きさの設定

座標変数

g_X :格子点座標(X)を格納した 1 次元配列
g_X_Weight :重み座標を格納した 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 次元格子点データの積分および平均

補間

Interpolate_e :スペクトルデータからの特定の座標の値の補間(1 次元)
a_Interpolate_ae :スペクトルデータからの特定の座標の値の補間(2 次元)

Methods

Avr_g   Int_g   Interpolate_e   a_Avr_ag   a_Int_ag   a_Interpolate_ae   ae_Dx_ae   ae_Initial   ae_ag   ag_ae   e_Dx_e   e_g   g_e   g_x   g_x_weight  

Included Modules

dc_message dc_types

Public Instance methods

Avr_g( g ) result(Avr_g)
Function :
Avr_g :real(DP)
g :real(DP), dimension(0:im-1), intent(in)

1 次元格子点データの平均

[Source]

    function Avr_g(g)
      !
      ! 1 次元格子点データの平均
      !
      real(DP), dimension(0:im-1), intent(in)   :: g
      real(DP)                                  :: Avr_g

      Avr_g = Int_g(g)/sum(g_X_Weight)

    end function Avr_g
Int_g( g ) result(Int_g)
Function :
Int_g :real(DP)
: (out) 積分結果
g :real(DP), dimension(0:im-1), intent(in)
: (in) 格子点データ

1 次元格子点データの積分

[Source]

    function Int_g(g)
      !
      ! 1 次元格子点データの積分
      !
      real(DP), dimension(0:im-1), intent(in)   :: g      !(in) 格子点データ
      real(DP)                                  :: Int_g  !(out) 積分結果

      Int_g = sum(g*g_X_Weight)

    end function Int_g
Interpolate_e( e_data, xval ) result(Interpolate_e)
Function :
Interpolate_e :real(DP)
: 補間した結果の値
e_data :real(DP), dimension(-km:km), intent(IN)
: (in) 入力フーリエデータ
xval :real(DP), intent(in)
: 補間する点の座標

[Source]

    function Interpolate_e(e_data,xval)
      real(DP), dimension(-km:km), intent(IN) :: e_data
      !(in) 入力フーリエデータ
      real(DP), intent(in)                    :: xval
      ! 補間する点の座標
      real(DP)                                :: Interpolate_e
      ! 補間した結果の値
      integer :: k
      ! DO 文変数

      Interpolate_e = e_data(0)

      do k=1,km
        Interpolate_e = Interpolate_e + 2.0D0 * (   e_data(k)  * cos(2*PI*k/xl*xval) - e_data(-k) * sin(2*PI*k/xl*xval) )
      enddo

    end function Interpolate_e
a_Avr_ag( ag ) result(a_Avr_ag)
Function :
a_Avr_ag :real(DP), dimension(size(ag,1))
: (out) 平均結果
ag :real(DP), dimension(:,0:), intent(in)
: (in) 格子点データ

1 次元格子点データが並んだ 2 次元配列の平均

[Source]

    function a_Avr_ag(ag)
      !
      ! 1 次元格子点データが並んだ 2 次元配列の平均
      !
      real(DP), dimension(:,0:), intent(in)     :: ag        !(in) 格子点データ
      real(DP), dimension(size(ag,1))           :: a_Avr_ag  !(out) 平均結果

      a_Avr_ag = a_Int_ag(ag)/sum(g_X_Weight)

    end function a_Avr_ag
a_Int_ag( ag ) result(a_Int_ag)
Function :
a_Int_ag :real(DP), dimension(size(ag,1))
: (out) 積分結果
ag :real(DP), dimension(:,0:), intent(in)
: (in) 格子点データ

1 次元格子点データが並んだ 2 次元配列の積分

[Source]

    function a_Int_ag(ag)
      !
      ! 1 次元格子点データが並んだ 2 次元配列の積分
      !
      real(DP), dimension(:,0:), intent(in)     :: ag        !(in) 格子点データ
      real(DP), dimension(size(ag,1))           :: a_Int_ag  !(out) 積分結果
      integer :: i

      if ( size(ag,2) < im ) then
         call MessageNotify('E','ae_Int_ag', 'The Grid points of input data too small.')
      elseif ( size(ag,2) > im ) then
         call MessageNotify('W','ae_Int_ag', 'The Grid points of input data too large.')
      endif

      a_Int_ag = 0.0D0
      do i=0,im-1
         a_Int_ag(:) = a_Int_ag(:) + ag(:,i)*g_X_Weight(i)
      enddo
    end function a_Int_ag
a_Interpolate_ae( ae_data, xval ) result(a_Interpolate_ae)
Function :
a_Interpolate_ae :real(DP), dimension(size(ae_data,1))
: 補間した結果の値
ae_data :real(DP), dimension(:,-km:), intent(IN)
: (in) 入力フーリエデータ
xval :real(DP), intent(in)
: 補間する点の座標

[Source]

    function a_Interpolate_ae(ae_data,xval)
      real(DP), dimension(:,-km:), intent(IN)   :: ae_data
      !(in) 入力フーリエデータ
      real(DP), intent(in)                      :: xval
      ! 補間する点の座標
      real(DP), dimension(size(ae_data,1))      :: a_Interpolate_ae
      ! 補間した結果の値
      integer :: k
      ! DO 文変数

      a_Interpolate_ae = ae_data(:,0)

      do k=1,km
         a_Interpolate_ae = a_Interpolate_ae + 2.0D0 * (   ae_data(:,k)  * cos(2*PI*k/xl*xval) - ae_data(:,-k) * sin(2*PI*k/xl*xval) )
      enddo

    end function a_Interpolate_ae
ae_Dx_ae( ae ) result(ae_Dx_ae)
Function :
ae :real(DP), dimension(:,-km:), intent(in)
: (in) 入力スペクトルデータ

入力スペクトルデータに X 微分を作用する(2 次元データ).

スペクトルデータの X 微分とは, 対応する格子点データに X 微分を 作用させたデータのスペクトル変換のことである.

[Source]

    function ae_Dx_ae(ae)
      !
      ! 入力スペクトルデータに X 微分を作用する(2 次元データ).
      !
      ! スペクトルデータの X 微分とは, 対応する格子点データに X 微分を
      ! 作用させたデータのスペクトル変換のことである.
      !
      !
      real(DP), dimension(:,-km:), intent(in)  :: ae
      !(in)  入力スペクトルデータ
      real(DP), dimension(size(ae,1),-km:km)   :: ae_dx_ae
      !(out) 入力スペクトルデータの X 微分

      integer ::  k

      if ( size(ae,2) < 2*km+1 ) then
         call MessageNotify('W','ae_Dx_ae', 'The Fourier dimension of input data too small.')
      elseif ( size(ae,2) > 2*km+1 ) then
         call MessageNotify('W','ae_Dx_ae', 'The Fourier dimension of input data too large.')
      endif

      do k=-km,km
         ae_Dx_ae(:,k) = -(2*pi*k/xl)*ae(:,-k)
      enddo
    end function ae_dx_ae
ae_Initial( i, k, xmin, xmax )
Subroutine :
i :integer,intent(in)
: 格子点の数
k :integer,intent(in)
: 切断波数
xmin :real(DP),intent(in)
: X 座標範囲
xmax :real(DP),intent(in)
: X 座標範囲

スペクトル変換の格子点数, 波数, 領域の大きさを設定する.

他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで 初期設定をしなければならない.

[Source]

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

      integer :: ii

      im = i
      km = k
      xl = xmax-xmin

      if ( im <= 0 .or. km <= 0 ) then
         call MessageNotify('E','ae_Initial', 'Number of grid points and waves should be positive')
      elseif ( mod(im,2) /= 0 ) then
         call MessageNotify('E','ae_Initial', 'Number of grid points should be even')
      elseif ( km >= im/2 ) then
         call MessageNotify('E','ae_Initial', 'KM shoud be less than IM/2')
      endif

      allocate(ti(im*2))

      call fttrui(im,iti,ti)

      allocate(g_x(0:im-1))
      do ii=0,im-1
         g_X(ii) = xmin + xl/im*ii
      enddo

      allocate(g_x_weight(0:im-1))
      g_X_Weight = xl/im

      call MessageNotify( 'M','ae_initial','ae_module (2012/12/14) is initialized')

    end subroutine ae_Initial
ae_ag( ag ) result(ae_ag)
Function :
ae_ag :real(DP), dimension(size(ag,1),-km:km)
: (out) スペクトルデータ
ag :real(DP), dimension(:,:), intent(in)
: (in) 格子点データ

格子点データからスペクトルデータへ正変換する(2 次元データ用)

スペクトル正変換の定義は以下のとおり. 

  格子点データ     g_j (j=-0,...,im-1)
  スペクトルデータ e_k (k=-km,...,km)

  e_0 = (1/im)\sum_{j=0}^{im-1} g_j
  e_k = (1/im)\sum_{j=0}^{im-1} g_j \cos(2\pi jk/im)       (k=1,2,...,km)
  e_{-k} = - (1/im)\sum_{j=0}^{im-1} g_j \sin(2\pi jk/im)  (k=1,2,...,km)

[Source]

    function ae_ag(ag)
      !
      ! 格子点データからスペクトルデータへ正変換する(2 次元データ用)
      !
      ! スペクトル正変換の定義は以下のとおり.
      !
      !   格子点データ     g_j (j=-0,...,im-1)
      !   スペクトルデータ e_k (k=-km,...,km)
      !
      !   e_0 = (1/im)\sum_{j=0}^{im-1} g_j
      !   e_k = (1/im)\sum_{j=0}^{im-1} g_j \cos(2\pi jk/im)       (k=1,2,...,km)
      !   e_{-k} = - (1/im)\sum_{j=0}^{im-1} g_j \sin(2\pi jk/im)  (k=1,2,...,km)
      !
      real(DP), dimension(:,:), intent(in)   :: ag     !(in)  格子点データ
      real(DP), dimension(size(ag,1),-km:km) :: ae_ag  !(out) スペクトルデータ

      real(DP), dimension(size(ag,1)*im)     :: y
      real(DP), dimension(size(ag,1),0:im-1) :: ag_work
      integer :: m, k

      m = size(ag,1)
      if ( size(ag,2) < im ) then
         call MessageNotify('E','ae_ag', 'The Grid points of input data too small.')
      elseif ( size(ag,2) > im ) then
         call MessageNotify('W','ae_ag', 'The Grid points of input data too large.')
      endif
      ag_work = ag

      call fttruf(m,im,ag_work,y,iti,ti)

      do k=1,km
         ae_ag(:,k) = ag_work(:,2*k)
         ae_ag(:,-k) = ag_work(:,2*k+1)
      enddo
      ae_ag(:,0) = ag_work(:,0)

    end function ae_ag
ag_ae( ae ) result(ag_ae)
Function :
ag_ae :real(DP), dimension(size(ae,1),0:im-1)
: (out) 格子点データ
ae :real(DP), dimension(:,-km:), intent(in)
: (in) スペクトルデータ

スペクトルデータから格子点データへ逆変換する(2 次元データ用)

スペクトル逆変換の定義は以下のとおり. 

  スペクトルデータ e_k (k=-km,...,km)
  格子点データ     g_j (j=-0,...,im-1)

 g_j = e_0
     + 2\sum_{k=1}^{km}(e_k\cos(2\pi jk/im) - e_{-k}\sin(2\pi jk/im))
                                                   (j=0,1,...,im-1).

[Source]

    function ag_ae(ae)
      !
      ! スペクトルデータから格子点データへ逆変換する(2 次元データ用)
      !
      ! スペクトル逆変換の定義は以下のとおり.
      !
      !   スペクトルデータ e_k (k=-km,...,km)
      !   格子点データ     g_j (j=-0,...,im-1)
      !
      !  g_j = e_0
      !      + 2\sum_{k=1}^{km}(e_k\cos(2\pi jk/im) - e_{-k}\sin(2\pi jk/im))
      !                                                    (j=0,1,...,im-1).
      !
      real(DP), dimension(:,-km:), intent(in)  :: ae    !(in)  スペクトルデータ
      real(DP), dimension(size(ae,1),0:im-1)   :: ag_ae !(out) 格子点データ

      real(DP), dimension(size(ae,1)*im)       :: y
      integer :: m, k

      m=size(ae,1)
      if ( size(ae,2) < 2*km+1 ) then
         call MessageNotify('E','ag_ae', 'The Fourier dimension of input data too small.')
      elseif ( size(ae,2) > 2*km+1 ) then
         call MessageNotify('W','ag_ae', 'The Fourier dimension of input data too large.')
      endif

      ag_ae = 0.0D0
      ag_ae(:,0)=ae(:,0)
      ag_ae(:,1)=0.0D0
      do k=1,km
         ag_ae(:,2*k)   = ae(:,k)
         ag_ae(:,2*k+1) = ae(:,-k)
      enddo
      ag_ae(:,2*km+2:im-1)=0.0D0

      call fttrub(m,im,ag_ae,y,iti,ti)
    end function ag_ae
e_Dx_e( e ) result(e_Dx_e)
Function :
e_Dx_e :real(DP), dimension(-km:km)
: (out) 入力スペクトルデータの X 微分
e :real(DP), dimension(-km:km), intent(in)
: (in) 入力スペクトルデータ

入力スペクトルデータに X 微分を作用する(1 次元データ).

スペクトルデータの X 微分とは, 対応する格子点データに X 微分を 作用させたデータのスペクトル変換のことである.

[Source]

    function e_Dx_e(e)
      !
      ! 入力スペクトルデータに X 微分を作用する(1 次元データ).
      !
      ! スペクトルデータの X 微分とは, 対応する格子点データに X 微分を
      ! 作用させたデータのスペクトル変換のことである.
      !
      !
      real(DP), dimension(-km:km), intent(in)     :: e
      !(in)  入力スペクトルデータ
      real(DP), dimension(-km:km)                 :: e_Dx_e
      !(out) 入力スペクトルデータの X 微分
      real(DP), dimension(1,-km:km)               :: ae_work

      ae_work(1,:) = e
      ae_work = ae_Dx_ae(ae_work)
      e_Dx_e = ae_work(1,:)

    end function e_Dx_e
e_g( g ) result(e_g)
Function :
e_g :real(DP), dimension(-km:km)
: (out) スペクトルデータ
g :real(DP), dimension(0:im-1), intent(in)
: (in) 格子点データ

格子点データからスペクトルデータへ正変換する(1 次元データ用)

スペクトル正変換の定義は以下のとおり. 

  格子点データ     g_j (j=-0,...,im-1)
  スペクトルデータ e_k (k=-km,...,km)

  e_0 = (1/im)\sum_{j=0}^{im-1} g_j
  e_k = (1/im)\sum_{j=0}^{im-1} g_j \cos(2\pi jk/im)       (k=1,2,...,km)
  e_{-k} = - (1/im)\sum_{j=0}^{im-1} g_j \sin(2\pi jk/im)  (k=1,2,...,km)

[Source]

    function e_g(g)
      !
      ! 格子点データからスペクトルデータへ正変換する(1 次元データ用)
      !
      ! スペクトル正変換の定義は以下のとおり.
      !
      !   格子点データ     g_j (j=-0,...,im-1)
      !   スペクトルデータ e_k (k=-km,...,km)
      !
      !   e_0 = (1/im)\sum_{j=0}^{im-1} g_j
      !   e_k = (1/im)\sum_{j=0}^{im-1} g_j \cos(2\pi jk/im)       (k=1,2,...,km)
      !   e_{-k} = - (1/im)\sum_{j=0}^{im-1} g_j \sin(2\pi jk/im)  (k=1,2,...,km)
      !
      real(DP), dimension(-km:km)              :: e_g   !(out) スペクトルデータ
      real(DP), dimension(0:im-1), intent(in)  :: g     !(in)  格子点データ

      real(DP), dimension(1,size(g))        :: ag_work
      real(DP), dimension(1,-km:km)         :: ae_work

      ag_work(1,:) = g
      ae_work = ae_ag(ag_work)
      e_g = ae_work(1,:)

    end function e_g
g_e( e ) result(g_e)
Function :
g_e :real(DP), dimension(0:im-1)
: (out) 格子点データ
e :real(DP), dimension(-km:km), intent(in)
: (in) スペクトルデータ

スペクトルデータから格子点データへ逆変換する(1 次元データ用)

スペクトル逆変換の定義は以下のとおり. 

  スペクトルデータ e_k (k=-km,...,km)
  格子点データ     g_j (j=-0,...,im-1)

 g_j = e_0
     + 2\sum_{k=1}^{km}(e_k\cos(2\pi jk/im) - e_{-k}\sin(2\pi jk/im))
                                                   (j=0,1,...,im-1).

[Source]

    function g_e(e)
      !
      ! スペクトルデータから格子点データへ逆変換する(1 次元データ用)
      !
      ! スペクトル逆変換の定義は以下のとおり.
      !
      !   スペクトルデータ e_k (k=-km,...,km)
      !   格子点データ     g_j (j=-0,...,im-1)
      !
      !  g_j = e_0
      !      + 2\sum_{k=1}^{km}(e_k\cos(2\pi jk/im) - e_{-k}\sin(2\pi jk/im))
      !                                                    (j=0,1,...,im-1).
      !
      real(DP), dimension(0:im-1)             :: g_e  !(out) 格子点データ
      real(DP), dimension(-km:km), intent(in) :: e    !(in)  スペクトルデータ

      real(DP), dimension(1,size(e))  :: ae_work
      real(DP), dimension(1,0:im-1)   :: ag_work

      ae_work(1,:) = e
      ag_work = ag_ae(ae_work)
      g_e = ag_work(1,:)

    end function g_e
g_x
Variable :
g_x(:) :real(DP), allocatable
: 格子点座標(X)を格納した 1 次元配列.
g_x_weight
Variable :
g_x_weight(:) :real(DP), allocatable
: 重み座標を格納した 1 次元配列. X 方向の格子点の間隔が格納してある.