Class w_base_module
In: src/w_base_module.f90

w_base_module

 spml/w_base_module モジュールは球面上での 2 次元流体運動を
 球面調和函数を用いたスペクトル法によって数値計算するための
 モジュール w_module の下部モジュールであり, スペクトル法の
 基本的なな Fortran90 関数を提供する.

 内部で ISPACK の SPPACK と SNPACK の Fortran77 サブルーチンを呼んでいる.
 スペクトルデータおよび格子点データの格納方法や変換の詳しい計算法に
 ついては ISPACK/SNPACK,SPPACK のマニュアルを参照されたい.

Methods

a   ia   im   ip   it   jm   l_nm   l_nm   l_nm   l_nm   nm   nm_l   nm_l   np   openmp   p   r   t   w_base_Initial   w_xy   x_Lon   x_Lon_Weight   xy_Lat   xy_Lon   xy_w   y   y_Lat   y_Lat_Weight  

Included Modules

dc_message

Public Instance methods

a
Variable :
a(:) :real(8), allocatable
: 変換用配列
ia
Variable :
ia(:) :integer, allocatable
: 変換用配列
im
Variable :
im =64 :integer
: 格子点の設定(東西)
ip
Variable :
ip(:) :integer, allocatable
: 変換用配列
it
Variable :
it(6) :integer
: 変換用配列
jm
Variable :
jm =32 :integer
: 格子点の設定(南北)
Function :
l_nm_array00 :integer
: (out) スペクトルデータの格納位置
n :integer, intent(in)
: (in) 全波数
m :integer, intent(in)
: (in) 帯状波数

全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.

引数 n,m がともに整数値の場合, 整数値を返す.

[Source]

    function l_nm_array00(n,m)
      !
      ! 全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.
      ! 
      ! 引数 n,m がともに整数値の場合, 整数値を返す. 
      !
      integer               :: l_nm_array00   
      !(out) スペクトルデータの格納位置 

      integer, intent(in)   :: n     !(in) 全波数
      integer, intent(in)   :: m     !(in) 帯状波数           

      call snnm2l(n,m,l_nm_array00)
    end function l_nm_array00
Function :
l_nm_array01(size(marray)) :integer
: (out) スペクトルデータ位置
n :integer, intent(in)
: (in) 全波数
marray(:) :integer, intent(in)
: (in) 帯状波数

スペクトルデータの格納位置

全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.

第 1 引数 n が整数, 第 2 引数 marray が整数 1 次元配列の場合, marray と同じ大きさの 1 次元整数配列を返す.

[Source]

    function l_nm_array01(n,marray)           ! スペクトルデータの格納位置 
      !
      ! 全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.
      ! 
      ! 第 1 引数 n が整数, 第 2 引数 marray が整数 1 次元配列の場合, 
      ! marray と同じ大きさの 1 次元整数配列を返す. 
      !
      integer, intent(in)  :: n               !(in) 全波数
      integer, intent(in)  :: marray(:)       !(in) 帯状波数
      integer              :: l_nm_array01(size(marray))
      !(out) スペクトルデータ位置

      integer              :: i 

      do i=1, size(marray)
         l_nm_array01(i) = l_nm_array00(n,marray(i))
      enddo
    end function l_nm_array01
Function :
l_nm_array10(size(narray)) :integer
: (out) スペクトルデータ位置
narray(:) :integer, intent(in)
: (in) 全波数
m :integer, intent(in)
: (in) 帯状波数

全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.

第 1 引数 narray が整数 1 次元配列, 第 2 引数 m が整数の場合, narray と同じ大きさの 1 次元整数配列を返す.

[Source]

    function l_nm_array10(narray,m)
      !
      ! 全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.
      ! 
      ! 第 1 引数 narray が整数 1 次元配列, 第 2 引数  m が整数の場合, 
      ! narray と同じ大きさの 1 次元整数配列を返す. 
      !
      integer, intent(in)  :: narray(:)           !(in) 全波数  
      integer, intent(in)  :: m                   !(in) 帯状波数
      integer              :: l_nm_array10(size(narray))
      !(out) スペクトルデータ位置

      integer              :: i 

      do i=1, size(narray)
         l_nm_array10(i) = l_nm_array00(narray(i),m)
      enddo
    end function l_nm_array10
Function :
l_nm_array11(size(narray)) :integer
: (out) スペクトルデータ位置
narray(:) :integer, intent(in)
: (in) 全波数
marray(:) :integer, intent(in)
: (in) 帯状波数

全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.

第 1,2 引数 narray, marray がともに整数 1 次元配列の場合, narray, marray と同じ大きさの 1 次元整数配列を返す. narray, marray は同じ大きさでなければならない.

[Source]

    function l_nm_array11(narray,marray)
      !
      ! 全波数(n)と東西波数(m)からそのスペクトルデータの格納位置を返す.
      ! 
      ! 第 1,2 引数 narray, marray がともに整数 1 次元配列の場合, 
      ! narray, marray と同じ大きさの 1 次元整数配列を返す. 
      ! narray, marray は同じ大きさでなければならない. 
      !
      integer, intent(in)  :: narray(:)          !(in) 全波数  
      integer, intent(in)  :: marray(:)          !(in) 帯状波数
      integer              :: l_nm_array11(size(narray))
      !(out) スペクトルデータ位置

      integer              :: i 

      if ( size(narray) .ne. size(marray) ) then
         call MessageNotify('E','l_nm_array11', 'dimensions of input arrays  n and m are different.')
      endif

      do i=1, size(narray)
         l_nm_array11(i) = l_nm_array00(narray(i),marray(i))
      enddo
    end function l_nm_array11
nm
Variable :
nm =21 :integer
: 切断波数の設定
Function :
nm_l_int(2) :integer
: (out) 全波数, 帯状波数
l :integer, intent(in)
: (in) スペクトルデータの格納位置

スペクトルデータの格納位置(l)から全波数(n)と東西波数(m)を返す.

引数 l が整数値の場合, 対応する全波数と帯状波数を 長さ 2 の 1 次元整数値を返す. nm_l(1) が全波数, nm_l(2) が帯状波数である.

[Source]

    function nm_l_int(l)
      ! 
      ! スペクトルデータの格納位置(l)から全波数(n)と東西波数(m)を返す.
      !
      ! 引数 l が整数値の場合, 対応する全波数と帯状波数を
      ! 長さ 2 の 1 次元整数値を返す. 
      ! nm_l(1) が全波数, nm_l(2) が帯状波数である. 
      !
      integer               :: nm_l_int(2)  !(out) 全波数, 帯状波数
      integer, intent(in)   :: l            !(in) スペクトルデータの格納位置
      
      call snl2nm(l,nm_l_int(1),nm_l_int(2))
    end function nm_l_int
Function :
nm_l_array(size(larray),2) :integer
: (in) スペクトルデータの格納位置
larray(:) :integer, intent(in)
: (out) 全波数, 帯状波数

スペクトルデータの格納位置(l)から全波数(n)と東西波数(m)を返す.

引数 larray が整数 1 次元配列の場合, larray に対応する n, m を格納した 2 次元整数配列を返す. nm_l_array(:,1) が全波数, nm_l_array(:,2) が帯状波数である.

[Source]

    function nm_l_array(larray)
      ! 
      ! スペクトルデータの格納位置(l)から全波数(n)と東西波数(m)を返す.
      !
      ! 引数 larray が整数 1 次元配列の場合, 
      ! larray に対応する n, m を格納した 2 次元整数配列を返す. 
      ! nm_l_array(:,1) が全波数, nm_l_array(:,2) が帯状波数である. 
      !
      integer, intent(in)  :: larray(:)
      !(out) 全波数, 帯状波数

      integer              :: nm_l_array(size(larray),2)
      !(in) スペクトルデータの格納位置

      integer              :: i

      do i=1, size(larray)
         nm_l_array(i,:) = nm_l_int(larray(i))
      enddo
    end function nm_l_array
np
Variable :
np =1 :integer
: OPENMP 最大スレッド数
openmp
Variable :
openmp =.false. :logical
: OPENMP スイッチ
p
Variable :
p(:) :real(8), allocatable
: 変換用配列
r
Variable :
r(:) :real(8), allocatable
: 変換用配列
t
Variable :
t(:) :real(8), allocatable
: 変換用配列
Subroutine :
n_in :integer,intent(in)
: (in) 切断全波数
i_in :integer,intent(in)
: (in) 格子点数(東西)
j_in :integer,intent(in)
: (in) 格子点数(南北)
np_in :integer,intent(in), optional
: (in) OPENMP での最大スレッド数

スペクトル変換の格子点数, 波数および OPENMP 使用時の 最大スレッド数を設定する.

実際の使用には上位サブルーチン w_Initial を用いること.

[Source]

    subroutine w_base_Initial(n_in,i_in,j_in,np_in)
      !
      ! スペクトル変換の格子点数, 波数および OPENMP 使用時の
      ! 最大スレッド数を設定する.
      !
      ! 実際の使用には上位サブルーチン w_Initial を用いること.
      !
      integer,intent(in) :: i_in              !(in) 格子点数(東西)
      integer,intent(in) :: j_in              !(in) 格子点数(南北)
      integer,intent(in) :: n_in              !(in) 切断全波数
      integer,intent(in), optional :: np_in   !(in) OPENMP での最大スレッド数

      integer :: iw, i, j

      im = i_in  
       jm = j_in  
       nm = n_in

      if ( present(np_in) )then
         np = np_in

         if ( np .gt. 1 ) then
            openmp = .true. 
            allocate(wv((nm+4)*(nm+3)*np))
            call MessageNotify('M','w_base_Initial', 'OpenMP computation was set up.')
         else
            openmp = .false. 
         endif

      else
         openmp = .false. 
      endif

      if ( im/2*2 .eq. im ) then
         id = im+1 
      else
         id = im
      endif
      if ( openmp ) then
         jd = jm
      else if ( jm/2*2 .eq. jm ) then
         jd = jm+1
      else
         jd = jm
      endif
      allocate(xy_work(id,jd))                ! 変換用配列

      allocate(t(im*2))                       ! 変換用配列
      allocate(ip(((nm+1)/2+nm+1)*2))         ! 変換用配列
      allocate(p(((nm+1)/2+nm+1)*jm))         ! 変換用配列
      allocate(r(((nm+1)/2*2+3)*(nm/2+1)))    ! 変換用配列
      allocate(ia((nm+1)*(nm+1)*4))           ! 変換用配列
      allocate(a((nm+1)*(nm+1)*6))            ! 変換用配列
      allocate(y(jm/2,4))                     ! 変換用配列

      allocate(q(((nm+1)/2+nm+1)*jm))         ! 作業配列
      if ( openmp ) then
         iw=(im+nm+1)*3*jm/2
      else
         iw=max((nm+4)*(nm+3),jd*3*(nm+1),jd*im)
      endif
      allocate(ws(iw),ww(iw))                 ! 作業用配列

      allocate(x_Lon(im), y_Lat(jm))              ! 格子点座標格納配列
      allocate(x_Lon_Weight(im),y_Lat_Weight(jm)) ! 格子点座標格納配列
      allocate(xy_Lon(im,jm),xy_Lat(im,jm))       ! 格子点座標格納配列

      call sninit(nm,im,jm,it,t,y,ip,p,r,ia,a)

      do i=1,im
         x_Lon(i)  = 2*pi/im*(i-1)              ! 経度座標
         x_Lon_Weight(i) = 2*pi/im              ! 経度座標重み
      enddo

      do j=1,jm/2
         y_Lat(jm/2+j)   =  asin(y(j,1))        ! 緯度座標
         y_Lat(jm/2-j+1) = -asin(y(j,1))        ! 緯度座標
         y_Lat_Weight(jm/2+j)   = 2*y(j,2)      ! 緯度重み(Gauss grid)
         y_Lat_Weight(jm/2-j+1) = 2*y(j,2)      ! 緯度重み(Gauss grid)
      enddo

      do j=1,jm
         xy_Lon(:,j) = x_Lon
      enddo

      do i=1,im
         xy_Lat(i,:) = y_Lat
      enddo

    end subroutine w_base_Initial
Function :
w_xy((nm+1)*(nm+1)) :real(8)
: (out) スペクトルデータ
xy_data(im,jm) :real(8), intent(in)
: (in) 格子点データ
ipow :integer, intent(in), optional
: (in) 変換時に同時に作用させる 1/cosφ の次数. 省略時は 0.
iflag :integer, intent(in), optional
: 変換の種類
    0 : 通常の正変換
    1 : 経度微分を作用させた正変換
   -1 : 緯度微分を作用させた正変換
    2 : sinφを作用させた正変換
  省略時は 0.

格子データからスペクトルデータへ(正)変換する(1 層用).

[Source]

    function w_xy(xy_data,ipow,iflag)
      !
      ! 格子データからスペクトルデータへ(正)変換する(1 層用).
      !
      real(8)               :: w_xy((nm+1)*(nm+1))
      !(out) スペクトルデータ

      real(8), intent(in)   :: xy_data(im,jm)
      !(in) 格子点データ

      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

      logical :: first=.true.                     ! 初回判定スイッチ
      save first

      if (present(ipow)) then
         ipval = ipow
      else
         ipval = ipow_default
      endif

      if (present(iflag)) then
         ifval = iflag
      else
         ifval = iflag_default
      endif

      xy_work(1:im,1:jm)=xy_data

      if ( openmp ) then
         if ( first ) then
            call MessageNotify('M','w_xy', 'OpenMP routine SNTGOS/SNPACK is used for spherical harmonic transformation.')
         endif
         call sntgos(nm,im,id,jm,1,xy_work,w_xy, it,t,y,ip,p,r,ia,a,q,ws,ww,wv,ipval,ifval)
      else
         call sntg2s(nm,im,id,jm,jd,1,xy_work,w_xy, it,t,y,ip,p,r,ia,a,q,ws,ww,ipval,ifval)
      endif
      first = .false.

    end function w_xy
x_Lon
Variable :
x_Lon(:) :real(8), allocatable
: 緯度経度
x_Lon_Weight
Variable :
x_Lon_Weight(:) :real(8), allocatable
: 座標重み
xy_Lat
Variable :
xy_Lat(:,:) :real(8), allocatable
xy_Lon
Variable :
xy_Lon(:,:) :real(8), allocatable
Function :
xy_w(im,jm) :real(8)
: (out) 格子点データ
w_data((nm+1)*(nm+1)) :real(8), intent(in)
: (in) スペクトルデータ
ipow :integer, intent(in), optional
: (in) 作用させる 1/cosφ の次数. 省略時は 0.
iflag :integer, intent(in), optional
: (in) 変換の種類
    0 : 通常の正変換
    1 : 経度微分を作用させた正変換
   -1 : 緯度微分を作用させた正変換
    2 : sinφを作用させた正変換
    省略時は 0.

スペクトルデータから格子データへ変換する(1 層用).

[Source]

    function xy_w(w_data,ipow,iflag)
      !
      ! スペクトルデータから格子データへ変換する(1 層用).
      !
      real(8)               :: xy_w(im,jm)
      !(out) 格子点データ

      real(8), intent(in)   :: w_data((nm+1)*(nm+1))
      !(in) スペクトルデータ

      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

      logical :: first=.true.                    ! 初回判定スイッチ
      save first

      if (present(ipow)) then
         ipval = ipow
      else
         ipval = ipow_default
      endif

      if (present(iflag)) then
         ifval = iflag
      else
         ifval = iflag_default
      endif

      if ( openmp ) then
         if ( first ) then
            call MessageNotify('M','xy_w', 'OpenMP routine SNTSOG/SNPACK is used for spherical harmonic transformation.')
         endif
         call sntsog(nm,im,id,jm,1,w_data,xy_work, it,t,y,ip,p,r,ia,a,q,ws,ww,wv,ipval,ifval)
      else
         call snts2g(nm,im,id,jm,jd,1,w_data,xy_work, it,t,y,ip,p,r,ia,a,q,ws,ww,ipval,ifval)
      endif
      xy_w=xy_work(1:im,1:jm)
      first = .false.

    end function xy_w
y
Variable :
y(:,:) :real(8), allocatable
: 変換用配列
y_Lat
Variable :
y_Lat(:) :real(8), allocatable
: 緯度経度
y_Lat_Weight
Variable :
y_Lat_Weight(:) :real(8), allocatable
: 座標重み

[Validate]