Class | wa_base_mpi_module_sjpack |
In: |
libsrc/wa_mpi_module_sjpack/wa_base_mpi_module_sjpack.f90
|
spml/wa_base_mpi_module_sjpack モジュールは球面上での流体運動を 球面調和函数を用いたスペクトル法と MPI 並列化によって 数値計算するための モジュール wa_mpi_module_sjpack の下部モジュールであり, スペクトル計算の基本的な Fortran90 関数を提供する. 球面上の 1 層モデル用 w_base_mpi_module_sjpack モジュールを多層モデル用に 拡張したものであり, 同時に複数個のスペクトルデータ, 格子点データに 対する変換が行える. 内部で ISPACK の SJPACK-MPI の Fortran77 サブルーチンを呼んでいる. スペクトルデータおよび格子点データの格納方法や変換の詳しい計算法に ついては ISPACK/SJPACK のマニュアルを参照されたい.
Variable : | |||
im =64 : | integer
|
Original external subprogram is w_base_mpi_module_sjpack#im
Variable : | |||
nn =22 : | integer
|
Original external subprogram is w_base_mpi_module_sjpack#nn
Variable : | |||
nn =22 : | integer
|
Original external subprogram is w_base_module_sjpack#nn
Subroutine : | |||
wa_Psi(:,:) : | real(8), intent(in)
| ||
wa_Chi(:,:) : | real(8), intent(in)
| ||
xva_U(0:im-1,1:jc,size(wa_Psi,2)) : | real(8), intent(out)
| ||
xva_V(0:im-1,1:jc,size(wa_Psi,2)) : | real(8), intent(out)
|
流線・ポテンシャル(スペクトルデータ)から速度場(格子データ)に (逆)変換する(多層用)
スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ.
u cosφ = ∂χ/∂λ - cosφ∂ψ/∂φ, v cosφ = cosφ∂χ/∂φ + ∂ψ/∂λ
subroutine wa_StreamPotential2VectorMPI(wa_Psi, wa_Chi, xva_U, xva_V) ! ! 流線・ポテンシャル(スペクトルデータ)から速度場(格子データ)に ! (逆)変換する(多層用) ! ! スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ. ! ! u cosφ = ∂χ/∂λ - cosφ∂ψ/∂φ, ! v cosφ = cosφ∂χ/∂φ + ∂ψ/∂λ ! real(8), intent(in) :: wa_Psi(:,:) !(in) 流線関数((nm+1)*(nm+1),:) real(8), intent(in) :: wa_Chi(:,:) !(in) 速度ポテンシャル((nm+1)*(nm+1),:) real(8), intent(out) :: xva_U(0:im-1,1:jc,size(wa_Psi,2)) !(out) 速度経度成分(0:im-1,1:jc,:) real(8), intent(out) :: xva_V(0:im-1,1:jc,size(wa_Psi,2)) !(out) 速度緯度成分(0:im-1,1:jc,:) integer :: k do k=1,size(wa_Psi,2) call w_StreamPotential2VectorMPI ( wa_Psi(:,k), wa_Chi(:,k), xva_U(:,:,k), xva_V(:,:,k) ) enddo end subroutine wa_StreamPotential2VectorMPI
Subroutine : | |||
xva_U(0:,:,:) : | real(8), intent(in)
| ||
xva_V(0:,:,:) : | real(8), intent(in)
| ||
wa_Vor((nm+1)*(nm+1),size(xva_U,3)) : | real(8), intent(out)
| ||
wa_Div((nm+1)*(nm+1),size(xva_U,3)) : | real(8), intent(out)
|
速度場(格子データ)から渦度・発散(スペクトルデータ)に (正)変換する(多層用)
スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ.
ζ = 1/cosφ∂v/∂λ - 1/cosφ ∂(u cosφ)/∂φ D = 1/cosφ∂u/∂λ + 1/cosφ ∂(v cosφ)/∂φ
subroutine wa_Vector2VorDivMPI(xva_U, xva_V, wa_Vor, wa_Div) ! ! 速度場(格子データ)から渦度・発散(スペクトルデータ)に ! (正)変換する(多層用) ! ! スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ. ! ! ζ = 1/cosφ∂v/∂λ - 1/cosφ ∂(u cosφ)/∂φ ! D = 1/cosφ∂u/∂λ + 1/cosφ ∂(v cosφ)/∂φ ! real(8), intent(in) :: xva_U(0:,:,:) !(in) 速度経度成分(0:im-1,1:jm,:) real(8), intent(in) :: xva_V(0:,:,:) !(in) 速度緯度成分(0:im-1,1:jm,:) real(8), intent(out) :: wa_Vor((nm+1)*(nm+1),size(xva_U,3)) !(out) 流線関数 real(8), intent(out) :: wa_Div((nm+1)*(nm+1),size(xva_U,3)) !(out) 速度ポテンシャル integer :: k do k=1,size(xva_U,3) call w_Vector2VorDivMPI (xva_U(:,:,k), xva_V(:,:,k), wa_Vor(:,k), wa_Div(:,k)) enddo end subroutine wa_Vector2VorDivMPI
Subroutine : | |||
xva_UCosLat(0:,:,:) : | real(8), intent(in)
| ||
xva_VCosLat(0:,:,:) : | real(8), intent(in)
| ||
wa_Vor((nm+1)*(nm+1),size(xva_UCosLat,3)) : | real(8), intent(out)
| ||
wa_Div((nm+1)*(nm+1),size(xva_UCosLat,3)) : | real(8), intent(out)
|
速度場(格子データ)から渦度・発散(スペクトルデータ)に (正)変換する(多層用)
スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ.
ζ = 1/cosφ∂v/∂λ - 1/cosφ ∂(u cosφ)/∂φ D = 1/cosφ∂u/∂λ + 1/cosφ ∂(v cosφ)/∂φ
subroutine wa_VectorCosLat2VorDivMPI(xva_UCosLat, xva_VCosLat, wa_Vor, wa_Div) ! ! 速度場(格子データ)から渦度・発散(スペクトルデータ)に ! (正)変換する(多層用) ! ! スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ. ! ! ζ = 1/cosφ∂v/∂λ - 1/cosφ ∂(u cosφ)/∂φ ! D = 1/cosφ∂u/∂λ + 1/cosφ ∂(v cosφ)/∂φ ! real(8), intent(in) :: xva_UCosLat(0:,:,:) !(in) 速度経度成分 * cos(lat) (0:im-1,1:jm,:) real(8), intent(in) :: xva_VCosLat(0:,:,:) !(in) 速度緯度成分 * cos(lat) (0:im-1,1:jm,:) real(8), intent(out) :: wa_Vor((nm+1)*(nm+1),size(xva_UCosLat,3)) !(out) 流線関数 real(8), intent(out) :: wa_Div((nm+1)*(nm+1),size(xva_UCosLat,3)) !(out) 速度ポテンシャル integer :: k do k=1,size(xva_UCosLat,3) call w_VectorCosLat2VorDivMPI ( xva_UCosLat(:,:,k), xva_VCosLat(:,:,k), wa_Vor(:,k), wa_Div(:,k) ) enddo end subroutine wa_VectorCosLat2VorDivMPI
Subroutine : |
モジュールの終了処理(割り付け配列の解放)をおこなう.
このサブルーチンを単独で用いるのでなく, 上位サブルーチン wa_mpi_Finalize を使用すること.
subroutine wa_base_mpi_Finalize ! ! モジュールの終了処理(割り付け配列の解放)をおこなう. ! ! このサブルーチンを単独で用いるのでなく, ! 上位サブルーチン wa_mpi_Finalize を使用すること. ! call MessageNotify('M','wa_base_mpi_finalize', 'wa_base_mpi_module_sjpack is finalized (dummy, 2013/02/23)') end subroutine wa_base_mpi_Finalize
Subroutine : |
SNPACK 用 wa_base_mpi_initial の互換のためのダミーサブルーチン
subroutine wa_base_mpi_initial ! ! SNPACK 用 wa_base_mpi_initial の互換のためのダミーサブルーチン ! call MessageNotify('M','wa_base_mpi_initial', 'wa_base_mpi_module_sjpack is initialized (dummy, 2013/02/23)') end subroutine wa_base_mpi_initial
Function : | |||
wa_xva((nm+1)*(nm+1),size(xva_data,3)) : | real(8)
| ||
xva_data(:,:,:) : | real(8), intent(in)
| ||
ipow : | integer, intent(in), optional
| ||
iflag : | integer, intent(in), optional
|
格子点 -> 球面調和関数スペクトル
格子データからスペクトルデータへ(正)変換する(多層用).
function wa_xva(xva_data,ipow,iflag) ! 格子点 -> 球面調和関数スペクトル ! ! 格子データからスペクトルデータへ(正)変換する(多層用). ! real(8), intent(in) :: xva_data(:,:,:) !(in) 格子点データ(im,jm,*) real(8) :: wa_xva((nm+1)*(nm+1),size(xva_data,3)) !(out) スペクトルデータ integer, intent(in), optional :: ipow !(in) 変換時に同時に作用させる 1/cosφ の次数. 省略時は 0. integer, intent(in), optional :: iflag ! 変換の種類 ! 0 : 通常の正変換 ! 1 : 経度微分を作用させた正変換 ! -1 : 緯度微分を作用させた正変換 ! 2 : sinφを作用させた正変換 ! 省略時は 0. integer, parameter :: ipow_default = 0 ! スイッチデフォルト値 integer, parameter :: iflag_default = 0 ! スイッチデフォルト値 integer ipval, ifval integer k if (present(ipow)) then ipval = ipow else ipval = ipow_default endif if (present(iflag)) then ifval = iflag else ifval = iflag_default endif do k=1,size(xva_data,3) wa_xva(:,k) = w_xv(xva_data(:,:,k),iflag=ifval,ipow=ipval) enddo end function wa_xva
Function : | |||
xva_wa(im,jc,size(wa_data,2)) : | real(8)
| ||
wa_data(:,:) : | real(8), intent(in)
| ||
ipow : | integer, intent(in), optional
| ||
iflag : | integer, intent(in), optional
|
球面調和関数スペクトル -> 格子点
スペクトルデータから格子データへ変換する(多層用).
function xva_wa(wa_data,ipow,iflag) ! 球面調和関数スペクトル -> 格子点 ! ! スペクトルデータから格子データへ変換する(多層用). ! real(8), intent(in) :: wa_data(:,:) !(in) スペクトルデータ real(8) :: xva_wa(im,jc,size(wa_data,2)) !(out) 格子点データ integer, intent(in), optional :: ipow !(in) 作用させる 1/cosφ の次数. 省略時は 0. integer, intent(in), optional :: iflag !(in) 変換の種類 ! 0 : 通常の正変換 ! 1 : 経度微分を作用させた正変換 ! -1 : 緯度微分を作用させた正変換 ! 2 : sinφを作用させた正変換 ! 省略時は 0. ! integer, parameter :: ipow_default = 0 integer, parameter :: iflag_default = 0 integer ipval, ifval integer k if (present(ipow)) then ipval = ipow else ipval = ipow_default endif if (present(iflag)) then ifval = iflag else ifval = iflag_default endif do k=1,size(wa_data,2) xva_wa(:,:,k) = xv_w(wa_data(:,k),iflag=ifval,ipow=ipval) enddo end function xva_wa