| Class | eee_module |
| In: |
libsrc/eee_module/eee_module.f90
|
| Authors: | Shin-ichi Takehiro, Youhei SASAKI |
| Version: | $Id: eee_module.f90 590 2013-08-19 08:48:21Z uwabami $ |
| Copyright&License: | See COPYRIGHT |
spml/eee_module モジュールは周期境界条件の下での 3 次元矩形領域の 流体運動をスペクトル法により数値計算するための Fortran90 関数を 提供する.
内部で ISPACK/P3PACK の Fortran77 サブルーチンを呼んでいる. スペクトルデータおよび格子点データの格納方法については ISPACK/P3PACK のマニュアルを参照されたい.
| eee_ : | スペクトルデータ(第 1,2 3 次元がそれぞれ Z, Y,X 方向波数) |
| eee2 : | 2 つのスペクトルデータの並んだもの |
zyx_ :; 3 次元格子点データ(第 1,2 3 次元がそれぞれ Z, Y,X 方向の格子点)
| yx_ : | XY 方向 2 次元格子点データ |
| zy_ : | YZ 方向 2 次元格子点データ |
| zx_ : | XZ 方向 2 次元格子点データ |
| x_ : | X 方向 1 次元格子点データ |
| y_ : | Y 方向 1 次元格子点データ |
| z_ : | Z 方向 1 次元格子点データ |
| _eee : | スペクトルデータ |
| _eee_eee : | 2 つのスペクトルデータ |
| _eee2 : | 2 つのスペクトルデータの並んだもの |
| _zyx : | 3 次元格子点データ |
| _x : | X 方向 1 次元格子点データ |
| _y : | Y 方向 1 次元格子点データ. |
| eee_Initial : | スペクトル変換の格子点数, 切断波数の設定 |
| eee_ChangeResolution : | 解像度設定の変更 |
| x_X, y_Y, z_Z : | 格子点座標(X,Y座標)を格納した 1 次元配列 |
| x_X_Weight, y_Y_Weight, z_Z_Weight : | 重み座標を格納した 1 次元配列 |
| zyx_X, zyx_Y, zyx_Z : | 格子点データの XYZ 座標(X,Y,Z) (格子点データ型 3 次元配列) |
| zyx_eee : | スペクトルデータから格子データへの変換 |
| eee_zyx : | 格子データからスペクトルデータへの変換 |
| xyz_zyx : | 格子点データの添え字順序変換 |
| eee_Lapla_eee : | スペクトルデータにラプラシアンを作用させる |
| eee_LaplaInv_eee : | スペクトルデータにラプラシアンの逆変換を作用させる |
| eee_Dx_eee : | スペクトルデータに X 微分を作用させる |
| eee_Dy_eee : | スペクトルデータに Y 微分を作用させる |
| eee_Dz_eee : | スペクトルデータに Z 微分を作用させる |
| eee2_RotVelxVor_eee2 : | Euler 方程式の非線形項を計算する |
| eee_VorFromZeta_eee2 : | 渦度 2 成分(ζ_1, ζ_2)ら渦度の 1 成分を計算する |
| eee_VelFromZeta_eee2 : | 渦度 2 成分(ζ_1, ζ_2)ら速度の 1 成分を計算する |
| eee2_ZetaFromVor_eee_eee_eee : | 渦度から渦度 2 成分(ζ_1, ζ_2)を計算する |
| IntZYX_zyx, AvrZYX_zyx : | 3 次元格子点データの全領域積分および平均 |
| z_IntYX_zyx, z_AvrYX_zyx : | 3 次元格子点データの X,Y 方向積分および平均 |
| x_IntZY_zyx, x_AvrZY_zyx : | 3 次元格子点データの Y,Z 方向積分および平均 |
| y_IntZX_zyx, y_AvrZX_zyx : | 3 次元格子点データの Z,X 方向積分および平均 |
| zy_IntX_zyx, zy_AvrX_zyx : | 3 次元格子点データの X 方向積分および平均 |
| zx_IntY_zyx, zx_AvrY_zyx : | 3 次元格子点データの Y 方向積分および平均 |
| yx_IntZ_zyx, yx_AvrZ_zyx : | 3 次元格子点データの Z 方向積分および平均 |
| IntYX_yx, AvrYX_yx : | 2 次元(XY)格子点データの X,Y 方向積分および平均 |
| y_IntX_yx, y_AvrX_yx : | 2 次元(XY)格子点データの X 方向積分および平均 |
| x_IntY_yx, x_AvrY_yx : | 2 次元(XY)格子点データの Y 方向積分および平均 |
| IntZX_zx, AvrZX_zx : | 2 次元(ZX)格子点データの Z,X 方向積分および平均 |
| z_IntX_zx, z_AvrX_zx : | 2 次元(ZX)格子点データの X 方向積分および平均 |
| x_IntZ_zx, x_AvrZ_zx : | 2 次元(ZX)格子点データの Z 方向積分および平均 |
| IntZY_zy, AvrZY_zy : | 2 次元(YZ)格子点データの Y,Z 方向積分および平均 |
| y_IntZ_zy, y_AvrZ_zy : | 2 次元(YZ)格子点データの Z 方向積分および平均 |
| z_IntY_zy, z_AvrY_zy : | 2 次元(YZ)格子点データの Y 方向積分および平均 |
| IntX_x, AvrX_x : | 1 次元(X)格子点データの X 方向積分および平均 |
| IntY_y, AvrY_y : | 1 次元(Y)格子点データの Y 方向積分および平均 |
| IntZ_z, AvrZ_z : | 1 次元(Z)格子点データの Z 方向積分および平均 |
| EnergyHelicityFromZeta_eee2 : | 全エネルギーと全ヘリシティーを計算する. |
| ESpectralFromZeta : | エネルギースペクトルを計算する. |
| Function : | |||
| AvrX_x : | real(8)
| ||
| x : | real(8), dimension(0:im-1), intent(IN)
|
1 次元(X)格子点データの X 方向平均
実際には格子点データ各点毎に x_X_Weight をかけた総和を計算し, x_X_Weight の総和で割ることで平均している.
function AvrX_x(x)
!
! 1 次元(X)格子点データの X 方向平均
!
! 実際には格子点データ各点毎に x_X_Weight をかけた総和を計算し,
! x_X_Weight の総和で割ることで平均している.
!
real(8), dimension(0:im-1), intent(IN) :: x !(in) 1 次元格子点データ
real(8) :: AvrX_x !(out) 平均値
AvrX_x = IntX_x(x)/sum(x_X_weight)
end function AvrX_x
| Function : | |||
| AvrYX_yx : | real(8)
| ||
| yx : | real(8), dimension(0:jm-1,0:im-1), intent(IN)
|
2 次元格子点データの全領域平均
実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた 総和を計算し, x_X_Weight*y_Y_Weight の総和で割ることで平均している.
function AvrYX_yx(yx)
!
! 2 次元格子点データの全領域平均
!
! 実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた
! 総和を計算し, x_X_Weight*y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:jm-1,0:im-1), intent(IN) :: yx
!(in) 2 次元格子点データ
real(8) :: AvrYX_yx
!(out) 平均値
AvrYX_yx = IntYX_yx(yx)/(sum(x_X_weight)*sum(y_Y_weight))
end function AvrYX_yx
| Function : | |||
| AvrY_y : | real(8)
| ||
| y : | real(8), dimension(0:jm-1), intent(IN)
|
1 次元(Y)格子点データの Y 方向平均
実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し, y_Y_Weight の総和で割ることで平均している.
function AvrY_y(y)
!
! 1 次元(Y)格子点データの Y 方向平均
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し,
! y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:jm-1), intent(IN) :: y !(in) 1 次元格子点データ
real(8) :: AvrY_y !(out) 平均値
AvrY_y = IntY_y(y)/sum(y_Y_weight)
end function AvrY_y
| Function : | |||
| AvrZX_zx : | real(8)
| ||
| zx : | real(8), dimension(0:km-1,0:im-1), intent(IN)
|
2 次元(ZX)格子点データの全領域平均
実際には格子点データ各点毎に x_X_Weight, z_Z_Weight をかけた 総和を計算し, x_X_Weight*z_Z_Weight の総和で割ることで平均している.
function AvrZX_zx(zx)
!
! 2 次元(ZX)格子点データの全領域平均
!
! 実際には格子点データ各点毎に x_X_Weight, z_Z_Weight をかけた
! 総和を計算し, x_X_Weight*z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:im-1), intent(IN) :: zx
!(in) 2 次元(ZX)格子点データ
real(8) :: AvrZX_zx
!(out) 平均値
AvrZX_zx = IntZX_zx(zx)/(sum(x_X_weight)*sum(z_Z_weight))
end function AvrZX_zx
| Function : | |||
| AvrZYX_zyx : | real(8)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
2 次元格子点データの全領域平均.
実際には格子点データ各点毎に x_X_Weight, y_Y_Weight, z_Z_Weight を かけた総和を計算し, x_X_Weight*y_Y_Weight*z_Z_Weight の総和で割る ことで平均している.
function AvrZYX_zyx(zyx)
!
! 2 次元格子点データの全領域平均.
!
! 実際には格子点データ各点毎に x_X_Weight, y_Y_Weight, z_Z_Weight を
! かけた総和を計算し, x_X_Weight*y_Y_Weight*z_Z_Weight の総和で割る
! ことで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 3 次元格子点データ
real(8) :: AvrZYX_zyx
!(out) 積分値
AvrZYX_zyx = IntZYX_zyx(zyx)/(sum(x_X_weight)*sum(y_Y_weight)*sum(z_Z_weight))
end function AvrZYX_zyx
| Function : | |||
| AvrZY_zy : | real(8)
| ||
| zy : | real(8), dimension(0:km-1,0:jm-1), intent(IN)
|
2 次元(ZY)格子点データの全領域平均
実際には格子点データ各点毎に z_Z_Weight, y_Y_Weight をかけた 総和を計算し, z_Z_Weight*y_Y_Weight の総和で割ることで平均している.
function AvrZY_zy(zy)
!
! 2 次元(ZY)格子点データの全領域平均
!
! 実際には格子点データ各点毎に z_Z_Weight, y_Y_Weight をかけた
! 総和を計算し, z_Z_Weight*y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1), intent(IN) :: zy
!(in) 2 次元(ZY)格子点データ
real(8) :: AvrZY_zy
!(out) 平均値
AvrZY_zy = IntZY_zy(zy)/(sum(y_Y_weight)*sum(z_Z_weight))
end function AvrZY_zy
| Function : | |||
| AvrZ_z : | real(8)
| ||
| z : | real(8), dimension(0:km-1), intent(IN)
|
1 次元(Z)格子点データの Z 方向平均
実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算し, z_Z_Weight の総和で割ることで平均している.
function AvrZ_z(z)
!
! 1 次元(Z)格子点データの Z 方向平均
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算し,
! z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1), intent(IN) :: z !(in) 1 次元格子点データ
real(8) :: AvrZ_z !(out) 平均値
AvrZ_z = IntZ_z(z)/sum(z_Z_weight)
end function AvrZ_z
| Subroutine : | |||
| esp : | real(8), dimension(:), intent(OUT)
| ||
| eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in)
|
渦度成分(ζ_1, ζ_2)からエネルギースペクトルを計算する.
* esp のサイズで求められる波数範囲が定められる
* esp の総和が EFFromZeta で求められるエネルギー
(全領域平均値)に等しい
subroutine ESpectralFromZeta(esp,eee2)
!
! 渦度成分(ζ_1, ζ_2)からエネルギースペクトルを計算する.
!
! * esp のサイズで求められる波数範囲が定められる
! * esp の総和が EFFromZeta で求められるエネルギー
! (全領域平均値)に等しい
!
real(8), dimension(:), intent(OUT) :: esp
! エネルギースペクトル
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in) :: eee2
! 渦度成分(ζ_1, ζ_2)
integer kmax
! 作業変数
kmax=size(esp)
call p3espt(nm,mm,lm,kmax,eee2,esp)
end subroutine ESpectralFromZeta
| Function : | |||
| EnergyHelicityFromZeta_eee2 : | real(8), dimension(2)
| ||
| eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in)
|
渦度成分(ζ_1, ζ_2)から全領域平均エネルギーとヘリシティーを計算する.
function EnergyHelicityFromZeta_eee2(eee2)
!
! 渦度成分(ζ_1, ζ_2)から全領域平均エネルギーとヘリシティーを計算する.
!
real(8), dimension(2) :: EnergyHelicityFromZeta_eee2
! エネルギーヘリシティー
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in) :: eee2
! 渦度成分(ζ_1, ζ_2)
real(8) :: E, H ! エネルギー, ヘリシティー
call p3cnsv(nm,mm,lm,eee2,E,H)
EnergyHelicityFromZeta_eee2(1) = E
EnergyHelicityFromZeta_eee2(2) = H
end function EnergyHelicityFromZeta_eee2
| Function : | |||
| IntX_x : | real(8)
| ||
| x : | real(8), dimension(0:im-1)
|
1 次元(X)格子点データの X 方向積分
実際には格子点データ各点毎に x_X_Weight をかけた総和を計算している.
function IntX_x(x)
!
! 1 次元(X)格子点データの X 方向積分
!
! 実際には格子点データ各点毎に x_X_Weight をかけた総和を計算している.
!
real(8), dimension(0:im-1) :: x !(in) 1 次元格子点データ
real(8) :: IntX_x !(out) 積分値
IntX_x = sum(x*x_X_Weight)
end function IntX_x
| Function : | |||
| IntYX_yx : | real(8)
| ||
| yx : | real(8), dimension(0:jm-1,0:im-1), intent(IN)
|
2 次元(YX)格子点データの全領域積分
実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた 総和を計算している.
function IntYX_yx(yx)
!
! 2 次元(YX)格子点データの全領域積分
!
! 実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:jm-1,0:im-1), intent(IN) :: yx
!(in) 2 次元(YX)格子点データ
real(8) :: IntYX_yx
!(out) 積分値
integer :: i, j
! 作業変数
IntYX_yx = 0.0d0
do i=0,im-1
do j=0,jm-1
IntYX_yx = IntYX_yx + yx(j,i) * y_Y_Weight(j) * x_X_Weight(i)
enddo
enddo
end function IntYX_yx
| Function : | |||
| IntY_y : | real(8)
| ||
| y : | real(8), dimension(0:jm-1)
|
1 次元(Y)格子点データの Y 方向積分
実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.
function IntY_y(y)
!
! 1 次元(Y)格子点データの Y 方向積分
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.
!
real(8), dimension(0:jm-1) :: y !(in) 1 次元格子点データ
real(8) :: IntY_y !(out) 積分値
IntY_y = sum(y*y_Y_Weight)
end function IntY_y
| Function : | |||
| IntZX_zx : | real(8)
| ||
| zx : | real(8), dimension(0:km-1,0:im-1), intent(IN)
|
2 次元(ZX)格子点データの全領域積分および平均.
実際には格子点データ各点毎に z_Z_Weight, x_X_Weight をかけた 総和を計算している.
function IntZX_zx(zx)
!
! 2 次元(ZX)格子点データの全領域積分および平均.
!
! 実際には格子点データ各点毎に z_Z_Weight, x_X_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:im-1), intent(IN) :: zx
!(in) 2 次元(ZX)格子点データ
real(8) :: IntZX_zx
!(out) 積分値
integer :: i, k
! 作業変数
IntZX_zx = 0.0d0
do i=0,im-1
do k=0,km-1
IntZX_zx = IntZX_zx + zx(k,i) * x_X_Weight(i) * z_Z_Weight(k)
enddo
enddo
end function IntZX_zx
| Function : | |||
| IntZYX_zyx : | real(8)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
2 次元格子点データの全領域積分
実際には格子点データ各点毎に x_X_Weight, y_Y_Weight, z_Z_Weight をかけた 総和を計算している.
function IntZYX_zyx(zyx)
!
! 2 次元格子点データの全領域積分
!
! 実際には格子点データ各点毎に x_X_Weight, y_Y_Weight, z_Z_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 3 次元格子点データ
real(8) :: IntZYX_zyx
!(out) 積分値
integer :: i, j, k
! 作業変数
IntZYX_zyx = 0.0d0
do i=0,im-1
do j=0,jm-1
do k=0,km-1
IntZYX_zyx = IntZYX_zyx + zyx(k,j,i) * z_Z_Weight(k) * y_Y_Weight(j) * x_X_Weight(i)
enddo
enddo
end do
end function IntZYX_zyx
| Function : | |||
| IntZY_zy : | real(8)
| ||
| zy : | real(8), dimension(0:km-1,0:jm-1), intent(IN)
|
2 次元(ZY)格子点データの全領域積分および平均.
実際には格子点データ各点毎に z_Z_Weight, y_Y_Weight をかけた 総和を計算している.
function IntZY_zy(zy)
!
! 2 次元(ZY)格子点データの全領域積分および平均.
!
! 実際には格子点データ各点毎に z_Z_Weight, y_Y_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1), intent(IN) :: zy
!(in) 2 次元(ZY)格子点データ
real(8) :: IntZY_zy
!(out) 積分値
integer :: j, k
! 作業変数
IntZY_zy = 0.0d0
do j=0,jm-1
do k=0,km-1
IntZY_zy = IntZY_zy + zy(k,j) * y_Y_Weight(j) * z_Z_Weight(k)
enddo
enddo
end function IntZY_zy
| Function : | |||
| IntZ_z : | real(8)
| ||
| z : | real(8), dimension(0:km-1)
|
1 次元(Z)格子点データの Z 方向積分
実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算している.
function IntZ_z(z)
!
! 1 次元(Z)格子点データの Z 方向積分
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算している.
!
real(8), dimension(0:km-1) :: z !(in) 1 次元格子点データ
real(8) :: IntZ_z !(out) 積分値
IntZ_z = sum(z*z_Z_Weight)
end function IntZ_z
| Function : | |||
| eee2_RotVelxVor_eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2)
| ||
| eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in)
|
渦度 2 成分(ζ_1, ζ_2)のスペクトルデータから Euler 方程式の非線形項
▽x(u x ω)
の 2 成分を計算する.
(ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
function eee2_RotVelxVor_eee2(eee2)
!
! 渦度 2 成分(ζ_1, ζ_2)のスペクトルデータから Euler 方程式の非線形項
!
! ▽x(u x ω)
!
! の 2 成分を計算する.
!
! (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2) :: eee2_RotVelxVor_eee2
!(out) 非線形項のスペクトルデータの 2 つの成分
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in) :: eee2
!(in) 入力スペクトルデータ. 渦度の 2 成分(ζ_1, ζ_2)
real(8), dimension((2*lm+1)*(2*mm+1)*(2*nm+1)) :: ws
real(8), dimension(km*jm*im*4) :: wg
call p3elnl(nm,mm,lm,km,jm,im,eee2,eee2_RotVelxVor_eee2, ws,wg,itk,tk,itj,tj,iti,ti)
end function eee2_RotVelxVor_eee2
| Function : | |||
| eee2_ZetaFromVor_eee_eee_eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2)
| ||
| eee_1 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
| ||
| eee_2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
| ||
| eee_3 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
渦度のスペクトルから渦度 2 成分(ζ_1, ζ_2)を計算する. (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
function eee2_ZetaFromVor_eee_eee_eee(eee_1,eee_2,eee_3)
!
! 渦度のスペクトルから渦度 2 成分(ζ_1, ζ_2)を計算する.
!
! (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2) :: eee2_ZetaFromVor_eee_eee_eee
!(out) 渦度の 2 成分(ζ_1, ζ_2)
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee_1, eee_2, eee_3
!(in) 渦度スペクトルデータの各成分
integer :: l,m,n
do l=-lm,lm
do m=-mm,mm
do n=-nm,nm
if ( l /= 0 ) then
eee2_ZetaFromVor_eee_eee_eee(n,m,l,1) = eee_2(n,m,l)
eee2_ZetaFromVor_eee_eee_eee(n,m,l,2) = eee_3(n,m,l)
elseif( m /= 0 ) then
eee2_ZetaFromVor_eee_eee_eee(n,m,l,1) = eee_3(n,m,l)
eee2_ZetaFromVor_eee_eee_eee(n,m,l,2) = eee_1(n,m,l)
else
eee2_ZetaFromVor_eee_eee_eee(n,m,l,1) = eee_1(n,m,l)
eee2_ZetaFromVor_eee_eee_eee(n,m,l,2) = eee_2(n,m,l)
endif
enddo
enddo
enddo
end function eee2_ZetaFromVor_eee_eee_eee
| Subroutine : | |
| id : | integer, intent(in) |
解像度設定の変更. eee_Initial で設定する際に かえってくるオプショナル引数 id の値を用いる.
subroutine eee_ChangeResolution(id)
!
! 解像度設定の変更. eee_Initial で設定する際に
! かえってくるオプショナル引数 id の値を用いる.
!
integer, intent(in) :: id
if (id .gt. nparams .or. id .lt. 1) then
write(*,*)"id is invalid"
end if
im = params(id)%im
jm = params(id)%jm
km = params(id)%km
lm = params(id)%lm
mm = params(id)%mm
nm = params(id)%nm
itk => params(id)%itk
tk => params(id)%tk
itj => params(id)%itj
tj => params(id)%tj
iti => params(id)%iti
ti => params(id)%ti
x_X => params(id)%x_X
y_Y => params(id)%y_Y
z_Z => params(id)%z_Z
x_X_Weight => params(id)%x_X_Weight
y_Y_Weight => params(id)%y_Y_Weight
z_Z_Weight => params(id)%z_Z_Weight
zyx_X => params(id)%zyx_X
zyx_Y => params(id)%zyx_Y
zyx_Z => params(id)%zyx_Z
end subroutine eee_ChangeResolution
| Function : | |||
| eee_Dx_eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
入力スペクトルデータに X 微分(∂x)を作用する.
スペクトルデータの X 微分とは, 対応する格子点データに X 微分を 作用させたデータのスペクトル変換のことである.
実際にはスペクトルデータに X 方向波数 l をかけて 実部 <-> 虚部成分に入れ換える計算を行っている.
function eee_Dx_eee(eee)
!
! 入力スペクトルデータに X 微分(∂x)を作用する.
!
! スペクトルデータの X 微分とは, 対応する格子点データに X 微分を
! 作用させたデータのスペクトル変換のことである.
!
! 実際にはスペクトルデータに X 方向波数 l をかけて
! 実部 <-> 虚部成分に入れ換える計算を行っている.
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_Dx_eee
!(out) スペクトルデータの X 微分
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee
!(in) 入力スペクトルデータ
integer l,m,n
! 作業変数
do l=-lm,lm
do m=-mm,mm
do n=-nm,nm
eee_Dx_eee(n,m,l) = -l*eee(-n,-m,-l)
enddo
enddo
enddo
end function eee_Dx_eee
| Function : | |||
| eee_Dy_eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
入力スペクトルデータに Y 微分(∂y)を作用する.
スペクトルデータの Y 微分とは, 対応する格子点データに Y 微分を 作用させたデータのスペクトル変換のことである.
実際にはスペクトルデータに Y 方向波数 m をかけて 実部 <-> 虚部成分に入れ換える計算を行っている.
function eee_Dy_eee(eee)
!
! 入力スペクトルデータに Y 微分(∂y)を作用する.
!
! スペクトルデータの Y 微分とは, 対応する格子点データに Y 微分を
! 作用させたデータのスペクトル変換のことである.
!
! 実際にはスペクトルデータに Y 方向波数 m をかけて
! 実部 <-> 虚部成分に入れ換える計算を行っている.
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_Dy_eee
!(out) スペクトルデータの X 微分
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee
!(in) 入力スペクトルデータ
integer l,m,n
! 作業変数
do l=-lm,lm
do m=-mm,mm
do n=-nm,nm
eee_Dy_eee(n,m,l) = -m*eee(-n,-m,-l)
enddo
enddo
enddo
end function eee_Dy_eee
| Function : | |||
| eee_Dz_eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
入力スペクトルデータに Z 微分(∂z)を作用する.
スペクトルデータの Z 微分とは, 対応する格子点データに Z 微分を 作用させたデータのスペクトル変換のことである.
実際にはスペクトルデータに Z 方向波数 n をかけて 実部 <-> 虚部成分に入れ換える計算を行っている.
function eee_Dz_eee(eee)
!
! 入力スペクトルデータに Z 微分(∂z)を作用する.
!
! スペクトルデータの Z 微分とは, 対応する格子点データに Z 微分を
! 作用させたデータのスペクトル変換のことである.
!
! 実際にはスペクトルデータに Z 方向波数 n をかけて
! 実部 <-> 虚部成分に入れ換える計算を行っている.
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_Dz_eee
!(out) スペクトルデータの X 微分
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee
!(in) 入力スペクトルデータ
integer l,m,n
! 作業変数
do l=-lm,lm
do m=-mm,mm
do n=-nm,nm
eee_Dz_eee(n,m,l) = -n*eee(-n,-m,-l)
enddo
enddo
enddo
end function eee_Dz_eee
| Subroutine : | |||
| i : | integer,intent(in)
| ||
| j : | integer,intent(in)
| ||
| k : | integer,intent(in)
| ||
| l : | integer,intent(in)
| ||
| m : | integer,intent(in)
| ||
| n : | integer,intent(in)
| ||
| id : | integer, intent(out), optional
|
スペクトル変換の格子点数, 波数を設定する. 領域の範囲は [0,2pi]x[0,2pi]x[0,2pi] に決め打ち
他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで 初期設定をしなければならない.
オプショナル引数 id を用いて異なる解像度, 領域を同時に 扱うことができる. eee_Initial を解像度, 領域ごとに呼んで 引数 id をキープし, eee_ChangeResolution で切替える.
subroutine eee_Initial(i,j,k,l,m,n,id)
!
! スペクトル変換の格子点数, 波数を設定する.
! 領域の範囲は [0,2pi]x[0,2pi]x[0,2pi] に決め打ち
!
! 他の関数や変数を呼ぶ前に, 最初にこのサブルーチンを呼んで
! 初期設定をしなければならない.
!
! オプショナル引数 id を用いて異なる解像度, 領域を同時に
! 扱うことができる. eee_Initial を解像度, 領域ごとに呼んで
! 引数 id をキープし, eee_ChangeResolution で切替える.
!
integer,intent(in) :: i ! 格子点数(X)
integer,intent(in) :: j ! 格子点数(Y)
integer,intent(in) :: k ! 格子点数(Z)
integer,intent(in) :: l ! 切断波数(X)
integer,intent(in) :: m ! 切断波数(Y)
integer,intent(in) :: n ! 切断波数(X)
integer, intent(out), optional :: id ! 解像度領域情報番号
character(len=3) cid
integer :: ii, jj, kk
im = i
jm = j
km = k
lm = l
mm = m
nm = n
if ( nparams .ge. nparams_max ) then
call MessageNotify('W','eee_initial', 'too many call of eee_Initial, nothing was done.')
if ( present(id) ) id = -1
return
end if
nparams = nparams + 1
params(nparams)%im = im
params(nparams)%jm = jm
params(nparams)%km = km
params(nparams)%lm = lm
params(nparams)%mm = mm
params(nparams)%nm = nm
allocate(params(nparams)%itk(5))
allocate(params(nparams)%itj(5))
allocate(params(nparams)%iti(5))
allocate(params(nparams)%tk(km*2))
allocate(params(nparams)%tj(jm*2))
allocate(params(nparams)%ti(im*2))
allocate(params(nparams)%x_X(0:im-1))
allocate(params(nparams)%x_X_Weight(0:im-1))
allocate(params(nparams)%y_Y(0:jm-1))
allocate(params(nparams)%y_Y_Weight(0:jm-1))
allocate(params(nparams)%z_Z(0:km-1))
allocate(params(nparams)%z_Z_Weight(0:km-1))
allocate(params(nparams)%zyx_X(0:km-1,0:jm-1,0:im-1))
allocate(params(nparams)%zyx_Y(0:km-1,0:jm-1,0:im-1))
allocate(params(nparams)%zyx_Z(0:km-1,0:jm-1,0:im-1))
call eee_ChangeResolution(nparams)
call p3init(km,jm,im,itk,tk,itj,tj,iti,ti)
do ii=0,im-1
x_X(ii) = 2*pi/im*ii
enddo
x_X_Weight = 2*pi/im
do jj=0,jm-1
y_Y(jj) = 2*pi/jm*jj
enddo
y_Y_Weight = 2*pi/jm
do kk=0,km-1
z_Z(kk) = 2*pi/km*kk
enddo
z_Z_Weight = 2*pi/km
zyx_X = spread(spread(x_X,1,jm),1,km)
zyx_Y = spread(spread(y_Y,2,im),1,km)
zyx_Z = spread(spread(z_Z,2,jm),3,im)
if ( present(id) ) id = nparams
write(cid,'(I3)') nparams
call MessageNotify('M','eee_initial','eee_module (2009/07/31) is initialized')
call MessageNotify('M','eee_initial', 'Resolution ID is '//trim(adjustl(cid)))
end subroutine eee_initial
| Function : | |||
| eee_LaplaInv_eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
入力スペクトルデータに逆ラプラシアン(∂xx+∂yy+∂zz)**(-1)を作用する.
スペクトルデータの逆ラプラシアンとは, 対応する格子点データに 逆ラプラシアンを作用させたデータのスペクトル変換のことである.
実際にはスペクトルデータに全波数 (l**2 + m**2 + n**2) で割る 計算を行っている. l=m=n=0 成分には 0 を代入している.
function eee_LaplaInv_eee(eee)
!
! 入力スペクトルデータに逆ラプラシアン(∂xx+∂yy+∂zz)**(-1)を作用する.
!
! スペクトルデータの逆ラプラシアンとは, 対応する格子点データに
! 逆ラプラシアンを作用させたデータのスペクトル変換のことである.
!
! 実際にはスペクトルデータに全波数 (l**2 + m**2 + n**2) で割る
! 計算を行っている. l=m=n=0 成分には 0 を代入している.
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_LaplaInv_eee
!(out) スペクトルデータの逆ラプラシアン
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee
!(in) スペクトルデータ
integer l,m,n
do l=-lm,lm
do m=-mm,mm
do n=-nm,nm
if ( l.ne.0 .or. m.ne.0 .or. n.ne.0 ) then
eee_LaplaInv_eee(n,m,l) = -eee(n,m,l)/(l**2+m**2+n**2)
else
eee_LaplaInv_eee(n,m,l) = 0.0
endif
enddo
enddo
enddo
end function eee_LaplaInv_eee
| Function : | |||
| eee_Lapla_eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
入力スペクトルデータにラプラシアン(∂xx+∂yy+∂zz)を作用する.
スペクトルデータのラプラシアンとは, 対応する格子点データに ラプラシアンを作用させたデータのスペクトル変換のことである.
実際にはスペクトルデータに全波数 (l**2 + m**2 + n**2) をかける 計算を行っている.
function eee_Lapla_eee(eee)
!
! 入力スペクトルデータにラプラシアン(∂xx+∂yy+∂zz)を作用する.
!
! スペクトルデータのラプラシアンとは, 対応する格子点データに
! ラプラシアンを作用させたデータのスペクトル変換のことである.
!
! 実際にはスペクトルデータに全波数 (l**2 + m**2 + n**2) をかける
! 計算を行っている.
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_Lapla_eee
!(out) スペクトルデータのラプラシアン
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee
!(in) 入力スペクトルデータ
integer l,m,n
! 作業変数
do l=-lm,lm
do m=-mm,mm
do n=-nm,nm
eee_Lapla_eee(n,m,l) = -(l**2+m**2+n**2)*eee(n,m,l)
enddo
enddo
enddo
end function eee_Lapla_eee
| Function : | |||
| eee_VelFromZeta_eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in)
| ||
| isw : | integer, intent(IN)
|
渦度 2 成分(ζ_1, ζ_2)のスペクトルデータから速度スペクトルの 1 成分 を計算する. (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
function eee_VelFromZeta_eee2(eee2,isw)
!
! 渦度 2 成分(ζ_1, ζ_2)のスペクトルデータから速度スペクトルの 1 成分
! を計算する.
!
! (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_VelFromZeta_eee2
!(out) 非線形項のスペクトルデータの 2 つの成分
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in) :: eee2
!(in) 入力スペクトルデータ. 渦度の 2 成分(ζ_1, ζ_2)
integer, intent(IN) :: isw
!(in) 出力する渦度の成分のインデックス(1,2,3)
call p3getu(nm,mm,lm,eee2,eee_VelFromZeta_eee2,isw)
end function eee_VelFromZeta_eee2
| Function : | |||
| eee_VorFromZeta_eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| eee2 : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in)
| ||
| isw : | integer, intent(IN)
|
渦度 2 成分(ζ_1, ζ_2)のスペクトルデータから渦度のスペクトルの 1 成分 を計算する. (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
function eee_VorFromZeta_eee2(eee2,isw)
!
! 渦度 2 成分(ζ_1, ζ_2)のスペクトルデータから渦度のスペクトルの 1 成分
! を計算する.
!
! (ζ_1, ζ_2) と渦度 ω との関係は ISPACK/P3PACK のマニュアルを参照
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_VorFromZeta_eee2
!(out) 非線形項のスペクトルデータの 2 つの成分
real(8), dimension(-nm:nm,-mm:mm,-lm:lm,2), intent(in) :: eee2
!(in) 入力スペクトルデータ. 渦度の 2 成分(ζ_1, ζ_2)
integer, intent(IN) :: isw
!(in) 出力する渦度の成分のインデックス(1,2,3)
call p3geto(nm,mm,lm,eee2,eee_VorFromZeta_eee2,isw)
end function eee_VorFromZeta_eee2
| Function : | |||
| eee_zyx : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(in)
|
格子データからスペクトルデータへ変換する.
function eee_zyx(zyx)
!
! 格子データからスペクトルデータへ変換する.
!
real(8), dimension(-nm:nm,-mm:mm,-lm:lm) :: eee_zyx
!(out) スペクトルデータ
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(in) :: zyx
!(in) 格子点データ
real(8), dimension(km*jm*im) :: w
real(8), dimension(0:km-1,0:jm-1,0:im-1) :: zyx_tmp
! 作業領域
zyx_tmp = zyx
call p3g2sa(nm,mm,lm,km,jm,im,zyx_tmp,eee_zyx,w,itk,tk,itj,tj,iti,ti)
end function eee_zyx
| Function : | |||
| x_AvrY_yx : | real(8), dimension(0:im-1)
| ||
| yx : | real(8), dimension(0:jm-1,0:im-1), intent(IN)
|
2 次元格子点データの Y 方向平均
実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し, y_Y_Weight の総和で割ることで平均している.
function x_AvrY_yx(yx)
!
! 2 次元格子点データの Y 方向平均
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し,
! y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:jm-1,0:im-1), intent(IN) :: yx
!(in) 2 次元格子点データ
real(8), dimension(0:im-1) :: x_AvrY_yx
!(out) 平均された 1 次元(X)格子点データ
x_AvrY_yx = x_IntY_yx(yx)/sum(y_Y_weight)
end function x_AvrY_yx
| Function : | |||
| x_AvrZY_zyx : | real(8), dimension(0:im-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの Y,Z 方向平均
実際には格子点データ各点毎に y_Y_Weight, z_Z_Weight をかけた 総和を計算し, y_Y_Weight*z_Z_Weight の総和で割ることで平均している.
function x_AvrZY_zyx(zyx)
!
! 3 次元格子点データの Y,Z 方向平均
!
! 実際には格子点データ各点毎に y_Y_Weight, z_Z_Weight をかけた
! 総和を計算し, y_Y_Weight*z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:im-1) :: x_AvrZY_zyx
!(out) 積分された 1 次元(Y)格子点データ
x_AvrZY_zyx = x_IntZY_zyx(zyx)/(sum(y_Y_weight)*sum(z_Z_weight))
end function x_AvrZY_zyx
| Function : | |||
| x_AvrZ_zx : | real(8), dimension(0:im-1)
| ||
| zx : | real(8), dimension(0:km-1,0:im-1), intent(IN)
|
2 次元(ZX)格子点データの Z 方向平均
実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算し, z_Z_Weight の総和で割ることで平均している.
function x_AvrZ_zx(zx)
!
! 2 次元(ZX)格子点データの Z 方向平均
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算し,
! z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:im-1), intent(IN) :: zx
!(in) 2 次元(ZX)格子点データ
real(8), dimension(0:im-1) :: x_AvrZ_zx
!(out) 平均された 1 次元(X)格子点データ
x_AvrZ_zx = x_IntZ_zx(zx)/sum(z_Z_weight)
end function x_AvrZ_zx
| Function : | |||
| x_IntY_yx : | real(8), dimension(0:im-1)
| ||
| yx : | real(8), dimension(0:jm-1,0:im-1), intent(IN)
|
2 次元(YX)格子点データの Y 方向積分
実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.
function x_IntY_yx(yx)
!
! 2 次元(YX)格子点データの Y 方向積分
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.
!
real(8), dimension(0:jm-1,0:im-1), intent(IN) :: yx
!(in) 2 次元(YX)格子点データ
real(8), dimension(0:im-1) :: x_IntY_yx
!(out) 積分された 1 次元(X)格子点データ
integer :: j
! 作業変数
x_IntY_yx = 0.0d0
do j=0,jm-1
x_IntY_yx(:) = x_IntY_yx(:) + yx(j,:) * y_Y_Weight(j)
enddo
end function x_IntY_yx
| Function : | |||
| x_IntZY_zyx : | real(8), dimension(0:im-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの Y,Z 方向積分
実際には格子点データ各点毎に y_Y_Weight, z_Z_Weight をかけた 総和を計算している.
function x_IntZY_zyx(zyx)
!
! 3 次元格子点データの Y,Z 方向積分
!
! 実際には格子点データ各点毎に y_Y_Weight, z_Z_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:im-1) :: x_IntZY_zyx
!(out) 積分された 1 次元(Y)格子点データ
integer :: j, k
! 作業変数
x_IntZY_zyx = 0.0d0
do j=0,jm-1
do k=0,km-1
x_IntZY_zyx(:) = x_IntZY_zyx(:) + zyx(k,j,:) * y_Y_Weight(j)* z_Z_Weight(k)
enddo
enddo
end function x_IntZY_zyx
| Function : | |||
| x_IntZ_zx : | real(8), dimension(0:im-1)
| ||
| zx : | real(8), dimension(0:km-1,0:im-1), intent(IN)
|
2 次元(ZX)格子点データの Z 方向積分
実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算している.
function x_IntZ_zx(zx)
!
! 2 次元(ZX)格子点データの Z 方向積分
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算している.
!
real(8), dimension(0:km-1,0:im-1), intent(IN) :: zx
!(in) 2 次元(ZX)格子点データ
real(8), dimension(0:im-1) :: x_IntZ_zx
!(out) 積分された 1 次元(Y)格子点データ
integer :: k
! 作業変数
x_IntZ_zx = 0.0d0
do k=0,km-1
x_IntZ_zx(:) = x_IntZ_zx(:) + zx(k,:) * z_Z_Weight(k)
enddo
end function x_IntZ_zx
| Variable : | |||
| x_X_Weight => null() : | real(8), dimension(:), pointer
|
| Function : | |||
| xyz_zyx : | real(8), dimension(0:im-1,0:jm-1,0:km-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(in)
|
格子データの添え字順序変更
function xyz_zyx(zyx)
!
! 格子データの添え字順序変更
!
real(8), dimension(0:im-1,0:jm-1,0:km-1) :: xyz_zyx
!(out) 格子点データ
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(in) :: zyx
!(in) 格子点データ
integer :: i,j,k
do k=0,km-1
do j=0,jm-1
do i=0,im-1
xyz_zyx(i,j,k) = zyx(k,j,i)
enddo
enddo
enddo
end function xyz_zyx
| Function : | |||
| y_AvrX_yx : | real(8), dimension(0:jm-1)
| ||
| yx : | real(8), dimension(0:jm-1,0:im-1), intent(IN)
|
2 次元格子点データの X 方向平均
実際には格子点データ各点毎に x_X_Weight をかけた総和を計算し, x_X_Weight の総和で割ることで平均している.
function y_AvrX_yx(yx)
!
! 2 次元格子点データの X 方向平均
!
! 実際には格子点データ各点毎に x_X_Weight をかけた総和を計算し,
! x_X_Weight の総和で割ることで平均している.
!
real(8), dimension(0:jm-1,0:im-1), intent(IN) :: yx
!(in) 2 次元格子点データ
real(8), dimension(0:jm-1) :: y_AvrX_yx
!(out) 平均された 1 次元(Y)格子点データ
y_AvrX_yx = y_IntX_yx(yx)/sum(x_X_weight)
end function y_AvrX_yx
| Function : | |||
| y_AvrZX_zyx : | real(8), dimension(0:jm-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの X,Z 方向平均
実際には格子点データ各点毎に x_X_Weight, z_Z_Weight をかけた 総和を計算し, x_X_Weight*z_Z_Weight の総和で割ることで平均している.
function y_AvrZX_zyx(zyx)
!
! 3 次元格子点データの X,Z 方向平均
!
! 実際には格子点データ各点毎に x_X_Weight, z_Z_Weight をかけた
! 総和を計算し, x_X_Weight*z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:jm-1) :: y_AvrZX_zyx
!(out) 積分された 1 次元(Y)格子点データ
y_AvrZX_zyx = y_IntZX_zyx(zyx)/(sum(x_X_weight)*sum(z_Z_weight))
end function y_AvrZX_zyx
| Function : | |||
| y_AvrZ_zy : | real(8), dimension(0:jm-1)
| ||
| zy : | real(8), dimension(0:km-1,0:jm-1), intent(IN)
|
2 次元(ZY)格子点データの Z 方向平均
実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算し, z_Z_Weight の総和で割ることで平均している.
function y_AvrZ_zy(zy)
!
! 2 次元(ZY)格子点データの Z 方向平均
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算し,
! z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1), intent(IN) :: zy
!(in) 2 次元(ZY)格子点データ
real(8), dimension(0:jm-1) :: y_AvrZ_zy
!(out) 平均された 1 次元(X)格子点データ
y_AvrZ_zy = y_IntZ_zy(zy)/sum(z_Z_weight)
end function y_AvrZ_zy
| Function : | |||
| y_IntX_yx : | real(8), dimension(0:jm-1)
| ||
| yx : | real(8), dimension(0:jm-1,0:im-1), intent(IN)
|
2 次元(YX)格子点データの X 方向積分
実際には格子点データ各点毎に x_X_Weight をかけた総和を計算している.
function y_IntX_yx(yx)
!
! 2 次元(YX)格子点データの X 方向積分
!
! 実際には格子点データ各点毎に x_X_Weight をかけた総和を計算している.
!
real(8), dimension(0:jm-1,0:im-1), intent(IN) :: yx
!(in) 2 次元(YX)格子点データ
real(8), dimension(0:jm-1) :: y_IntX_yx
!(out) 積分された 1 次元(Y)格子点データ
integer :: i
! 作業変数
y_IntX_yx = 0.0d0
do i=0,im-1
y_IntX_yx(:) = y_IntX_yx(:) + yx(:,i) * x_X_Weight(i)
enddo
end function y_IntX_yx
| Function : | |||
| y_IntZX_zyx : | real(8), dimension(0:jm-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの X,Z 方向積分
実際には格子点データ各点毎に x_X_Weight, z_Z_Weight をかけた 総和を計算している.
function y_IntZX_zyx(zyx)
!
! 3 次元格子点データの X,Z 方向積分
!
! 実際には格子点データ各点毎に x_X_Weight, z_Z_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:jm-1) :: y_IntZX_zyx
!(out) 積分された 1 次元(Y)格子点データ
integer :: i, k
! 作業変数
y_IntZX_zyx = 0.0d0
do i=0,im-1
do k=0,km-1
y_IntZX_zyx(:) = y_IntZX_zyx(:) + zyx(k,:,i) * x_X_Weight(i)* z_Z_Weight(k)
enddo
enddo
end function y_IntZX_zyx
| Function : | |||
| y_IntZ_zy : | real(8), dimension(0:jm-1)
| ||
| zy : | real(8), dimension(0:km-1,0:jm-1), intent(IN)
|
2 次元(ZY)格子点データの Z 方向積分
実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算している.
function y_IntZ_zy(zy)
!
! 2 次元(ZY)格子点データの Z 方向積分
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1), intent(IN) :: zy
!(in) 2 次元(ZY)格子点データ
real(8), dimension(0:jm-1) :: y_IntZ_zy
!(out) 積分された 1 次元(Y)格子点データ
integer :: k
! 作業変数
y_IntZ_zy = 0.0d0
do k=0,km-1
y_IntZ_zy(:) = y_IntZ_zy(:) + zy(k,:) * z_Z_Weight(k)
enddo
end function y_IntZ_zy
| Variable : | |||
| y_Y_Weight => null() : | real(8), dimension(:), pointer
|
| Function : | |||
| yx_AvrZ_zyx : | real(8), dimension(0:jm-1,0:im-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの Z 方向平均
実際には格子点データ各点毎に z_Z_Weight をかけた 総和を計算し, z_Z_Weight の総和で割ることで平均している.
function yx_AvrZ_zyx(zyx)
!
! 3 次元格子点データの Z 方向平均
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた
! 総和を計算し, z_Z_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:jm-1,0:im-1) :: yx_AvrZ_zyx
!(out) 積分された 1 次元(YX)格子点データ
yx_AvrZ_zyx = yx_IntZ_zyx(zyx)/sum(z_Z_weight)
end function yx_AvrZ_zyx
| Function : | |||
| yx_IntZ_zyx : | real(8), dimension(0:jm-1,0:im-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの Z 方向積分
実際には格子点データ各点毎に z_Z_Weight をかけた 総和を計算している.
function yx_IntZ_zyx(zyx)
!
! 3 次元格子点データの Z 方向積分
!
! 実際には格子点データ各点毎に z_Z_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:jm-1,0:im-1) :: yx_IntZ_zyx
!(out) 積分された 1 次元(YX)格子点データ
integer :: k
! 作業変数
yx_IntZ_zyx = 0.0d0
do k=0,km-1
yx_IntZ_zyx(:,:) = yx_IntZ_zyx(:,:) + zyx(k,:,:) * z_Z_Weight(k)
enddo
end function yx_IntZ_zyx
| Function : | |||
| z_AvrX_zx : | real(8), dimension(0:km-1)
| ||
| zx : | real(8), dimension(0:km-1,0:im-1), intent(IN)
|
2 次元(ZX)格子点データの X 方向平均
実際には格子点データ各点毎に x_X_Weight をかけた総和を計算し, x_X_Weight の総和で割ることで平均している.
function z_AvrX_zx(zx)
!
! 2 次元(ZX)格子点データの X 方向平均
!
! 実際には格子点データ各点毎に x_X_Weight をかけた総和を計算し,
! x_X_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:im-1), intent(IN) :: zx
!(in) 2 次元(ZX)格子点データ
real(8), dimension(0:km-1) :: z_AvrX_zx
!(out) 平均された 1 次元(Z)格子点データ
z_AvrX_zx = z_IntX_zx(zx)/sum(x_X_weight)
end function z_AvrX_zx
| Function : | |||
| z_AvrYX_zyx : | real(8), dimension(0:km-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの X,Y 方向平均
実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた 総和を計算し, x_X_Weight*y_Y_Weight の総和で割ることで平均している.
function z_AvrYX_zyx(zyx)
!
! 3 次元格子点データの X,Y 方向平均
!
! 実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた
! 総和を計算し, x_X_Weight*y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1) :: z_AvrYX_zyx
!(out) 積分された 1 次元(Y)格子点データ
z_AvrYX_zyx = z_IntYX_zyx(zyx)/(sum(x_X_weight)*sum(y_Y_weight))
end function z_AvrYX_zyx
| Function : | |||
| z_AvrY_zy : | real(8), dimension(0:km-1)
| ||
| zy : | real(8), dimension(0:km-1,0:jm-1), intent(IN)
|
2 次元(ZY)格子点データの Y 方向平均
実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し, y_Y_Weight の総和で割ることで平均している.
function z_AvrY_zy(zy)
!
! 2 次元(ZY)格子点データの Y 方向平均
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し,
! y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1), intent(IN) :: zy
!(in) 2 次元(ZY)格子点データ
real(8), dimension(0:km-1) :: z_AvrY_zy
!(out) 平均された 1 次元(Z)格子点データ
z_AvrY_zy = z_IntY_zy(zy)/sum(y_Y_weight)
end function z_AvrY_zy
| Function : | |||
| z_IntX_zx : | real(8), dimension(0:km-1)
| ||
| zx : | real(8), dimension(0:km-1,0:im-1), intent(IN)
|
2 次元(ZX)格子点データの X 方向積分
実際には格子点データ各点毎に x_X_Weight をかけた総和を計算している.
function z_IntX_zx(zx)
!
! 2 次元(ZX)格子点データの X 方向積分
!
! 実際には格子点データ各点毎に x_X_Weight をかけた総和を計算している.
!
real(8), dimension(0:km-1,0:im-1), intent(IN) :: zx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1) :: z_IntX_zx
!(out) 積分された 1 次元(X)格子点データ
integer :: i
! 作業変数
z_IntX_zx = 0.0d0
do i=0,im-1
z_IntX_zx(:) = z_IntX_zx(:) + zx(:,i) * x_X_Weight(i)
enddo
end function z_IntX_zx
| Function : | |||
| z_IntYX_zyx : | real(8), dimension(0:km-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの X,Y 方向積分
実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた 総和を計算している.
function z_IntYX_zyx(zyx)
!
! 3 次元格子点データの X,Y 方向積分
!
! 実際には格子点データ各点毎に x_X_Weight, y_Y_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1) :: z_IntYX_zyx
!(out) 積分された 1 次元(Y)格子点データ
integer :: i, j
! 作業変数
z_IntYX_zyx = 0.0d0
do i=0,im-1
do j=0,jm-1
z_IntYX_zyx(:) = z_IntYX_zyx(:) + zyx(:,j,i) * x_X_Weight(i)* y_Y_Weight(j)
enddo
enddo
end function z_IntYX_zyx
| Function : | |||
| z_IntY_zy : | real(8), dimension(0:km-1)
| ||
| zy : | real(8), dimension(0:km-1,0:jm-1), intent(IN)
|
2 次元(ZY)格子点データの Y 方向積分
実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.
function z_IntY_zy(zy)
!
! 2 次元(ZY)格子点データの Y 方向積分
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1), intent(IN) :: zy
!(in) 2 次元格子点データ
real(8), dimension(0:km-1) :: z_IntY_zy
!(out) 積分された 1 次元(X)格子点データ
integer :: j
! 作業変数
z_IntY_zy = 0.0d0
do j=0,jm-1
z_IntY_zy(:) = z_IntY_zy(:) + zy(:,j) * y_Y_Weight(j)
enddo
end function z_IntY_zy
| Variable : | |||
| z_Z_Weight => null() : | real(8), dimension(:), pointer
|
| Function : | |||
| zx_AvrY_zyx : | real(8), dimension(0:km-1,0:im-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの Y 方向平均
実際には格子点データ各点毎に y_Y_Weight をかけた 総和を計算し, y_Y_Weight の総和で割ることで平均している.
function zx_AvrY_zyx(zyx)
!
! 3 次元格子点データの Y 方向平均
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた
! 総和を計算し, y_Y_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1,0:im-1) :: zx_AvrY_zyx
!(out) 積分された 1 次元(Y)格子点データ
zx_AvrY_zyx = zx_IntY_zyx(zyx)/sum(y_Y_weight)
end function zx_AvrY_zyx
| Function : | |||
| zx_IntY_zyx : | real(8), dimension(0:km-1,0:im-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの Y 方向積分
実際には格子点データ各点毎に y_Y_Weight をかけた 総和を計算している.
function zx_IntY_zyx(zyx)
!
! 3 次元格子点データの Y 方向積分
!
! 実際には格子点データ各点毎に y_Y_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1,0:im-1) :: zx_IntY_zyx
!(out) 積分された 1 次元(Y)格子点データ
integer :: j
! 作業変数
zx_IntY_zyx = 0.0d0
do j=0,jm-1
zx_IntY_zyx(:,:) = zx_IntY_zyx(:,:) + zyx(:,j,:) * y_Y_Weight(j)
enddo
end function zx_IntY_zyx
| Function : | |||
| zy_AvrX_zyx : | real(8), dimension(0:km-1,0:jm-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの X 方向平均
実際には格子点データ各点毎に x_X_Weight をかけた 総和を計算し, x_X_Weight の総和で割ることで平均している.
function zy_AvrX_zyx(zyx)
!
! 3 次元格子点データの X 方向平均
!
! 実際には格子点データ各点毎に x_X_Weight をかけた
! 総和を計算し, x_X_Weight の総和で割ることで平均している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1,0:jm-1) :: zy_AvrX_zyx
!(out) 積分された 2 次元(ZY)格子点データ
zy_AvrX_zyx = zy_IntX_zyx(zyx)/sum(x_X_weight)
end function zy_AvrX_zyx
| Function : | |||
| zy_IntX_zyx : | real(8), dimension(0:km-1,0:jm-1)
| ||
| zyx : | real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN)
|
3 次元格子点データの X 方向積分
実際には格子点データ各点毎に x_X_Weight をかけた 総和を計算している.
function zy_IntX_zyx(zyx)
!
! 3 次元格子点データの X 方向積分
!
! 実際には格子点データ各点毎に x_X_Weight をかけた
! 総和を計算している.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1), intent(IN) :: zyx
!(in) 2 次元格子点データ
real(8), dimension(0:km-1,0:jm-1) :: zy_IntX_zyx
!(out) 積分された 2 次元(ZY)格子点データ
integer :: i
! 作業変数
zy_IntX_zyx = 0.0d0
do i=0,im-1
zy_IntX_zyx(:,:) = zy_IntX_zyx(:,:) + zyx(:,:,i) * x_X_Weight(i)
enddo
end function zy_IntX_zyx
| Variable : | |||
| zyx_X => null() : | real(8), dimension(:,:,:), pointer
|
| Variable : | |||
| zyx_Y => null() : | real(8), dimension(:,:,:), pointer
|
| Variable : | |||
| zyx_Z => null() : | real(8), dimension(:,:,:), pointer
|
| Function : | |||
| zyx_eee : | real(8), dimension(0:km-1,0:jm-1,0:im-1)
| ||
| eee : | real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in)
|
スペクトルデータから格子データへ変換する.
function zyx_eee(eee)
!
! スペクトルデータから格子データへ変換する.
!
real(8), dimension(0:km-1,0:jm-1,0:im-1) :: zyx_eee
!(out) 格子点データ
real(8), dimension(-nm:nm,-mm:mm,-lm:lm), intent(in) :: eee
!(in) スペクトルデータ
real(8), dimension(km*jm*im) :: w
! 作業領域
call p3s2ga(nm,mm,lm,km,jm,im,eee,zyx_eee,w,itk,tk,itj,tj,iti,ti)
end function zyx_eee