Class | z_base_module |
In: |
setup/z_base_module.f90
|
z_base_module は, 1 次元 (z 方向) 等間隔交互格子を用いた有限差分 法に基づく数値モデルのための基本的な Fortran90 副プログラムおよび 関数を提供する. このモジュールは z_module の下位モジュールである.
Function : | |||
AvrZ_r : | real(DBKIND)
| ||
r_Var(kmin:kmax) : | real(DBKIND), intent(in)
|
整数格子上の配列に対し z 方向の平均操作を行う
function AvrZ_r(r_Var) ! 整数格子上の配列に対し z 方向の平均操作を行う real(DBKIND), intent(in) :: r_Var(kmin:kmax) ! 入力 real(DBKIND) :: AvrZ_r ! 出力 AvrZ_r = IntZ_r(r_Var)/sum(r_dz(1:km)) end function AvrZ_r
Function : | |||
AvrZ_z : | real(DBKIND)
| ||
z_Var(kmin:kmax) : | real(DBKIND), intent(in)
|
半整数格子上の配列に対し z 方向の平均操作を行う
function AvrZ_z(z_Var) ! 半整数格子上の配列に対し z 方向の平均操作を行う real(DBKIND), intent(in) :: z_Var(kmin:kmax) ! 入力 real(DBKIND) :: AvrZ_z ! 出力 AvrZ_z = IntZ_z(z_Var)/sum(z_dz(1:km)) end function AvrZ_z
Function : | |||
IntZ_r : | real(DBKIND)
| ||
r_Var(kmin:kmax) : | real(DBKIND), intent(in)
|
整数格子上の配列に対し z 方向に重み付きの積分を行う
function IntZ_r(r_Var) ! 整数格子上の配列に対し z 方向に重み付きの積分を行う real(DBKIND), intent(in) :: r_Var(kmin:kmax) ! 入力 real(DBKIND) :: IntZ_r ! 出力 ! 初期化 IntZ_r = 0.0d0 ! 積分 IntZ_r = sum(r_Var(1:km)*r_dz(1:km)) end function IntZ_r
Function : | |||
IntZ_z : | real(DBKIND)
| ||
z_Var(kmin:kmax) : | real(DBKIND), intent(in)
|
半整数格子上の配列に対し z 方向に重み付きの積分を行う
function IntZ_z(z_Var) ! 半整数格子上の配列に対し z 方向に重み付きの積分を行う real(DBKIND), intent(in) :: z_Var(kmin:kmax) ! 入力 real(DBKIND) :: IntZ_z ! 出力 ! 初期化 IntZ_z = 0.0d0 ! 積分 IntZ_z = sum(z_Var(1:km)*z_dz(1:km)) end function IntZ_z
Function : | |||
r_avr_z(kmin:kmax) : | real(DBKIND)
| ||
z_Var(kmin:kmax) : | real(DBKIND),intent(in)
|
平均操作を行い半整数格子点の配列値を整数格子点上へ返す
function r_avr_z(z_Var) ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す real(DBKIND),intent(in) :: z_Var(kmin:kmax) ! 入力 real(DBKIND) :: r_avr_z(kmin:kmax) ! 出力 integer :: kz ! ループ添字 ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! r_avr_z = 0.0d0 ! 平均操作 ! * 平均がとれない kmax 格子上の値は計算しない. ! * 不等間隔格子なので重みをつけて平均する. ! 点 P, Q を a:b に内分する点の X の値は, a*X(Q) + b*X(P) ! (a+b=1 を仮定) であることから, ! ! r_Var(kz) = [0.5*z_dz(kz) /r_dz(kz)]*z_Var(kz+1) + ! [0.5*z_dz(kz+1)/r_dz(kz)]*z_Var(kz) ! ! 変形すると ! ! r_Var(kz) = [z_dz(kz)*z_Var(kz+1) + z_dz(kz+1)*z_Var(kz)] ! *0.5/r_dz(kz) ! do kz = kmin, kmax-1 r_avr_z(kz) = (z_dz(kz)*z_Var(kz+1) + z_dz(kz+1)*z_Var(kz))*0.5d0/r_dz(kz) end do ! kmax 格子の値 r_avr_z(kmax) = r_avr_z(kmax-1) end function r_avr_z
Function : | |||
z_avr_r(kmin:kmax) : | real(DBKIND)
| ||
r_Var(kmin:kmax) : | real(DBKIND),intent(in)
|
平均操作を行い整数格子点の配列値を半整数格子点上へ返す
function z_avr_r(r_Var) ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す real(DBKIND),intent(in) :: r_Var(kmin:kmax) ! 入力 real(DBKIND) :: z_avr_r(kmin:kmax) ! 出力 integer :: kz ! ループ添字 ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! z_avr_r = 0.0d0 ! 平均操作 ! * 平均がとれない kmin 格子上の値は計算しない. ! do kz = kmin+1, kmax z_avr_r(kz) = (r_Var(kz) + r_Var(kz-1))*0.5d0 end do ! kmin 格子の値 z_avr_r(kmin) = z_avr_r(kmin+1) end function z_avr_r
Subroutine : | |||
k : | integer,intent(in)
| ||
zmg : | integer,intent(in)
| ||
zmin : | real(DBKIND),intent(in)
| ||
zmax : | real(DBKIND),intent(in)
|
z 方向の座標値と格子点間隔を設定する
subroutine z_axis_init(k, zmg, zmin, zmax) ! z 方向の座標値と格子点間隔を設定する integer,intent(in) :: k ! z 方向格子点数 integer,intent(in) :: zmg ! z 方向糊代格子点数 real(DBKIND),intent(in) :: zmin ! z 座標最小値 real(DBKIND),intent(in) :: zmax ! z 座標最大値 integer :: kz ! ループ添字 real(DBKIND) :: dz km = k zmargin = zmg kmin = 1 - zmargin kmax = km + zmargin allocate(z_Z(kmin:kmax)) allocate(r_Z(kmin:kmax)) allocate(z_dz(kmin:kmax)) allocate(r_dz(kmin:kmax)) dz = (zmax - zmin)/km do kz = kmin, kmax r_Z(kz) = dz * kz z_Z(kz) = dz * (kz - 0.5) r_dz(kz) = dz z_dz(kz) = dz end do end subroutine z_axis_init