! -*- mode: f90; coding: utf-8 -*- !---------------------------------------------------------------------- ! Copyright (c) 2019 SPMODEL Development Group. All rights reserved. ! ! 履歴 2016/03/09 竹広真一 ! 2019/10/04 佐々木洋平 ! 2019/10/09 佐々木洋平 ! 2020/11/06 竹広真一 セクター計算オプション導入 !---------------------------------------------------------------------- !> @author Shin-ichi Takehiro, Youhei SASAKI !> @copyright Copyright (C) SPMODEL Development Group, 2019. !> License is MIT/X11. see [COPYRIGHT](@ref COPYRIGHT) in detail !> @brief 2次元球面: 球面調和関数展開 !> !> @details !> spml/w_module モジュールは球面上での 2 次元流体運動を !> 球面調和函数を用いたスペクトル法によって数値計算するための !> Fortran90 関数を提供する. !> w_module は実際には以下の下部モジュールからなる． !> * w_module_base: 基本変換 !> * w_module_deriv: 微分計算 !> * w_module_integral: 積分・平均計算 !> * w_module_spectrum: スペクトル解析 !> * w_module_interpolate: 補間計算 !> !> 内部で ISPACK3 の SXPACKサブルーチンを呼んでいる. スペクトルデータお !> よび格子点データの格納方法や変換の詳しい計算法については !> ISPACK3/SXPACK のマニュアルを参照されたい. !> !> @warning !> 後方非互換である点に注意されたい !> * 旧版と異なり (ISPACK/SXPACKのレベルで) 球面調和関数の正規化が異なる． !> * 旧版と異なり, !> [l_nm](@ref w_module_base#l_nm), !> [nm_l](@ref w_module_base#nm_l) !> は初期化したのちにしか使うことができない. !> !> ## 関数・変数の名前と型について !> !> ### 命名法 !> !> * 関数名の先頭 (`w_`, `nm_`, `n_`, `xy_`, `x_`, `y_`) は, !> 返す値の形を示している. !> * `w_` : スペクトルデータ !> * `xy_` : 2 次元格子点データ !> * `nm_` : スペクトルデータの並んだ 3 次元配列 !> * `n_` : スペクトルデータの並んだ 2 次元配列 !> * `x_` : 経度方向 1 次元格子点データ !> * `y_` : 緯度方向 1 次元格子点データ !> !> * 関数名の間の文字列(`DLon`, `GradLat`, `GradLat`, `DivLon`, `DivLat`, !> `Lapla`, `LaplaInv`, `Jacobian`)は, その関数の作用を表している. !> !> * 関数名の最後 (`_w_w`, `_w`, `_xy`, `_x`, `_y`) は, !> 入力変数の形スペクトルデータおよび格子点データであることを示している. !> * `_w` :: スペクトルデータ !> * `_w_w` :: 2 つのスペクトルデータ !> * `_xy` :: 2 次元格子点データ !> * `_x` :: 経度方向 1 次元格子点データ !> * `_y` :: 緯度方向 1 次元格子点データ !> !> ### 各データの種類の説明 !> !> * 以下, `DP = kind(1.0D0)` !> * `xy` : 2 次元格子点データ. !> * 変数の種類と次元は `real(DP), dimension(0:im-1,1:jm)`. !> * `im`, `jm` はそれぞれ経度, 緯度座標の格子点数であり, サブルーチン !> [w_Initial](@ref w_initial) にてあらかじめ設定しておく. !> !> * `w` : スペクトルデータ. !> * 変数の種類と次元は `real(DP), dimension((nm+1)*(nm+1))`. !> * nm は球面調和函数の最大全波数であり, サブルーチン w_Initial にて !> あらかじめ設定しておく. !> * スペクトルデータの格納のされ方は関数 !> [l_nm](@ref w_module_base#l_nm), !> [nm_l](@ref w_module_base#nm_l) によって !> 調べることができる. !> !> * `nm` : スペクトルデータの並んだ 2 次元配列. !> * 変数の種類と次元は `real(DP), dimension(0:nm,-nm:nm)`. !> 第 1 次元が水平全波数, 第 2 次元が帯状波数を表す. !> * `nm` は球面調和函数の最大全波数であり, !> サブルーチン [w_Initial](@ref w_initial) !> にてあらかじめ設定しておく. !> !> * `n` : スペクトルデータの並んだ 1 次元配列. !> * 変数の種類と次元は `real(DP), dimension(0:nm)`. !> * 第 1 次元が水平全波数を表す. nm は球面調和函数の最大全波数であり, !> サブルーチン [w_Initial](@ref w_initial) !> にてあらかじめ設定しておく. !> !> * `x`, `y` : 経度, 緯度方向 1 次元格子点データ. !> * 変数の種類と次元はそれぞれ `real(DP), dimension(0:im-1)` !> および `real(DP), dimension(1:jm)`. !> !> * `w_` で始まる関数が返す値はスペクトルデータに同じ. !> !> * `xy_` で始まる関数が返す値は 2 次元格子点データに同じ. !> !> * `x_`, `y_` で始まる関数が返す値は 1 次元格子点データに同じ. !> !> * スペクトルデータに対する微分等の作用とは, 対応する格子点データに !> 微分などを作用させたデータをスペクトル変換したものことである. !> !> ## 変数・手続き群の要約 !> !> ### 初期化, 終了処理 !> !> w_module_base にて定義 !> !> * [w_Initial] !> (@ref w_initial) !> : スペクトル変換の格子点数, 波数, 領域の大きさの設定 !> * [w_Finalize] !> (@ref w_finalize) !> : モジュールの終了処理(割り付け配列の解放)をおこなう. !> !> ### 座標変数 !> !> * [x_Lon] !> (@ref w_module_base#x_lon) !> : 経度座標格子点を格納した 1 次元配列 !> * [y_Lat] !> (@ref w_module_base#y_lat) !> : 緯度座標格子点を格納した 1 次元配列 !> * [x_Lon_Weight] !> (@ref w_module_base#x_lon_weight) !> : 重み座標を格納した 1 次元配列 !> * [y_Lat_Weight] !> (@ref w_module_base#y_lat_weight) !> : 重み座標を格納した 1 次元配列 !> * [xy_Lon] !> (@ref w_module_base#xy_lon) !> : 経度座標(X,Y)の格子点データ(2 次元配列) !> * [xy_Lat] !> (@ref w_module_base#xy_lat) !> : 緯度座標(X,Y)の格子点データ(2 次元配列) !> !> ### 基本変換 !> !> * [xy_w] !> (@ref w_module_base#xy_w) !> : スペクトルデータから格子データへの変換 !> * [w_xy] !> (@ref w_module_base#w_xy) !> : 格子データからスペクトルデータへの変換 !> * [l_nm] !> (@ref w_module_base#l_nm), !> [nm_l] !> (@ref w_module_base#nm_l) !> : スペクトルデータの格納位置と全波数・帯状波数の変換 !> * [w_StreamPotential2Vector] !> (@ref w_module_base#w_streampotential2vector) !> : 流線ポテンシャルから速度場計算 !> * [w_Vector2VorDiv] !> (@ref w_module_base#w_vector2vordiv) !> : 速度場から渦度発散を計算 !> * [w_VectorCosLat2VorDiv] !> (@ref w_module_base#w_vectorcoslat2vordiv) !> : 速度場から渦度発散を計算 !> !> ### 微分 !> !> w_module_deriv にて定義 !> !> * [w_Lapla_w] !> (@ref w_module_deriv#w_lapla_w) !> : スペクトルデータにラプラシアンを作用させる !> * [rn] !> (@ref w_module_deriv#rn) !> : スペクトルデータのラプラシアンを計算するための係数. !> * [w_LaplaInv_w] !> (@ref w_module_deriv#w_laplainv_w) !> : スペクトルデータにラプラシアンの逆変換を作用させる !> * [w_DLon_w] !> (@ref w_module_deriv#w_dlon_w) !> : スペクトルデータに経度微分\f$∂/∂λ\f$を作用させる !> * [xy_GradLon_w] !> (@ref w_module_deriv#xy_gradlon_w) !> : スペクトルデータに勾配型経度微分 !> \f$1/\cosφ・∂/∂λ\f$を作用させる !> * [xy_GradLat_w] !> (@ref w_module_deriv#xy_gradlat_w) !> : スペクトルデータに勾配型緯度微分\f$∂/∂φ\f$を作用させる !> * [w_DivLon_xy] !> (@ref w_module_deriv#w_divlon_xy) !> : 格子データに発散型経度微分 !> \f$1/\cosφ・∂/∂λ\f$を作用させる !> * [w_DivLat_xy] !> (@ref w_module_deriv#w_divlat_xy) !> : 格子データに発散型緯度微分 !> \f$1/\cosφ・∂(g \cosφ)/∂φ\f$を作用させる !> * [w_Div_xy_xy] !> (@ref w_module_deriv#w_div_xy_xy) !> : ベクトル成分である 2 つの格子データに発散を作用させる !> * [w_Jacobian_w_w] !> (@ref w_module_deriv#w_jacobian_w_w) !> : 2 つのスペクトルデータからヤコビアンを計算する !> !> \f$(λ,μ)\f$ 座標 \f$(μ=\sinφ)\f$ での微分 !> !> * [xy_GradLambda_w] !> (@ref w_module_deriv#xy_gradlambda_w) !> : スペクトルデータに勾配型経度微分\f$∂/∂λ\f$を作用させる !> * [xy_GradMu_w] !> (@ref w_module_deriv#xy_gradmu_w) !> : スペクトルデータに勾配型緯度微分\f$(1-μ^2)∂/∂μ\f$を作用させる !> * [w_DivLambda_xy] !> (@ref w_module_deriv#w_divlambda_xy) !> : 格子データに発散型経度微分\f$1/(1-μ^2)・∂/∂λ\f$を作用させる !> * [w_DivMu_xy] !> (@ref w_module_deriv#w_divmu_xy) !> : 格子データに発散型緯度微分\f$∂/∂μ\f$を作用させる !> !> ### 補間 !> !> w_module_interpolate にて定義 !> !> * [Interpolate_w] !> (@ref w_module_interpolate#interpolate_w) !> : スペクトルデータから任意の点での値を求める. !> !> ### 積分・平均 !> !> w_module_integral にて定義 !> !> * [IntLonLat_xy] !> (@ref w_module_integral#intlonlat_xy) !> : 2 次元格子点データの全領域積分 !> * [AvrLonLat_xy] !> (@ref w_module_integral#avrlonlat_xy) !> : 2 次元格子点データの全領域平均 !> * [y_IntLon_xy] !> (@ref w_module_integral#y_intlon_xy) !> : 2 次元格子点データの経度方向積分 !> * [y_AvrLon_xy] !> (@ref w_module_integral#y_avrlon_xy) !> : 2 次元格子点データの経度方向平均 !> * [IntLon_x] !> (@ref w_module_integral#intlon_x) !> : 1 次元(X)格子点データの経度方向積分 !> * [AvrLon_x] !> (@ref w_module_integral#avrlon_x) !> : 1 次元(X)格子点データの経度方向平均 !> * [x_IntLat_xy] !> (@ref w_module_integral#x_intlat_xy) !> : 2 次元格子点データの緯度方向積分 !> * [x_AvrLat_xy] !> (@ref w_module_integral#x_avrlat_xy) !> : 2 次元格子点データの緯度方向平均 !> * [IntLat_y] !> (@ref w_module_integral#intlat_y) !> : 1 次元(Y)格子点データの緯度方向積分 !> * [AvrLat_y] !> (@ref w_module_integral#avrlat_y) !> : 1 次元(Y)格子点データの緯度方向平均 !> !> ### スペクトル解析 !> !> w_module_spectrum にて定義 !> !> * [nm_EnergyFromStreamfunc_w] !> (@ref w_module_spectrum#nm_energyfromstreamfunc_w) !> : 流線関数からエネルギースペクトルを計算する(水平全波数 n, 帯状波数 m 空間) !> * [n_EnergyFromStreamfunc_w] !> (@ref w_module_spectrum#n_energyfromstreamfunc_w) !> : 流線関数からエネルギースペクトルを計算する(水平全波数 n 空間) !> * [nm_EnstrophyFromStreamfunc_w] !> (@ref w_module_spectrum#nm_enstrophyfromstreamfunc_w) !> : 流線関数からエンストロフィースペクトルを計算する (水平全波数 n, 帯状波数 m 空間) !> * [n_EnstrophyFromStreamfunc_w] !> (@ref w_module_spectrum#n_enstrophyfromstreamfunc_w) !> : 流線関数からエンストロフィースペクトルを計算する (水平全波数 n 空間) !> * [w_spectrum_VMiss] !> (@ref w_module_spectrum#w_spectrum_vmiss) !> : 欠損値 !> module w_module_mint use dc_message, only : MessageNotify use w_module_base_mint, only : & & x_Lon, y_Lat, & & x_Lon_weight, y_Lat_Weight, & & xy_Lon, xy_Lat, & & xy_w, w_xy, l_nm, nm_l, & & w_StreamPotential2Vector, & & w_Vector2VorDiv, & & w_VectorCosLat2VorDiv, & & w_base_initial, & & w_base_finalize use w_module_deriv_mint, only: & & w_Lapla_w, w_LaplaInv_w, & & rn, & & w_DLon_w, & & xy_GradLon_w, xy_GradLat_w, & & w_DivLon_xy, w_DivLat_xy, & & w_Div_xy_xy, & & w_Jacobian_w_w, & & xy_GradLambda_w, xy_GradMu_w, & & w_DivLambda_xy, w_DivMu_xy, & & w_deriv_initial, & & w_deriv_Finalize use w_module_integral_mint, only: & & IntLonLat_xy, & & y_IntLon_xy, IntLon_x, & & x_IntLat_xy, IntLat_y, & & AvrLonLat_xy, & & y_AvrLon_xy, AvrLon_x, & & x_AvrLat_xy, AvrLat_y use w_module_spectrum_mint, only: & & nm_EnergyFromStreamfunc_w, & & n_EnergyFromStreamfunc_w, & & nm_EnstrophyFromStreamfunc_w, & & n_EnstrophyFromStreamfunc_w, & & w_spectrum_VMiss use w_module_interpolate_mint, only: & & Interpolate_w private public w_Initial public w_Finalize ! w_base_module public x_Lon public y_Lat public x_Lon_weight public y_Lat_Weight public xy_Lon public xy_Lat public xy_w public w_xy public l_nm public nm_l public w_StreamPotential2Vector public w_Vector2VorDiv public w_VectorCosLat2VorDiv ! w_derive_module public w_Lapla_w public w_LaplaInv_w public rn public w_DLon_w public xy_GradLon_w public xy_GradLat_w public w_DivLon_xy public w_DivLat_xy public w_Div_xy_xy public w_Jacobian_w_w public xy_GradLambda_w public xy_GradMu_w public w_DivLambda_xy public w_DivMu_xy ! w_integral_module public IntLonLat_xy public y_IntLon_xy public IntLon_x public x_IntLat_xy public IntLat_y public AvrLonLat_xy public y_AvrLon_xy public AvrLon_x public x_AvrLat_xy public AvrLat_y ! w_spectrum_module public nm_EnergyFromStreamfunc_w public n_EnergyFromStreamfunc_w public nm_EnstrophyFromStreamfunc_w public n_EnstrophyFromStreamfunc_w public w_spectrum_VMiss ! w_interpolate_module public Interpolate_w contains !--------------------------------------------------------------------- !> @brief !> スペクトル変換の格子点数, 波数および OPENMP 使用時の !> 最大スレッド数を設定する. !> !> @details !> 他の関数を呼ぶ前に, 最初にこのサブルーチンを呼んで初期設定を !> しなければならない. !> !> np_in に 1 より大きな値を指定すれば ISPACK の球面調和函数変換 !> OPENMP 並列計算ルーチンが用いられる. 並列計算を実行するには, !> 実行時に環境変数 OMP_NUM_THREADS を np_in 以下の数字に設定する等の !> システムに応じた準備が必要となる. !> !> np_in に 1 より大きな値を指定しなければ並列計算ルーチンは呼ばれない. !> subroutine w_initial(n_in,i_in,j_in,np,mint) !> 格子点数(東西) integer,intent(in) :: i_in !> 格子点数(南北) integer,intent(in) :: j_in !> 切断波数の設定 integer,intent(in) :: n_in !> OPENMP での最大スレッド数 integer,intent(in), optional :: np !> 経度方向対称性 integer,intent(in), optional :: mint if ( present (np) ) then call MessageNotify('M','w_Initial', & & 'Optional argnument np is dummy.') endif if ( present (mint) ) then call w_base_initial(n_in,i_in,j_in,mint=mint) else call w_base_initial(n_in,i_in,j_in) endif call w_deriv_initial call MessageNotify('M','w_initial',& 'w_module_mint (2020/11/18) is initialized') end subroutine w_initial !--------------------------------------------------------------------- !> @brief !> モジュールの終了処理(割り付け配列の解放)をおこなう. !> !> @warning !> 解像度を変更する際にはこのサブルーチンを呼んで終了処理を !> おこなったのちに再度 w_Initial で初期設定しなければ !> ならない. !> subroutine w_Finalize call w_base_Finalize call w_deriv_Finalize call MessageNotify('M','w_Finalize',& 'w_module_mint (2020/11/18) is finalized') end subroutine w_Finalize end module w_module_mint