!--------------------------------------------------------------------- ! Copyright (C) GFD Dennou Club, 2006. All rights reserved. !--------------------------------------------------------------------- ! != 3 次元 (xyz 方向) 等間隔交互格子 格子点設定モジュール ! !* 履歴 ! * 2007/07/15 (小高正嗣) : 3D deepconv へ移植, dc_types を Use. ! * 2006/12/26 (小高正嗣) : 平均関数の名前を元に戻す ! * 2006/06/05 (小高正嗣) : 新規作成 ! module xyz_base_module != 3 次元 (xyz 方向) 等間隔交互格子 格子点設定モジュール ! !== 概要 ! ! xyz_base_module は, 3 次元 (xyz 方向) 等間隔交互格子を用いた有限差分 ! 法に基づく数値モデルのための, 基本的な Fortran90 副プログラムおよび ! 関数を提供する. ! ! このモジュールは xyz_module の下位モジュールである. 下請けモジュール ! として data_type, x_base_module, y_base_module z_base_module ! モジュールを用いている. ! ! !== 備考 ! ! * 例外処理がない ! * その場合のメッセージ表示がない ! * 初期化の値はマシンイプシロン値としておくべき ! use dc_types, only : DBKIND => DP use x_base_module, only : im, imin, imax, xmargin, x_axis_init, & & x_X, p_X, x_dx, p_dx, x_avr_p, p_avr_x, & & IntX_p, IntX_x, AvrX_p, AvrX_x use y_base_module, only : jm, jmin, jmax, ymargin, y_axis_init, & & y_Y, q_Y, y_dy, q_dy, y_avr_q, q_avr_y, & & IntY_q, IntY_y, AvrY_q, AvrY_y use z_base_module, only : km, kmin, kmax, zmargin, z_axis_init, & & z_Z, r_Z, z_dz, r_dz, z_avr_r, r_avr_z, & & IntZ_r, IntZ_z, AvrZ_r, AvrZ_z use xy_base_module, only : y_IntX_xy, y_IntX_py, q_IntX_xq implicit none private public :: im, imin, imax, xmargin, x_X, p_X, x_dx, p_dx ! Cascaded public :: jm, jmin, jmax, ymargin, y_Y, q_Y, y_dy, q_dy ! Cascaded public :: km, kmin, kmax, zmargin, z_Z, r_Z, z_dz, r_dz ! Cascaded public :: x_avr_p, p_avr_x, IntX_p, IntX_x, AvrX_p, AvrX_x ! Cascaded public :: y_avr_q, q_avr_y, IntY_q, IntY_y, AvrY_q, AvrY_y ! Cascaded public :: z_avr_r, r_avr_z, IntZ_r, IntZ_z, AvrZ_r, AvrZ_z ! Cascaded public :: xyz_X, xyz_Y, xyz_Z, xyz_dX, xyz_dY, xyz_dZ public :: xyz_axis_init !! public :: xyz_avr_pyz, xyr_avr_pyr, xqz_avr_pqz public :: pyz_avr_xyz, pyr_avr_xyr, pqz_avr_xqz public :: xyz_avr_xqz, pyz_avr_pqz, xyr_avr_xqr public :: xqz_avr_xyz, pqz_avr_pyz, xqr_avr_xyr public :: xyz_avr_xyr, pyz_avr_pyr, xqz_avr_xqr public :: xyr_avr_xyz, pyr_avr_pyz, xqr_avr_xqz !! public :: pqz_avr_xyz, pyr_avr_xyz, xqr_avr_xyz public :: xyz_avr_pqz, xyz_avr_pyr, xyz_avr_xqr !! public :: yz_IntX_xyz, yz_IntX_pyz, qz_IntX_xqz, yr_IntX_xyr public :: xz_IntY_xyz, xz_IntY_xqz, pz_IntY_pyz, xr_IntY_xyr public :: xy_IntZ_xyz, xy_IntZ_xyr, py_IntZ_pyz, xq_IntZ_xqz public :: z_IntXY_xyz, z_IntXY_pyz, z_IntXY_xqz, r_IntXY_xyr public :: IntXYZ_xyz , IntXYZ_pyz , IntXYZ_xqz , IntXYZ_xyr public :: yz_AvrX_xyz, yz_AvrX_pyz, qz_AvrX_xqz, yr_AvrX_xyr public :: xz_AvrY_xyz, xz_AvrY_xqz, pz_AvrY_pyz, xr_AvrY_xyr public :: xy_AvrZ_xyz, xy_AvrZ_xyr, py_AvrZ_pyz, xq_AvrZ_xqz public :: z_AvrXY_xyz, z_AvrXY_pyz, z_AvrXY_xqz, r_AvrXY_xyr public :: AvrXYZ_xyz , AvrXYZ_pyz , AvrXYZ_xqz , AvrXYZ_xyr real(DBKIND),allocatable :: xyz_X(:,:,:) ! x 座標(半整数格子) real(DBKIND),allocatable :: xyz_Y(:,:,:) ! y 座標(半整数格子) real(DBKIND),allocatable :: xyz_Z(:,:,:) ! z 座標(半整数格子) real(DBKIND),allocatable :: xyz_dX(:,:,:) ! x 格子間隔(半整数格子) real(DBKIND),allocatable :: xyz_dY(:,:,:) ! y 格子間隔(半整数格子) real(DBKIND),allocatable :: xyz_dZ(:,:,:) ! z 格子間隔(半整数格子) save xyz_X, xyz_Y, xyz_Z interface xyz_avr_pyz module procedure xaa_avr_paa end interface interface xyr_avr_pyr module procedure xaa_avr_paa end interface interface xqz_avr_pqz module procedure xaa_avr_paa end interface interface pyz_avr_xyz module procedure paa_avr_xaa end interface interface pqz_avr_xqz module procedure paa_avr_xaa end interface interface pyr_avr_xyr module procedure paa_avr_xaa end interface interface xyz_avr_xqz module procedure aya_avr_aqa end interface interface pyz_avr_pqz module procedure aya_avr_aqa end interface interface xyr_avr_xqr module procedure aya_avr_aqa end interface interface xqz_avr_xyz module procedure aqa_avr_aya end interface interface pqz_avr_pyz module procedure aqa_avr_aya end interface interface xqr_avr_xyr module procedure aqa_avr_aya end interface interface xyz_avr_xyr module procedure aaz_avr_aar end interface interface pyz_avr_pyr module procedure aaz_avr_aar end interface interface xqz_avr_xqr module procedure aaz_avr_aar end interface interface xyr_avr_xyz module procedure aar_avr_aaz end interface interface pyr_avr_pyz module procedure aar_avr_aaz end interface interface xqr_avr_xqz module procedure aar_avr_aaz end interface interface yz_IntX_xyz module procedure aa_IntX_xaa end interface interface qz_IntX_xqz module procedure aa_IntX_xaa end interface interface yr_IntX_xyr module procedure aa_IntX_xaa end interface interface xz_IntY_xyz module procedure aa_IntY_aya end interface interface pz_IntY_pyz module procedure aa_IntY_aya end interface interface xr_IntY_xyr module procedure aa_IntY_aya end interface interface xy_IntZ_xyz module procedure aa_IntZ_aaz end interface interface py_IntZ_pyz module procedure aa_IntZ_aaz end interface interface xq_IntZ_xqz module procedure aa_IntZ_aaz end interface interface z_IntXY_xyz module procedure a_IntXY_xya end interface interface r_IntXY_xyr module procedure a_IntXY_xya end interface interface yz_AvrX_xyz module procedure aa_AvrX_xaa end interface interface qz_AvrX_xqz module procedure aa_AvrX_xaa end interface interface yr_AvrX_xyr module procedure aa_AvrX_xaa end interface interface xz_AvrY_xyz module procedure aa_AvrY_aya end interface interface pz_AvrY_pyz module procedure aa_AvrY_aya end interface interface xr_AvrY_xyr module procedure aa_AvrY_aya end interface interface xy_AvrZ_xyz module procedure aa_AvrZ_aaz end interface interface py_AvrZ_pyz module procedure aa_AvrZ_aaz end interface interface xq_AvrZ_xqz module procedure aa_AvrZ_aaz end interface interface z_AvrXY_xyz module procedure a_AvrXY_xya end interface interface r_AvrXY_xyr module procedure a_AvrXY_xya end interface contains !-------------------------------------------------------------------- subroutine xyz_axis_init(i, j, k, xmg, ymg, zmg, & & xmin, xmax, ymin, ymax, zmin, zmax) integer ,intent(in) :: i ! x 方向格子点数 integer ,intent(in) :: j ! y 方向格子点数 integer ,intent(in) :: k ! z 方向格子点数 integer ,intent(in) :: xmg ! x 方向糊代格子点数 integer ,intent(in) :: ymg ! y 方向糊代格子点数 integer ,intent(in) :: zmg ! z 方向糊代格子点数 real(DBKIND),intent(in) :: xmin ! x 座標最小値 real(DBKIND),intent(in) :: xmax ! x 座標最大値 real(DBKIND),intent(in) :: ymin ! y 座標最小値 real(DBKIND),intent(in) :: ymax ! y 座標最大値 real(DBKIND),intent(in) :: zmin ! z 座標最小値 real(DBKIND),intent(in) :: zmax ! z 座標最大値 real(DBKIND),allocatable :: xy_X(:,:)! x 座標(半整数格子, 作業配列) real(DBKIND),allocatable :: xy_Y(:,:)! y 座標(半整数格子, 作業配列) real(DBKIND),allocatable :: yz_Z(:,:)! z 座標(半整数格子, 作業配列) ! 配列の上下限の値, 座標値と格子点間隔を設定 ! * 1 次元用のサブルーチンを用いる ! call x_axis_init(i, xmg, xmin, xmax) call y_axis_init(j, ymg, ymin, ymax) call z_axis_init(k, zmg, zmin, zmax) ! 3 次元座標配列の設定 ! * 組み込み関数 spread を用いる. ! * 中間配列として 2 次元座標配列を作り, それを 3 次元に拡張する. ! allocate(xy_X(imin:imax,jmin:jmax)) allocate(xy_Y(imin:imax,jmin:jmax)) allocate(yz_Z(jmin:jmax,kmin:kmax)) allocate(xyz_X(imin:imax,jmin:jmax,kmin:kmax)) allocate(xyz_Y(imin:imax,jmin:jmax,kmin:kmax)) allocate(xyz_Z(imin:imax,jmin:jmax,kmin:kmax)) allocate(xyz_dX(imin:imax,jmin:jmax,kmin:kmax)) allocate(xyz_dY(imin:imax,jmin:jmax,kmin:kmax)) allocate(xyz_dZ(imin:imax,jmin:jmax,kmin:kmax)) xy_X = spread(x_X, 2,size(y_Y)) xyz_X = spread(xy_X,3,size(z_Z)) xy_X = spread(x_dX, 2,size(y_dY)) xyz_dX = spread(xy_X,3,size(z_dZ)) xy_Y = spread(y_Y, 1,size(x_X)) xyz_Y = spread(xy_Y,3,size(z_Z)) xy_Y = spread(y_dY, 1,size(x_dX)) xyz_dY = spread(xy_Y,3,size(z_dZ)) yz_Z = spread(z_Z, 1,size(y_Y)) xyz_Z = spread(yz_Z,1,size(x_X)) yz_Z = spread(z_dZ, 1,size(y_dY)) xyz_dZ = spread(yz_Z,1,size(x_dX)) deallocate(xy_X) deallocate(xy_Y) deallocate(yz_Z) end subroutine xyz_axis_init !-------------------------------------------------------------------- function xaa_avr_paa(paa_Var) ! 平均操作を行い x 方向整数格子点の配列値を半整数点上へ返す real(DBKIND),intent(in) :: paa_Var(imin:imax,jmin:jmax,kmin:kmax) real(DBKIND) :: xaa_avr_paa(imin:imax,jmin:jmax,kmin:kmax) ! integer :: jy, kz integer :: ix ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! xaa_avr_paa = 0.0d0 ! 平均操作 ! * 関数 x_avr_p を用いて計算 ! ! do kz = kmin, kmax ! do jy = jmin, jmax ! xaa_avr_paa(:,jy,kz) = x_avr_p(paa_Var(:,jy,kz)) ! end do ! end do do ix = imin+1, imax xaa_avr_paa(ix,:,:) = (paa_Var(ix,:,:) + paa_Var(ix-1,:,:))*0.5d0 end do ! imin 格子上の値 xaa_avr_paa(imin,:,:) = xaa_avr_paa(imin+1,:,:) end function xaa_avr_paa !-------------------------------------------------------------------- function paa_avr_xaa(xaa_Var) ! 半整数格子点の配列値を整数点上へ返す real(DBKIND),intent(in) :: xaa_Var(imin:imax,jmin:jmax,kmin:kmax) real(DBKIND) :: paa_avr_xaa(imin:imax,jmin:jmax,kmin:kmax) ! integer :: jy, kz integer :: ix ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! paa_avr_xaa = 0.0d0 ! 平均操作 ! * 関数 p_avr_x を用いて計算 ! ! do kz = kmin, kmax ! do jy = jmin, jmax ! paa_avr_xaa(:,jy,kz) = p_avr_x(xaa_Var(:,jy,kz)) ! end do ! end do do ix = imin, imax-1 paa_avr_xaa(ix,:,:) = & & (x_dx(ix)*xaa_Var(ix+1,:,:) + x_dx(ix+1)*xaa_Var(ix,:,:)) & & *0.5d0/p_dx(ix) end do paa_avr_xaa(imax,:,:) = paa_avr_xaa(imax-1,:,:) end function paa_avr_xaa !-------------------------------------------------------------------- function aya_avr_aqa(aqa_Var) ! 平均操作を行い y 方向半整数格子点の配列値を整数格子点上へ返す real(DBKIND),intent(in) :: aqa_Var(imin:imax,jmin:jmax,kmin:kmax) real(DBKIND) :: aya_avr_aqa(imin:imax,jmin:jmax,kmin:kmax) ! integer :: ix, kz integer :: jy ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! aya_avr_aqa = 0.0d0 ! 平均操作 ! * 関数 y_q を用いて計算 ! ! do kz = kmin, kmax ! do ix = imin, imax ! aya_avr_aqa(ix,:,kz) = y_avr_q(aqa_Var(ix,:,kz)) ! end do ! end do do jy = jmin+1, jmax aya_avr_aqa(:,jy,:) = (aqa_Var(:,jy,:) + aqa_Var(:,jy-1,:))*0.5d0 end do aya_avr_aqa(:,jmin,:) = aya_avr_aqa(:,jmin+1,:) end function aya_avr_aqa !-------------------------------------------------------------------- function aqa_avr_aya(aya_Var) ! 平均操作を行い y 方向半整数格子点の配列値を整数格子点上へ返す real(DBKIND),intent(in) :: aya_Var(imin:imax,jmin:jmax,kmin:kmax) real(DBKIND) :: aqa_avr_aya(imin:imax,jmin:jmax,kmin:kmax) ! integer :: ix, kz integer :: jy ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! aqa_avr_aya = 0.0d0 ! 平均操作 ! * 関数 q_avr_y を用いて計算 ! ! do kz = kmin, kmax ! do ix = imin, imax ! aqa_avr_aya(ix,:,kz) = q_avr_y(aya_Var(ix,:,kz)) ! end do ! end do do jy = jmin, jmax-1 aqa_avr_aya(:,jy,:) = & & (y_dy(jy)*aya_Var(:,jy+1,:) + y_dy(jy+1)*aya_Var(:,jy,:)) & & * 0.5d0/q_dy(jy) end do aqa_avr_aya(:,jmax,:) = aqa_avr_aya(:,jmax-1,:) end function aqa_avr_aya !-------------------------------------------------------------------- function aaz_avr_aar(aar_Var) ! 平均操作を行い z 方向整数格子点の配列値を半整数格子点上へ返す real(DBKIND),intent(in) :: aar_Var(imin:imax,jmin:jmax,kmin:kmax) real(DBKIND) :: aaz_avr_aar(imin:imax,jmin:jmax,kmin:kmax) ! integer :: ix, jy integer :: kz ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! aaz_avr_aar = 0.0d0 ! 平均操作 ! * 関数 z_avr_r を用いて計算 ! ! do jy = jmin, jmax ! do ix = imin, imax ! aaz_avr_aar(ix,jy,:) = z_avr_r(aar_Var(ix,jy,:)) ! end do ! end do do kz = kmin+1, kmax aaz_avr_aar(:,:,kz) = (aar_Var(:,:,kz) + aar_Var(:,:,kz-1))*0.5d0 end do aaz_avr_aar(:,:,kmin) = aaz_avr_aar(:,:,kmin+1) end function aaz_avr_aar !-------------------------------------------------------------------- function aar_avr_aaz(aaz_Var) ! 平均操作を行い z 方向半整数格子点の配列値を整数格子点上へ返す real(DBKIND),intent(in) :: aaz_Var(imin:imax,jmin:jmax,kmin:kmax) real(DBKIND) :: aar_avr_aaz(imin:imax,jmin:jmax,kmin:kmax) ! integer :: ix, jy integer :: kz ! 初期化 ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき ! aar_avr_aaz = 0.0d0 ! 平均操作 ! * 関数 r_avr_z を用いて計算 ! ! do jy = jmin, jmax ! do ix = imin, imax ! aar_avr_aaz(ix,jy,:) = r_avr_z(aaz_Var(ix,jy,:)) ! end do ! end do do kz = kmin, kmax-1 aar_avr_aaz(:,:,kz) = & & (z_dz(kz)*aaz_Var(:,:,kz+1) + z_dz(kz+1)*aaz_Var(:,:,kz)) & & *0.5d0/r_dz(kz) end do aar_avr_aaz(:,:,kmax) = aar_avr_aaz(:,:,kmax-1) end function aar_avr_aaz !-------------------------------------------------------------------- function pqz_avr_xyz(xyz_Var) real(DBKIND), intent(in) :: xyz_Var(imin:imax,jmin:jmax,kmin:kmax)! 入力 real(DBKIND) :: pqz_avr_xyz(imin:imax,jmin:jmax,kmin:kmax) ! 出力 integer :: ix, jy ! ループ添字 ! 初期化 pqz_avr_xyz = 0.0d0 do jy = jmin, jmax-1 do ix = imin, imax-1 pqz_avr_xyz(ix,jy,:) = & & ( xyz_Var(ix,jy,:) + xyz_Var(ix+1,jy,:) + & & xyz_Var(ix,jy+1,:) + xyz_Var(ix+1,jy+1,:) ) * 0.25d0 end do end do end function pqz_avr_xyz !-------------------------------------------------------------------- function pyr_avr_xyz(xyz_Var) real(DBKIND), intent(in) :: xyz_Var(imin:imax,jmin:jmax,kmin:kmax)! 入力 real(DBKIND) :: pyr_avr_xyz(imin:imax,jmin:jmax,kmin:kmax) ! 出力 integer :: ix, kz ! ループ添字 ! 初期化 pyr_avr_xyz = 0.0d0 do kz = kmin, kmax-1 do ix = imin, imax-1 pyr_avr_xyz(ix,:,kz) = & & ( xyz_Var(ix,:,kz) + xyz_Var(ix+1,:,kz) + & & xyz_Var(ix,:,kz+1) + xyz_Var(ix+1,:,kz+1) ) * 0.25d0 end do end do end function pyr_avr_xyz !-------------------------------------------------------------------- function xqr_avr_xyz(xyz_Var) real(DBKIND), intent(in) :: xyz_Var(imin:imax,jmin:jmax,kmin:kmax)! 入力 real(DBKIND) :: xqr_avr_xyz(imin:imax,jmin:jmax,kmin:kmax) ! 出力 integer :: jy, kz ! ループ添字 ! 初期化 xqr_avr_xyz = 0.0d0 do kz = kmin, kmax-1 do jy = jmin, jmax-1 xqr_avr_xyz(:,jy,kz) = & & ( xyz_Var(:,jy,kz) + xyz_Var(:,jy+1,kz) + & & xyz_Var(:,jy,kz+1) + xyz_Var(:,jy+1,kz+1) ) * 0.25d0 end do end do end function xqr_avr_xyz !-------------------------------------------------------------------- function xyz_avr_pqz(pqz_Var) real(DBKIND), intent(in) :: pqz_Var(imin:imax,jmin:jmax,kmin:kmax)! 入力 real(DBKIND) :: xyz_avr_pqz(imin:imax,jmin:jmax,kmin:kmax) ! 出力 integer :: ix, jy ! ループ添字 ! 初期化 xyz_avr_pqz = 0.0d0 do jy = jmin+1, jmax do ix = imin+1, imax xyz_avr_pqz(ix,jy,:) = & & ( pqz_Var(ix-1,jy-1,:) + pqz_Var(ix,jy-1,:) + & & pqz_Var(ix-1,jy,:) + pqz_Var(ix,jy,:) ) * 0.25d0 end do end do end function xyz_avr_pqz !-------------------------------------------------------------------- function xyz_avr_pyr(pyr_Var) real(DBKIND), intent(in) :: pyr_Var(imin:imax,jmin:jmax,kmin:kmax)! 入力 real(DBKIND) :: xyz_avr_pyr(imin:imax,jmin:jmax,kmin:kmax) ! 出力 integer :: ix, kz ! ループ添字 ! 初期化 xyz_avr_pyr = 0.0d0 do kz = kmin+1, kmax do ix = imin+1, imax xyz_avr_pyr(ix,:,kz) = & & ( pyr_Var(ix-1,:,kz-1) + pyr_Var(ix,:,kz-1) + & & pyr_Var(ix-1,:,kz) + pyr_Var(ix,:,kz) ) * 0.25d0 end do end do end function xyz_avr_pyr !-------------------------------------------------------------------- function xyz_avr_xqr(xqr_Var) real(DBKIND), intent(in) :: xqr_Var(imin:imax,jmin:jmax,kmin:kmax)! 入力 real(DBKIND) :: xyz_avr_xqr(imin:imax,jmin:jmax,kmin:kmax) ! 出力 integer :: jy, kz ! ループ添字 ! 初期化 xyz_avr_xqr = 0.0d0 do kz = kmin+1, kmax do jy = jmin+1, jmax xyz_avr_xqr(:,jy,kz) = & & ( xqr_Var(:,jy-1,kz-1) + xqr_Var(:,jy,kz-1) + & & xqr_Var(:,jy-1,kz) + xqr_Var(:,jy,kz) ) * 0.25d0 end do end do end function xyz_avr_xqr !-------------------------------------------------------------------- function aa_IntX_xaa(xaa_Var) ! xaa 格子上の配列に対し x 方向積分を行う real(DBKIND), intent(in) :: xaa_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: aa_IntX_xaa(jmin:jmax,kmin:kmax) ! 出力 integer :: jy, kz ! ループ添字 ! 初期化 aa_IntX_xaa = 0.0d0 ! 積分 do kz = kmin, kmax do jy = jmin, jmax aa_IntX_xaa(jy,kz) = IntX_x(xaa_Var(:,jy,kz)) end do end do end function aa_IntX_xaa !-------------------------------------------------------------------- function yz_IntX_pyz(pyz_Var) ! pyz 格子上の配列に対し x 方向積分を行う real(DBKIND), intent(in) :: pyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: yz_IntX_pyz(jmin:jmax,kmin:kmax) ! 出力 integer :: jy, kz ! ループ添字 ! 初期化 yz_IntX_pyz = 0.0d0 ! 積分 do kz = kmin, kmax do jy = jmin, jmax yz_IntX_pyz(jy,kz) = IntX_p(pyz_Var(:,jy,kz)) end do end do end function yz_IntX_pyz !-------------------------------------------------------------------- function aa_IntY_aya(aya_Var) ! aya 格子上の配列に対し y 方向積分を行う real(DBKIND), intent(in) :: aya_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: aa_IntY_aya(imin:imax,kmin:kmax) ! 出力 integer :: ix, kz ! ループ添字 ! 初期化 aa_IntY_aya = 0.0d0 ! 積分 do kz = kmin, kmax do ix = imin, imax aa_IntY_aya(ix,kz) = IntY_y(aya_Var(ix,:,kz)) end do end do end function aa_IntY_aya !-------------------------------------------------------------------- function xz_IntY_xqz(xqz_Var) ! xqz 格子上の配列に対し y 方向積分を行う real(DBKIND), intent(in) :: xqz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: xz_IntY_xqz(imin:imax,kmin:kmax) ! 出力 integer :: ix, kz ! ループ添字 ! 初期化 xz_IntY_xqz = 0.0d0 ! 積分 do kz = kmin, kmax do ix = imin, imax xz_IntY_xqz(ix,kz) = IntY_q(xqz_Var(ix,:,kz)) end do end do end function xz_IntY_xqz !-------------------------------------------------------------------- function aa_IntZ_aaz(aaz_Var) ! aaz 格子上の配列に対し z 方向積分を行う real(DBKIND), intent(in) :: aaz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: aa_IntZ_aaz(imin:imax,jmin:jmax) ! 出力 integer :: ix, jy ! ループ添字 ! 初期化 aa_IntZ_aaz = 0.0d0 ! 積分 do jy = jmin, jmax do ix = imin, imax aa_IntZ_aaz(ix,jy) = IntZ_z(aaz_Var(ix,jy,:)) end do end do end function aa_IntZ_aaz !-------------------------------------------------------------------- function xy_IntZ_xyr(xyr_Var) ! xyr 格子上の配列に対し z 方向積分を行う real(DBKIND), intent(in) :: xyr_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: xy_IntZ_xyr(imin:imax,jmin:jmax) ! 出力 integer :: ix, jy ! ループ添字 ! 初期化 xy_IntZ_xyr = 0.0d0 ! 積分 do jy = jmin, jmax do ix = imin, imax xy_IntZ_xyr(ix,jy) = IntZ_r(xyr_Var(ix,jy,:)) end do end do end function xy_IntZ_xyr !-------------------------------------------------------------------- function a_IntXY_xya(xya_Var) ! xya 格子上の配列に対し xy 方向積分を行う real(DBKIND), intent(in) :: xya_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: a_IntXY_xya(kmin:kmax) ! 出力 integer :: kz ! ループ添字 ! 初期化 a_IntXY_xya = 0.0d0 ! 積分 do kz = kmin, kmax a_IntXY_xya(kz) = IntY_y(y_IntX_xy(xya_Var(:,:,kz))) end do end function a_IntXY_xya !-------------------------------------------------------------------- function z_IntXY_pyz(pyz_Var) ! pyz 格子上の配列に対し xy 方向積分を行う real(DBKIND), intent(in) :: pyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: z_IntXY_pyz(kmin:kmax) ! 出力 integer :: kz ! ループ添字 ! 初期化 z_IntXY_pyz = 0.0d0 ! 積分 do kz = kmin, kmax z_IntXY_pyz(kz) = IntY_y(y_IntX_py(pyz_Var(:,:,kz))) end do end function z_IntXY_pyz !-------------------------------------------------------------------- function z_IntXY_xqz(xqz_Var) ! xqz 格子上の配列に対し xy 方向積分を行う real(DBKIND), intent(in) :: xqz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: z_IntXY_xqz(kmin:kmax) ! 出力 integer :: kz ! ループ添字 ! 初期化 z_IntXY_xqz = 0.0d0 ! 積分 do kz = kmin, kmax z_IntXY_xqz(kz) = IntY_q(q_IntX_xq(xqz_Var(:,:,kz))) end do end function z_IntXY_xqz !-------------------------------------------------------------------- function IntXYZ_xyz(xyz_Var) ! xyz 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: xyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: IntXYZ_xyz ! 出力 ! 初期化 IntXYZ_xyz = 0.0d0 IntXYZ_xyz = IntZ_z(a_IntXY_xya(xyz_Var)) end function IntXYZ_xyz !-------------------------------------------------------------------- function IntXYZ_pyz(pyz_Var) ! pyz 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: pyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: IntXYZ_pyz ! 出力 ! 初期化 IntXYZ_pyz = 0.0d0 IntXYZ_pyz = IntZ_z(z_IntXY_pyz(pyz_Var)) end function IntXYZ_pyz !-------------------------------------------------------------------- function IntXYZ_xqz(xqz_Var) ! xqz 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: xqz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: IntXYZ_xqz ! 出力 ! 初期化 IntXYZ_xqz = 0.0d0 IntXYZ_xqz = IntZ_z(z_IntXY_xqz(xqz_Var)) end function IntXYZ_xqz !-------------------------------------------------------------------- function IntXYZ_xyr(xyr_Var) ! xyr 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: xyr_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: IntXYZ_xyr ! 出力 ! 初期化 IntXYZ_xyr = 0.0d0 IntXYZ_xyr = IntZ_r(a_IntXY_xya(xyr_Var)) end function IntXYZ_xyr !-------------------------------------------------------------------- function aa_AvrX_xaa(xaa_Var) ! xaa 格子上の配列に対し x 方向平均を行う real(DBKIND), intent(in) :: xaa_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: aa_AvrX_xaa(jmin:jmax,kmin:kmax) ! 出力 ! 初期化 aa_AvrX_xaa = 0.0d0 ! 平均 aa_AvrX_xaa = aa_IntX_xaa(xaa_Var)/sum(x_dx(1:im)) end function aa_AvrX_xaa !-------------------------------------------------------------------- function yz_AvrX_pyz(pyz_Var) ! pyz 格子上の配列に対し x 方向平均を行う real(DBKIND), intent(in) :: pyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: yz_AvrX_pyz(jmin:jmax,kmin:kmax) ! 出力 ! 初期化 yz_AvrX_pyz = 0.0d0 ! 平均 yz_AvrX_pyz = yz_IntX_pyz(pyz_Var)/sum(x_dx(1:im)) end function yz_AvrX_pyz !-------------------------------------------------------------------- function aa_AvrY_aya(aya_Var) ! aya 格子上の配列に対し y 方向平均を行う real(DBKIND), intent(in) :: aya_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: aa_AvrY_aya(imin:imax,kmin:kmax) ! 出力 ! 初期化 aa_AvrY_aya = 0.0d0 ! 平均 aa_AvrY_aya = aa_IntY_aya(aya_Var)/sum(y_dy(1:jm)) end function aa_AvrY_aya !-------------------------------------------------------------------- function xz_AvrY_xqz(xqz_Var) ! xqz 格子上の配列に対し y 方向平均を行う real(DBKIND), intent(in) :: xqz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: xz_AvrY_xqz(imin:imax,kmin:kmax) ! 出力 ! 初期化 xz_AvrY_xqz = 0.0d0 ! 平均 xz_AvrY_xqz = xz_IntY_xqz(xqz_Var)/sum(y_dy(1:jm)) end function xz_AvrY_xqz !-------------------------------------------------------------------- function aa_AvrZ_aaz(aaz_Var) ! aaz 格子上の配列に対し z 方向平均を行う real(DBKIND), intent(in) :: aaz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: aa_AvrZ_aaz(imin:imax,jmin:jmax) ! 出力 ! 初期化 aa_AvrZ_aaz = 0.0d0 ! 平均 aa_AvrZ_aaz = aa_IntZ_aaz(aaz_Var)/sum(z_dz(1:km)) end function aa_AvrZ_aaz !-------------------------------------------------------------------- function xy_AvrZ_xyr(xyr_Var) ! xyr 格子上の配列に対し z 方向平均を行う real(DBKIND), intent(in) :: xyr_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: xy_AvrZ_xyr(imin:imax,jmin:jmax) ! 出力 ! 初期化 xy_AvrZ_xyr = 0.0d0 ! 平均 xy_AvrZ_xyr = xy_IntZ_xyr(xyr_Var)/sum(z_dz(1:km)) end function xy_AvrZ_xyr !-------------------------------------------------------------------- function AvrXYZ_xyz(xyz_Var) ! xyz 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: xyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: AvrXYZ_xyz ! 出力 ! 初期化 AvrXYZ_xyz = 0.0d0 AvrXYZ_xyz = IntXYZ_xyz(xyz_Var)/ & & (sum(x_dx(1:im))*sum(y_dy(1:jm))*sum(z_dz(1:km))) end function AvrXYZ_xyz !-------------------------------------------------------------------- function AvrXYZ_pyz(pyz_Var) ! pyz 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: pyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: AvrXYZ_pyz ! 出力 ! 初期化 AvrXYZ_pyz = 0.0d0 AvrXYZ_pyz = IntXYZ_pyz(pyz_Var)/ & & (sum(x_dx(1:im))*sum(y_dy(1:jm))*sum(z_dz(1:km))) end function AvrXYZ_pyz !-------------------------------------------------------------------- function AvrXYZ_xqz(xqz_Var) ! xqz 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: xqz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: AvrXYZ_xqz ! 出力 ! 初期化 AvrXYZ_xqz = 0.0d0 AvrXYZ_xqz = IntXYZ_xqz(xqz_Var)/ & & (sum(x_dx(1:im))*sum(y_dy(1:jm))*sum(z_dz(1:km))) end function AvrXYZ_xqz !-------------------------------------------------------------------- function AvrXYZ_xyr(xyr_Var) ! xyr 格子上の配列に対し領域積分を行う real(DBKIND), intent(in) :: xyr_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: AvrXYZ_xyr ! 出力 ! 初期化 AvrXYZ_xyr = 0.0d0 AvrXYZ_xyr = IntXYZ_xyr(xyr_Var)/ & & (sum(x_dx(1:im))*sum(y_dy(1:jm))*sum(z_dz(1:km))) end function AvrXYZ_xyr !-------------------------------------------------------------------- function a_AvrXY_xya(xya_Var) ! xya 格子上の配列に対し xy 方向平均を行う real(DBKIND), intent(in) :: xya_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: a_AvrXY_xya(kmin:kmax) ! 出力 ! 初期化 a_AvrXY_xya = 0.0d0 ! 平均 a_AvrXY_xya = a_IntXY_xya(xya_Var)/(sum(x_dx(1:im))*sum(y_dy(1:jm))) end function a_AvrXY_xya !-------------------------------------------------------------------- function z_AvrXY_pyz(pyz_Var) ! pyz 格子上の配列に対し xy 方向平均を行う real(DBKIND), intent(in) :: pyz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: z_AvrXY_pyz(kmin:kmax) ! 出力 ! 初期化 z_AvrXY_pyz = 0.0d0 ! 平均 z_AvrXY_pyz = z_IntXY_pyz(pyz_Var)/(sum(x_dx(1:im))*sum(y_dy(1:jm))) end function z_AvrXY_pyz !-------------------------------------------------------------------- function z_AvrXY_xqz(xqz_Var) ! xqz 格子上の配列に対し xy 方向積分を行う real(DBKIND), intent(in) :: xqz_Var(imin:imax,jmin:jmax,kmin:kmax) ! 入力 real(DBKIND) :: z_AvrXY_xqz(kmin:kmax) ! 出力 ! 初期化 z_AvrXY_xqz = 0.0d0 ! 平均 z_AvrXY_xqz = z_IntXY_xqz(xqz_Var)/(sum(x_dx(1:im))*sum(y_dy(1:jm))) end function z_AvrXY_xqz !-------------------------------------------------------------------- end module xyz_base_module