Class xy_base_module
In: setup/xy_base_module.f90

2 次元 (xy 方向) 等間隔交互格子 有限差分モデル用 基本モジュール

概要

xy_base_module は, 2 次元 (xy 方向) 等間隔交互格子を用いた有限差分 法に基づく数値モデルのための, 基本的な Fortran90 副プログラムおよび 関数を提供する.

このモジュールは xy_module の下位モジュールである. 下請けモジュール として data_type, x_base_module, y_base_module モジュールを用いている.

備考

  • 例外処理がない
  • その場合のメッセージ表示がない
  • 初期化の値はマシンイプシロン値としておくべき

Methods

AvrXY_py   AvrXY_xq   AvrXY_xy   AvrX_p   AvrX_x   AvrY_q   AvrY_y   IntXY_py   IntXY_xq   IntXY_xy   IntX_p   IntX_x   IntY_q   IntY_y   im   imax   imin   jm   jmax   jmin   p_AvrY_py   p_IntY_py   p_X   p_avr_x   p_dx   pq_avr_py   pq_avr_xq   py_avr_pq   py_avr_xy   py_dx   q_AvrX_xq   q_IntX_xq   q_Y   q_avr_y   q_dy   x_AvrY_xq   x_AvrY_xy   x_IntY_xq   x_IntY_xy   x_X   x_avr_p   x_dx   xmargin   xq_avr_pq   xq_avr_xy   xq_dy   xy_X   xy_Y   xy_avr_py   xy_avr_xq   xy_axis_init   xy_dx   xy_dy   y_AvrX_py   y_AvrX_xy   y_IntX_py   y_IntX_xy   y_Y   y_avr_q   y_dy   ymargin  

Included Modules

dc_types x_base_module y_base_module

Public Instance methods

Function :
AvrXY_py :real(DBKIND)
: 出力
py_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

py 格子上の配列に対し領域平均を行う

[Source]

    function AvrXY_py(py_Var)
      ! py 格子上の配列に対し領域平均を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: AvrXY_py                    ! 出力


      AvrXY_py = IntXY_py(py_Var)/(sum(p_dx(1:im))*sum(y_dy(1:jm)))

    end function AvrXY_py
Function :
AvrXY_xq :real(DBKIND)
: 出力
xq_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xq 格子上の配列に対し領域平均を行う

[Source]

    function AvrXY_xq(xq_Var)
      ! xq 格子上の配列に対し領域平均を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: AvrXY_xq                    ! 出力


      AvrXY_xq = IntXY_xq(xq_Var)/(sum(x_dx(1:im))*sum(q_dy(1:jm)))

    end function AvrXY_xq
Function :
AvrXY_xy :real(DBKIND)
: 出力
xy_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xy 格子上の配列に対し領域平均を行う

[Source]

    function AvrXY_xy(xy_Var)
      ! xy 格子上の配列に対し領域平均を行う

      real(DBKIND), intent(in) :: xy_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: AvrXY_xy                    ! 出力


      AvrXY_xy = IntXY_xy(xy_Var)/(sum(x_dx(1:im))*sum(y_dy(1:jm)))

    end function AvrXY_xy
AvrX_p( p_Var ) result(AvrX_p)
Function :
AvrX_p :real(DBKIND)
: 出力
p_Var(imin:imax) :real(DBKIND), intent(in)
: 入力

整数格子上の配列に対し x 方向の平均操作を行う

Original external subprogram is x_base_module#AvrX_p

AvrX_x( x_Var ) result(AvrX_x)
Function :
AvrX_x :real(DBKIND)
: 出力
x_Var(imin:imax) :real(DBKIND), intent(in)
: 入力

半整数格子上の配列に対し x 方向の平均操作を行う

Original external subprogram is x_base_module#AvrX_x

AvrY_q( q_Var ) result(AvrY_q)
Function :
AvrY_q :real(DBKIND)
: 出力
q_Var(jmin:jmax) :real(DBKIND), intent(in)
: 入力

整数格子上の配列に対し y 方向の平均操作を行う

Original external subprogram is y_base_module#AvrY_q

AvrY_y( y_Var ) result(AvrY_y)
Function :
AvrY_y :real(DBKIND)
: 出力
y_Var(jmin:jmax) :real(DBKIND), intent(in)
: 入力

半整数格子上の配列に対し y 方向の平均操作を行う

Original external subprogram is y_base_module#AvrY_y

Function :
IntXY_py :real(DBKIND)
: 出力
py_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

py 格子上の配列に対し領域積分を行う

[Source]

    function IntXY_py(py_Var)
      ! py 格子上の配列に対し領域積分を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: IntXY_py                    ! 出力


      IntXY_py = IntY_y(y_IntX_py(py_Var))

    end function IntXY_py
Function :
IntXY_xq :real(DBKIND)
: 出力
xq_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xq 格子上の配列に対し領域積分を行う

[Source]

    function IntXY_xq(xq_Var)
      ! xq 格子上の配列に対し領域積分を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: IntXY_xq                    ! 出力


      IntXY_xq = IntY_q(a_IntX_xa(xq_Var))

    end function IntXY_xq
Function :
IntXY_xy :real(DBKIND)
: 出力
xy_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xy 格子上の配列に対し領域積分を行う

[Source]

    function IntXY_xy(xy_Var)
      ! xy 格子上の配列に対し領域積分を行う

      real(DBKIND), intent(in) :: xy_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: IntXY_xy                    ! 出力


      IntXY_xy = IntY_y(a_IntX_xa(xy_Var))

    end function IntXY_xy
IntX_p( p_Var ) result(IntX_p)
Function :
IntX_p :real(DBKIND)
: 出力
p_Var(imin:imax) :real(DBKIND), intent(in)
: 入力

整数格子上の配列に対し x 方向に重み付きの積分を行う

Original external subprogram is x_base_module#IntX_p

IntX_x( x_Var ) result(IntX_x)
Function :
IntX_x :real(DBKIND)
: 出力
x_Var(imin:imax) :real(DBKIND), intent(in)
: 入力

半整数格子上の配列に対し x 方向に重み付きの積分を行う

Original external subprogram is x_base_module#IntX_x

IntY_q( q_Var ) result(IntY_q)
Function :
IntY_q :real(DBKIND)
: 出力
q_Var(jmin:jmax) :real(DBKIND), intent(in)
: 入力

整数格子上の配列に対し y 方向に重み付きの積分を行う

Original external subprogram is y_base_module#IntY_q

IntY_y( y_Var ) result(IntY_y)
Function :
IntY_y :real(DBKIND)
: 出力
y_Var(jmin:jmax) :real(DBKIND), intent(in)
: 入力

半整数格子上の配列に対し y 方向に重み付きの積分を行う

Original external subprogram is y_base_module#IntY_y

im
Variable :
im = 10 :integer
: 格子点数

Original external subprogram is x_base_module#im

imax
Variable :
imax :integer
: 配列添字の上限

Original external subprogram is x_base_module#imax

imin
Variable :
imin :integer
: 配列添字の下限

Original external subprogram is x_base_module#imin

jm
Variable :
jm = 10 :integer
: 格子点数

Original external subprogram is y_base_module#jm

jmax
Variable :
jmax :integer
: 配列添字の上限

Original external subprogram is y_base_module#jmax

jmin
Variable :
jmin :integer
: 配列添字の下限

Original external subprogram is y_base_module#jmin

Function :
a_AvrY_ay(imin:imax) :real(DBKIND)
: 出力
ay_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

ay 格子上の配列に対し y 方向に平均を行う

[Source]

    function a_AvrY_ay(ay_Var)
      ! ay 格子上の配列に対し y 方向に平均を行う

      real(DBKIND), intent(in) :: ay_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_AvrY_ay(imin:imax)        ! 出力


      a_AvrY_ay = a_IntY_ay(ay_Var)/sum(y_dy(1:jm))
      
    end function a_AvrY_ay
Function :
a_IntY_ay(imin:imax) :real(DBKIND)
: 出力
ay_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

ay 格子上の配列に対し y 方向に重み付きの積分を行う

[Source]

    function a_IntY_ay(ay_Var)
      ! ay 格子上の配列に対し y 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: ay_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_IntY_ay(imin:imax)        ! 出力
      integer                  :: ix                          ! ループ添字

      ! 初期化
      a_IntY_ay = 0.0d0

      ! 積分
      do ix = imin, imax
        a_IntY_ay(ix) = IntY_y(ay_Var(ix,:))
      end do
      
    end function a_IntY_ay
p_X
Variable :
p_X(:) :real(DBKIND),allocatable
: 整数格子点座標

Original external subprogram is x_base_module#p_X

p_avr_x( x_Var ) result(p_avr_x)
Function :
p_avr_x(imin:imax) :real(DBKIND)
: 出力
x_Var(imin:imax) :real(DBKIND),intent(in)
: 入力

平均操作を行い半整数格子点の配列値を整数格子点上へ返す

Original external subprogram is x_base_module#p_avr_x

p_dx
Variable :
p_dx(:) :real(DBKIND),allocatable
: 整数格子点座標

Original external subprogram is x_base_module#p_dx

Function :
aq_avr_ay(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
ay_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い半整数格子点の配列値を整数格子点上へ返す

[Source]

    function aq_avr_ay(ay_Var)
      ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す
  
      real(DBKIND),intent(in) :: ay_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: aq_avr_ay(imin:imax,jmin:jmax) ! 出力
      integer                 :: ix                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      aq_avr_ay = 0.0d0

      ! 平均操作
      ! * 関数 q_y を用いて計算.
      ! 
      do ix = imin, imax
        aq_avr_ay(ix,:) = q_avr_y(ay_Var(ix,:))
      end do

    end function aq_avr_ay
Function :
pa_avr_xa(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
xa_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い半整数格子点の配列値を整数格子点上へ返す

[Source]

    function pa_avr_xa(xa_Var)
      ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す
  
      real(DBKIND),intent(in) :: xa_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: pa_avr_xa(imin:imax,jmin:jmax) ! 出力
      integer                 :: jy                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      pa_avr_xa = 0.0d0

      ! 平均操作
      ! * 関数 p_x を用いて計算
      !
      do jy = jmin, jmax
        pa_avr_xa(:,jy) = p_avr_x(xa_Var(:,jy))
      end do

    end function pa_avr_xa
Function :
ay_avr_aq(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
aq_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い整数格子点の配列値を半整数格子点上へ返す

[Source]

    function ay_avr_aq(aq_Var)
      ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す
  
      real(DBKIND),intent(in) :: aq_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: ay_avr_aq(imin:imax,jmin:jmax) ! 出力
      integer                 :: ix                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      ay_avr_aq = 0.0d0

      ! 平均操作
      ! * 関数 y_q を用いて計算
      !
      do ix = imin, imax
        ay_avr_aq(ix,:) = y_avr_q(aq_Var(ix,:))
      end do

    end function ay_avr_aq
Function :
pa_avr_xa(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
xa_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い半整数格子点の配列値を整数格子点上へ返す

[Source]

    function pa_avr_xa(xa_Var)
      ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す
  
      real(DBKIND),intent(in) :: xa_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: pa_avr_xa(imin:imax,jmin:jmax) ! 出力
      integer                 :: jy                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      pa_avr_xa = 0.0d0

      ! 平均操作
      ! * 関数 p_x を用いて計算
      !
      do jy = jmin, jmax
        pa_avr_xa(:,jy) = p_avr_x(xa_Var(:,jy))
      end do

    end function pa_avr_xa
py_dx
Variable :
py_dx(:,:) :real(DBKIND),allocatable
: x 方向格子間隔(整数格子)
Function :
a_AvrX_xa(jmin:jmax) :real(DBKIND)
: 出力
xa_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xa 格子上の配列に対し x 方向に平均を行う

[Source]

    function a_AvrX_xa(xa_Var)
      ! xa 格子上の配列に対し x 方向に平均を行う

      real(DBKIND), intent(in) :: xa_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_AvrX_xa(jmin:jmax)        ! 出力


      a_AvrX_xa = a_IntX_xa(xa_Var)/sum(x_dx(1:im))
      
    end function a_AvrX_xa
Function :
a_IntX_xa(jmin:jmax) :real(DBKIND)
: 出力
xa_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xa 格子上の配列に対し x 方向に重み付きの積分を行う

[Source]

    function a_IntX_xa(xa_Var)
      ! xa 格子上の配列に対し x 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: xa_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_IntX_xa(jmin:jmax)        ! 出力
      integer                  :: jy                          ! ループ添字

      ! 初期化
      a_IntX_xa = 0.0d0

      ! 積分
      do jy = jmin, jmax
        a_IntX_xa(jy) = IntX_x(xa_Var(:,jy))
      end do
      
    end function a_IntX_xa
q_Y
Variable :
q_Y(:) :real(DBKIND),allocatable
: 整数格子点座標

Original external subprogram is y_base_module#q_Y

q_avr_y( y_Var ) result(q_avr_y)
Function :
q_avr_y(jmin:jmax) :real(DBKIND)
: 出力
y_Var(jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い半整数格子点の配列値を整数格子点上へ返す

Original external subprogram is y_base_module#q_avr_y

q_dy
Variable :
q_dy(:) :real(DBKIND),allocatable
: 整数格子点座標

Original external subprogram is y_base_module#q_dy

Function :
x_AvrY_xq(imin:imax) :real(DBKIND)
: 出力
xq_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xq 格子上の配列に対し y 方向に平均を行う

[Source]

    function x_AvrY_xq(xq_Var)
      ! xq 格子上の配列に対し y 方向に平均を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: x_AvrY_xq(imin:imax)        ! 出力


      x_AvrY_xq = x_IntY_xq(xq_Var)/sum(q_dy(1:jm))
      
    end function x_AvrY_xq
Function :
a_AvrY_ay(imin:imax) :real(DBKIND)
: 出力
ay_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

ay 格子上の配列に対し y 方向に平均を行う

[Source]

    function a_AvrY_ay(ay_Var)
      ! ay 格子上の配列に対し y 方向に平均を行う

      real(DBKIND), intent(in) :: ay_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_AvrY_ay(imin:imax)        ! 出力


      a_AvrY_ay = a_IntY_ay(ay_Var)/sum(y_dy(1:jm))
      
    end function a_AvrY_ay
Function :
x_IntY_xq(imin:imax) :real(DBKIND)
: 出力
xq_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

ay 格子上の配列に対し y 方向に重み付きの積分を行う

[Source]

    function x_IntY_xq(xq_Var)
      ! ay 格子上の配列に対し y 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: x_IntY_xq(imin:imax)        ! 出力
      integer                  :: ix                          ! ループ添字

      ! 初期化
      x_IntY_xq = 0.0d0

      ! 積分
      do ix = imin, imax
        x_IntY_xq(ix) = IntY_q(xq_Var(ix,:))
      end do
      
    end function x_IntY_xq
Function :
a_IntY_ay(imin:imax) :real(DBKIND)
: 出力
ay_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

ay 格子上の配列に対し y 方向に重み付きの積分を行う

[Source]

    function a_IntY_ay(ay_Var)
      ! ay 格子上の配列に対し y 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: ay_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_IntY_ay(imin:imax)        ! 出力
      integer                  :: ix                          ! ループ添字

      ! 初期化
      a_IntY_ay = 0.0d0

      ! 積分
      do ix = imin, imax
        a_IntY_ay(ix) = IntY_y(ay_Var(ix,:))
      end do
      
    end function a_IntY_ay
x_X
Variable :
x_X(:) :real(DBKIND),allocatable
: 半整数格子点座標

Original external subprogram is x_base_module#x_X

x_avr_p( p_Var ) result(x_avr_p)
Function :
x_avr_p(imin:imax) :real(DBKIND)
: 出力
p_Var(imin:imax) :real(DBKIND),intent(in)
: 入力

平均操作を行い整数格子点の配列値を半整数格子点上へ返す

Original external subprogram is x_base_module#x_avr_p

x_dx
Variable :
x_dx(:) :real(DBKIND),allocatable
: 半整数格子点座標

Original external subprogram is x_base_module#x_dx

xmargin
Variable :
xmargin = 2 :integer
: 糊代格子点数

Original external subprogram is x_base_module#xmargin

Function :
xa_avr_pa(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
pa_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い整数格子点の配列値を半整数格子点上へ返す

[Source]

    function xa_avr_pa(pa_Var)
      ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す
  
      real(DBKIND),intent(in) :: pa_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: xa_avr_pa(imin:imax,jmin:jmax) ! 出力
      integer                 :: jy                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      xa_avr_pa = 0.0d0

      ! 平均操作
      ! * 関数 x_avr_p を用いて計算
      !
      do jy = jmin, jmax
        xa_avr_pa(:,jy) = x_avr_p(pa_Var(:,jy))
      end do

    end function xa_avr_pa
Function :
aq_avr_ay(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
ay_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い半整数格子点の配列値を整数格子点上へ返す

[Source]

    function aq_avr_ay(ay_Var)
      ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す
  
      real(DBKIND),intent(in) :: ay_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: aq_avr_ay(imin:imax,jmin:jmax) ! 出力
      integer                 :: ix                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      aq_avr_ay = 0.0d0

      ! 平均操作
      ! * 関数 q_y を用いて計算.
      ! 
      do ix = imin, imax
        aq_avr_ay(ix,:) = q_avr_y(ay_Var(ix,:))
      end do

    end function aq_avr_ay
xq_dy
Variable :
xq_dy(:,:) :real(DBKIND),allocatable
: y 方向格子間隔(整数格子)
xy_X
Variable :
xy_X(:,:) :real(DBKIND),allocatable
: x 座標(半整数格子)
xy_Y
Variable :
xy_Y(:,:) :real(DBKIND),allocatable
: y 座標(半整数格子)
Function :
xa_avr_pa(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
pa_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い整数格子点の配列値を半整数格子点上へ返す

[Source]

    function xa_avr_pa(pa_Var)
      ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す
  
      real(DBKIND),intent(in) :: pa_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: xa_avr_pa(imin:imax,jmin:jmax) ! 出力
      integer                 :: jy                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      xa_avr_pa = 0.0d0

      ! 平均操作
      ! * 関数 x_avr_p を用いて計算
      !
      do jy = jmin, jmax
        xa_avr_pa(:,jy) = x_avr_p(pa_Var(:,jy))
      end do

    end function xa_avr_pa
Function :
ay_avr_aq(imin:imax,jmin:jmax) :real(DBKIND)
: 出力
aq_Var(imin:imax,jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い整数格子点の配列値を半整数格子点上へ返す

[Source]

    function ay_avr_aq(aq_Var)
      ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す
  
      real(DBKIND),intent(in) :: aq_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: ay_avr_aq(imin:imax,jmin:jmax) ! 出力
      integer                 :: ix                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      ay_avr_aq = 0.0d0

      ! 平均操作
      ! * 関数 y_q を用いて計算
      !
      do ix = imin, imax
        ay_avr_aq(ix,:) = y_avr_q(aq_Var(ix,:))
      end do

    end function ay_avr_aq
Subroutine :
i :integer,intent(in)
: x 方向格子点数
j :integer,intent(in)
: y 方向格子点数
xmg :integer,intent(in)
: x 方向糊代格子点数
ymg :integer,intent(in)
: y 方向糊代格子点数
xmin :real(DBKIND),intent(in)
: x 座標最小値
xmax :real(DBKIND),intent(in)
: x 座標最大値
ymin :real(DBKIND),intent(in)
: y 座標最小値
ymax :real(DBKIND),intent(in)
: y 座標最大値

xy 方向の座標値と格子点間隔を設定する

[Source]

    subroutine xy_axis_init(i, j, xmg, ymg, xmin, xmax, ymin, ymax)
      ! xy 方向の座標値と格子点間隔を設定する

      integer,intent(in) :: i     ! x 方向格子点数
      integer,intent(in) :: j     ! y 方向格子点数
      integer,intent(in) :: xmg   ! x 方向糊代格子点数
      integer,intent(in) :: ymg   ! y 方向糊代格子点数
      real(DBKIND),intent(in)     :: xmin  ! x 座標最小値     
      real(DBKIND),intent(in)     :: xmax  ! x 座標最大値  
      real(DBKIND),intent(in)     :: ymin  ! y 座標最小値     
      real(DBKIND),intent(in)     :: ymax  ! y 座標最大値  

      ! 配列の上下限の値, 座標値と格子点間隔を設定
      ! * 1 次元用のサブルーチンを用いる
      !
      call x_axis_init(i, xmg, xmin, xmax)
      call y_axis_init(j, ymg, ymin, ymax)

      ! 2 次元座標配列の設定
      ! * 組み込み関数 spread を用いて 1 次元配列を 2 次元に拡張する.
      ! 
      allocate(xy_X(imin:imax,jmin:jmax))
      allocate(xy_Y(imin:imax,jmin:jmax))
      allocate(xy_dx(imin:imax,jmin:jmax))
      allocate(py_dx(imin:imax,jmin:jmax))
      allocate(xy_dy(imin:imax,jmin:jmax))
      allocate(xq_dy(imin:imax,jmin:jmax))

      xy_X  = spread(x_X,2,size(y_Y))
      xy_Y  = spread(y_Y,1,size(x_X))
      xy_dx = spread(x_dx,2,size(y_Y))
      py_dx = spread(p_dx,2,size(y_Y))
      xy_dy = spread(y_dy,1,size(x_X))
      xq_dy = spread(q_dy,1,size(x_X))

    end subroutine xy_axis_init  
xy_dx
Variable :
xy_dx(:,:) :real(DBKIND),allocatable
: x 方向格子間隔(半整数格子)
xy_dy
Variable :
xy_dy(:,:) :real(DBKIND),allocatable
: y 方向格子間隔(半整数格子)
Function :
y_AvrX_py(jmin:jmax) :real(DBKIND)
: 出力
py_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

py 格子上の配列に対し x 方向に平均を行う

[Source]

    function y_AvrX_py(py_Var)
      ! py 格子上の配列に対し x 方向に平均を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: y_AvrX_py(jmin:jmax)        ! 出力


      y_AvrX_py = y_IntX_py(py_Var)/sum(p_dx(1:im))
      
    end function y_AvrX_py
Function :
a_AvrX_xa(jmin:jmax) :real(DBKIND)
: 出力
xa_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xa 格子上の配列に対し x 方向に平均を行う

[Source]

    function a_AvrX_xa(xa_Var)
      ! xa 格子上の配列に対し x 方向に平均を行う

      real(DBKIND), intent(in) :: xa_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_AvrX_xa(jmin:jmax)        ! 出力


      a_AvrX_xa = a_IntX_xa(xa_Var)/sum(x_dx(1:im))
      
    end function a_AvrX_xa
Function :
y_IntX_py(jmin:jmax) :real(DBKIND)
: 出力
py_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

py 格子上の配列に対し x 方向に重み付きの積分を行う

[Source]

    function y_IntX_py(py_Var)
      ! py 格子上の配列に対し x 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: y_IntX_py(jmin:jmax)        ! 出力
      integer                  :: jy                          ! ループ添字

      ! 初期化
      y_IntX_py = 0.0d0

      ! 積分
      do jy = jmin, jmax
        y_IntX_py(jy) = IntX_p(py_Var(:,jy))
      end do
      
    end function y_IntX_py
Function :
a_IntX_xa(jmin:jmax) :real(DBKIND)
: 出力
xa_Var(imin:imax,jmin:jmax) :real(DBKIND), intent(in)
: 入力

xa 格子上の配列に対し x 方向に重み付きの積分を行う

[Source]

    function a_IntX_xa(xa_Var)
      ! xa 格子上の配列に対し x 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: xa_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_IntX_xa(jmin:jmax)        ! 出力
      integer                  :: jy                          ! ループ添字

      ! 初期化
      a_IntX_xa = 0.0d0

      ! 積分
      do jy = jmin, jmax
        a_IntX_xa(jy) = IntX_x(xa_Var(:,jy))
      end do
      
    end function a_IntX_xa
y_Y
Variable :
y_Y(:) :real(DBKIND),allocatable
: 半整数格子点座標

Original external subprogram is y_base_module#y_Y

y_avr_q( q_Var ) result(y_avr_q)
Function :
y_avr_q(jmin:jmax) :real(DBKIND)
: 出力
q_Var(jmin:jmax) :real(DBKIND),intent(in)
: 入力

平均操作を行い整数格子点の配列値を半整数格子点上へ返す

Original external subprogram is y_base_module#y_avr_q

y_dy
Variable :
y_dy(:) :real(DBKIND),allocatable
: 半整数格子点座標

Original external subprogram is y_base_module#y_dy

ymargin
Variable :
ymargin = 2 :integer
: 糊代格子点数

Original external subprogram is y_base_module#ymargin

[Validate]