Class | ae_module |
In: |
src/ae_module.f90
|
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 次元データ
* g : 1 次元格子点データ. 変数の種類と次元は real(8), dimension(0:im-1). im は X 座標の格子点数であり, サブルーチン ae_Initial にて あらかじめ設定しておく. * e : スペクトルデータ. 変数の種類と次元は real(8), dimension(-km:km). km は X 方向の最大波数であり, サブルーチン ae_Initial にて あらかじめ設定しておく. スペクトルデータの格納のされ方については... * 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 次元格子点データの積分および平均 |
Function : | |||
a_Avr_ag : | real(8), dimension(size(ag,1))
| ||
ag : | real(8), dimension(:,0:), intent(in)
|
1 次元格子点データが並んだ 2 次元配列の平均
function a_Avr_ag(ag) ! ! 1 次元格子点データが並んだ 2 次元配列の平均 ! real(8), dimension(:,0:), intent(in) :: ag !(in) 格子点データ real(8), dimension(size(ag,1)) :: a_Avr_ag !(out) 平均結果 a_Avr_ag = a_Int_ag(ag)/sum(g_X_Weight) end function a_Avr_ag
Function : | |||
a_Int_ag : | real(8), dimension(size(ag,1))
| ||
ag : | real(8), dimension(:,0:), intent(in)
|
1 次元格子点データが並んだ 2 次元配列の積分
function a_Int_ag(ag) ! ! 1 次元格子点データが並んだ 2 次元配列の積分 ! real(8), dimension(:,0:), intent(in) :: ag !(in) 格子点データ real(8), 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
Function : | |||
ae : | real(8), dimension(:,-km:), intent(in)
|
入力スペクトルデータに X 微分を作用する(2 次元データ).
スペクトルデータの X 微分とは, 対応する格子点データに X 微分を 作用させたデータのスペクトル変換のことである.
function ae_Dx_ae(ae) ! ! 入力スペクトルデータに X 微分を作用する(2 次元データ). ! ! スペクトルデータの X 微分とは, 対応する格子点データに X 微分を ! 作用させたデータのスペクトル変換のことである. ! ! real(8), dimension(:,-km:), intent(in) :: ae !(in) 入力スペクトルデータ real(8), 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
Subroutine : | |||
i : | integer,intent(in)
| ||
k : | integer,intent(in)
| ||
xmin : | real(8),intent(in)
| ||
xmax : | real(8),intent(in)
|
スペクトル変換の格子点数, 波数, 領域の大きさを設定する.
他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで 初期設定をしなければならない.
subroutine ae_Initial(i,k,xmin,xmax) ! ! スペクトル変換の格子点数, 波数, 領域の大きさを設定する. ! ! 他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで ! 初期設定をしなければならない. ! integer,intent(in) :: i ! 格子点の数 integer,intent(in) :: k ! 切断波数 real(8),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 end subroutine ae_Initial
Function : | |||
ae_ag : | real(8), dimension(size(ag,1),-km:km)
| ||
ag : | real(8), dimension(:,:), intent(in)
|
格子点データからスペクトルデータへ変換する(2 次元データ用)
function ae_ag(ag) ! ! 格子点データからスペクトルデータへ変換する(2 次元データ用) ! real(8), dimension(:,:), intent(in) :: ag !(in) 格子点データ real(8), dimension(size(ag,1),-km:km) :: ae_ag !(out) スペクトルデータ real(8), dimension(size(ag,1)*im) :: y real(8), 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
Function : | |||
ag_ae : | real(8), dimension(size(ae,1),0:im-1)
| ||
ae : | real(8), dimension(:,-km:), intent(in)
|
スペクトルデータから格子点データへ変換する(2 次元データ用)
function ag_ae(ae) ! ! スペクトルデータから格子点データへ変換する(2 次元データ用) ! real(8), dimension(:,-km:), intent(in) :: ae !(in) スペクトルデータ real(8), dimension(size(ae,1),0:im-1) :: ag_ae !(out) 格子点データ real(8), 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 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 call fttrub(m,im,ag_ae,y,iti,ti) end function ag_ae
Function : | |||
e_Dx_e : | real(8), dimension(-km:km)
| ||
e : | real(8), dimension(-km:km), intent(in)
|
入力スペクトルデータに X 微分を作用する(1 次元データ).
スペクトルデータの X 微分とは, 対応する格子点データに X 微分を 作用させたデータのスペクトル変換のことである.
function e_Dx_e(e) ! ! 入力スペクトルデータに X 微分を作用する(1 次元データ). ! ! スペクトルデータの X 微分とは, 対応する格子点データに X 微分を ! 作用させたデータのスペクトル変換のことである. ! ! real(8), dimension(-km:km), intent(in) :: e !(in) 入力スペクトルデータ real(8), dimension(-km:km) :: e_Dx_e !(out) 入力スペクトルデータの X 微分 real(8), 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
Function : | |||
e_g : | real(8), dimension(-km:km)
| ||
g : | real(8), dimension(0:im-1), intent(in)
|
格子 -> スペクトル
格子点データからスペクトルデータへ変換する(1 次元データ用)
function e_g(g) ! 格子 -> スペクトル ! ! 格子点データからスペクトルデータへ変換する(1 次元データ用) ! real(8), dimension(-km:km) :: e_g !(in) 格子点データ real(8), dimension(0:im-1), intent(in) :: g !(out) スペクトルデータ real(8), dimension(1,size(g)) :: ag_work real(8), 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
Function : | |||
g_e : | real(8), dimension(0:im-1)
| ||
e : | real(8), dimension(-km:km), intent(in)
|
スペクトルデータから格子点データへ変換する(1 次元データ用)
function g_e(e) ! ! スペクトルデータから格子点データへ変換する(1 次元データ用) ! real(8), dimension(0:im-1) :: g_e !(in) スペクトルデータ real(8), dimension(-km:km), intent(in) :: e !(out) 格子点データ real(8), dimension(1,size(e)) :: ae_work real(8), 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