Class l_module_sjpack
In: libsrc/l_module_sjpack/l_module_sjpack.f90

l_module_sjpack

Authors:Shin-ichi Takehiro, Youhei SASAKI
Version:$Id: l_module_sjpack.f90 590 2013-08-19 08:48:21Z uwabami $
Copyright&License:See COPYRIGHT

概要

  spml/l_module_sjpack モジュールは球面上での経度方向に一様な帯状的
  1 次元流体運動をルジャンドル多項式を用いたスペクトル法によって
  数値計算するための Fortran90 関数を提供する.

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

  関数, サブルーチンの名前と機能は l_module のものと同じに設計してある.
  したがって use 文を l_module から l_module_sjpack に変更するだけで
  SJPACK の機能が使えるようになる.

関数・変数の名前と型について

命名法

  • 関数名の先頭 (l_, y_) は, 返す値の形を示している.
    l_ :スペクトル(ルジャンドル多項式成分)データ
    y_ :1 次元緯度格子点データ
  • 関数名の間の文字列(GradLat, DivLat, Lapla, LaplaInv)は, その関数の作用を表している.
  • 関数名の最後 (_l, _y) は, 入力変数の形スペクトルデータ および格子点データであることを示している.
    _l :スペクト(ルジャンドル多項式成分)ルデータ
    _y :緯度方向 1 次元格子点データ

各データの種類の説明

  • p : スペクトルデータ.
    • 変数の種類と次元は real(8), dimension(0:nn).
    • nn はルジャンドル多項式の最大次数であり, サブルーチン l_Initial にて あらかじめ設定しておく.
  • y : 緯度方向 1 次元格子点データ.
    • 変数の種類と次元は real(8), dimension(1:jm).
  • l_ で始まる関数が返す値はスペクトルデータに同じ.
  • y_ で始まる関数が返す値は 1 次元格子点データに同じ.
  • スペクトルデータに対する微分等の作用とは, 対応する格子点データに 微分などを作用させたデータをスペクトル変換したものことである.

変数・手続き群の要約

初期化

l_Initial :スペクトル変換の格子点数, 波数, 領域の大きさの設定

終了処理

l_Finalize :モジュールの終了処理(割り付け配列の解放)をおこなう.

座標変数

y_Lat :格子点座標(緯度, 経度座標)を格納した 1 次元配列
y_Lat_Weight :重み座標を格納した 1 次元配列

基本変換

y_l :スペクトルデータから格子データへの変換
l_y :格子データからスペクトルデータへの変換

微分

l_Lapla_l :スペクトルデータにラプラシアンを作用させる
l_LaplaInv_l :スペクトルデータにラプラシアンの逆変換を作用させる
y_GradLat_l :スペクトルデータに勾配型緯度微分∂/∂φを作用させる
l_DivLat_y :格子データに発散型緯度微分 1 /cosφ・∂(g cosφ)/∂φを作用させる

微分(λ,μ=sinφ 座標)

y_GradMu_l :スペクトルデータに 勾配型緯度微分 (1-μ^2)∂/∂μを作用させる
l_DivMu_y :格子データに発散型緯度微分∂/∂μを作用させる

補間

Interpolate_l :スペクトルデータから任意の点での値を求める.

積分・平均

IntLat_y, AvrLat_y :1 次元(Y)格子点データの緯度方向積分および平均

スペクトル解析

作成中

Methods

AvrLat_y   IntLat_y   Interpolate_l   l_DivLat_y   l_DivMu_y   l_Finalize   l_LaplaInv_l   l_Lapla_l   l_initial   l_y   y_GradLat_l   y_GradMu_l   y_Lat   y_Lat_Weight   y_l  

Included Modules

dc_message

Public Instance methods

AvrLat_y( y_data ) result(AvrLat_y)
Function :
AvrLat_y :real(8)
: (out) 平均値
y_data(1:jm) :real(8), intent(in)
: (in) 1 次元緯度格子点データ

1 次元(Y)格子点データの緯度(Y)方向平均.

実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し, y_Y_Weight の総和で割ることで平均している.

[Source]

    function AvrLat_y(y_data)
      !
      ! 1 次元(Y)格子点データの緯度(Y)方向平均.
      !
      ! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算し, 
      ! y_Y_Weight の総和で割ることで平均している. 
      !
      real(8), intent(in) :: y_data(1:jm)    !(in)  1 次元緯度格子点データ
      real(8)             :: AvrLat_y        !(out) 平均値

      AvrLat_y = IntLat_y(y_data)/sum(y_Lat_weight)

    end function AvrLat_y
IntLat_y( y_data ) result(IntLat_y)
Function :
IntLat_y :real(8)
: (out) 積分値
y_data(1:jm) :real(8), intent(in)
: (in) 1 次元緯度(Y)格子点データ

1 次元緯度(Y)格子点データの Y 方向積分.

実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している.

[Source]

    function IntLat_y(y_data)
      !
      ! 1 次元緯度(Y)格子点データの Y 方向積分.
      !
      ! 実際には格子点データ各点毎に y_Y_Weight をかけた総和を計算している. 
      !
      real(8), intent(in) :: y_data(1:jm)    !(in)  1 次元緯度(Y)格子点データ
      real(8)             :: IntLat_y        !(out) 積分値

      IntLat_y = sum(y_data * y_Lat_weight)

    end function IntLat_y
Interpolate_l( l_data, alat ) result(Interpolate_l)
Function :
Interpolate_l :real(8)
: 補間した値
l_data(0:nn) :real(8), intent(IN)
: スペクトルデータ
alat :real(8), intent(IN)
: 補間する位置(緯度)

緯度 alat における関数値を そのルジャンドル変換係数 l_data から補間計算する

[Source]

    function Interpolate_l(l_data,alat)
      !
      ! 緯度 alat における関数値を
      ! そのルジャンドル変換係数 l_data から補間計算する
      !
      real(8), intent(IN) :: l_data(0:nn)           ! スペクトルデータ
      real(8), intent(IN) :: alat                   ! 補間する位置(緯度)
      real(8)             :: Interpolate_l          ! 補間した値
      
      real(8) :: mu
      real(8) :: y0, y1, y2
      integer :: k

      mu = sin(alat)
      Interpolate_l = 0.0D0

      !---- Σa_n^0 L_n^0 の計算
      y2 = 0 
       y1 = 0
      do k=nn,1,-1
         y0 = alpha(k,mu) * y1 + beta(k+1)*y2 + l_data(k)
         y2 = y1 
          y1 = y0
      enddo
      Interpolate_l = beta(1) * y2 + mu*sqrt(3.0D0) * y1 + l_data(0) 

    end function Interpolate_l
l_DivLat_y( y_data ) result(l_DivLat_y)
Function :
l_DivLat_y(0:nn) :real(8)
: (out) 格子点データを発散型緯度微分したスペクトルデータ
y_data(1:jm) :real(8), intent(in)
: (in) 入力格子点データ

格子点データに発散型緯度微分 1/cosφ・∂(f cosφ)/∂φ を作用させて スペクトルデータに変換して返す.

[Source]

    function l_DivLat_y(y_data)
      !
      ! 格子点データに発散型緯度微分 1/cosφ・∂(f cosφ)/∂φ を作用させて
      ! スペクトルデータに変換して返す. 
      !
      real(8)              :: l_DivLat_y(0:nn)
      !(out) 格子点データを発散型緯度微分したスペクトルデータ

      real(8), intent(in)  :: y_data(1:jm)
      !(in) 入力格子点データ

      real(8)              :: work(0:nn+1)
      ! 作業スペクトルデータ

      real(8)  :: q(jm/2*5)                ! 変換用作業配列
      real(8)  :: ws(nn+2)                 ! 変換用作業配列

      call ljtgzs(nm,nn+1,jm,work,y_data,p,q,r,ws,1)
      call ljcyzs(nn,work,l_DivLat_y,c)

    end function l_DivLat_y
l_DivMu_y( y_data ) result(l_DivMu_y)
Function :
l_DivMu_y(0:nn) :real(8)
: (out) 格子点データを発散型緯度微分したスペクトルデータ
y_data(1:jm) :real(8), intent(in)
: (in) 入力格子点データ

格子点データに発散型緯度微分 ∂/∂μ (μ=sinφ)を作用させて スペクトルデータに変換して返す(1 層用).

[Source]

    function l_DivMu_y(y_data)
      !
      ! 格子点データに発散型緯度微分 ∂/∂μ (μ=sinφ)を作用させて
      ! スペクトルデータに変換して返す(1 層用).
      !
      real(8)              :: l_DivMu_y(0:nn)
      !(out) 格子点データを発散型緯度微分したスペクトルデータ

      real(8), intent(in)  :: y_data(1:jm)
      !(in) 入力格子点データ

      l_DivMu_y = l_DivLat_y(y_data/cos(y_Lat))

    end function l_DivMu_y
l_Finalize( )
Subroutine :

モジュールの終了処理(割り付け配列の解放)をおこなう.

解像度を変更する際にはこのサブルーチンを呼んで終了処理を おこなったのちに再度 l_Initial で初期設定しなければ ならない.

[Source]

    subroutine l_Finalize
      !
      ! モジュールの終了処理(割り付け配列の解放)をおこなう. 
      !
      ! 解像度を変更する際にはこのサブルーチンを呼んで終了処理を
      ! おこなったのちに再度 l_Initial で初期設定しなければ
      ! ならない. 
      !
      if ( .not. l_initialize ) then
         call MessageNotify('W','l_finalize', 'l_module_sjpack not initialized yet')
         return
      endif

      deallocate(p)                     ! 変換用配列
      deallocate(r)                     ! 変換用配列
      deallocate(c)                     ! 微分用配列

      deallocate(y_Lat,y_Lat_Weight)    ! 座標変数配列

      l_initialize = .false.

      call MessageNotify('M','l_finalize', 'l_module_sjpack (2013/02/24) is finalized')

    end subroutine l_Finalize
l_LaplaInv_l( l_data ) result(l_LaplaInv_l)
Function :
l_LaplaInv_l(0:nn) :real(8)
: (out) スペクトルデータの逆ラプラシアン
l_data(0:nn) :real(8), intent(in)
: (in) 入力スペクトルデータ

入力スペクトルデータに逆ラプラシアン 

   ▽^{-2}
     =[1/cosφ・∂/∂φ(cosφ∂/∂φ)]^{-1}

を作用する.

スペクトルデータの逆ラプラシアンとは, 対応する格子点データに 逆ラプラシアンを作用させたデータのスペクトル変換のことである.

[Source]

    function l_LaplaInv_l(l_data)
      !
      ! 入力スペクトルデータに逆ラプラシアン
      !
      !    ▽^{-2}
      !      =[1/cosφ・∂/∂φ(cosφ∂/∂φ)]^{-1}
      !
      ! を作用する. 
      !
      ! スペクトルデータの逆ラプラシアンとは, 対応する格子点データに
      ! 逆ラプラシアンを作用させたデータのスペクトル変換のことである. 
      !
      real(8)              :: l_LaplaInv_l(0:nn)
      !(out) スペクトルデータの逆ラプラシアン

      real(8), intent(in)  :: l_data(0:nn)
      !(in) 入力スペクトルデータ

      integer ::  n

      l_LaplaInv_l(0) = 0.0D0
      do n=1,nn
         l_LaplaInv_l(n) = -l_data(n)/(n*(n+1))
      enddo

    end function l_LaplaInv_l
l_Lapla_l( l_data ) result(l_Lapla_l)
Function :
l_Lapla_l(0:nn) :real(8)
: (out) 入力スペクトルデータのラプラシアン
l_data(0:nn) :real(8), intent(in)
: (in) 入力スペクトルデータ

入力スペクトルデータにラプラシアン 

   ▽^2 = 1/cosφ・∂/∂φ(cosφ∂/∂φ)

を作用する.

スペクトルデータのラプラシアンとは, 対応する格子点データに ラプラシアンを作用させたデータのスペクトル変換のことである.

[Source]

    function l_Lapla_l(l_data)
      !
      ! 入力スペクトルデータにラプラシアン
      !
      !    ▽^2 = 1/cosφ・∂/∂φ(cosφ∂/∂φ)
      !
      ! を作用する. 
      !
      ! スペクトルデータのラプラシアンとは, 対応する格子点データに
      ! ラプラシアンを作用させたデータのスペクトル変換のことである. 
      !
      real(8)              :: l_Lapla_l(0:nn)
      !(out) 入力スペクトルデータのラプラシアン

      real(8), intent(in)  :: l_data(0:nn)
      !(in) 入力スペクトルデータ

      integer :: n

      do n=0,nn
         l_Lapla_l(n) = -n*(n+1)*l_data(n)
      enddo

    end function l_Lapla_l
l_initial( n_in, j_in )
Subroutine :
n_in :integer,intent(in)
: (in) 切断波数の設定
j_in :integer,intent(in)
: (in) 格子点数(南北)

スペクトル変換の格子点数, 波数を設定する.

他の関数を呼ぶ前に, 最初にこのサブルーチンを呼んで初期設定を しなければならない.

[Source]

    subroutine l_initial(n_in,j_in)
      !
      ! スペクトル変換の格子点数, 波数を設定する.
      !
      ! 他の関数を呼ぶ前に, 最初にこのサブルーチンを呼んで初期設定を
      ! しなければならない. 
      !
      integer,intent(in) :: j_in              !(in) 格子点数(南北)
      integer,intent(in) :: n_in              !(in) 切断波数の設定

      integer :: j

      jm = j_in 
       nn = n_in 
       nm = n_in+1 
       mm = 1

      allocate(p(jm/2,mm+4))                  ! 変換用配列
      allocate(r((mm+1)*(2*nm-mm-1)+1))       ! 変換用配列
      allocate(c(2*nn+1))                     ! 微分用配列

      allocate(y_Lat(jm),y_Lat_Weight(jm))      ! 座標変数配列

      call ljinit(mm,nm,jm,p,r)

      call ljinic(nn,c)

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

      l_initialize = .true.

      call MessageNotify( 'M','l_initial','l_module_sjpack (2013/02/24) is initialized')

    end subroutine l_initial
l_y( y_data ) result(l_y)
Function :
l_y(0:nn) :real(8)
: (in) スペクトルデータ
y_data(1:jm) :real(8), intent(in)
: (in) 格子点データ

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

[Source]

    function l_y(y_data)
      !
      ! スペクトルデータから格子データへ変換する(1 層用).
      !
      real(8)               :: l_y(0:nn)
      !(in) スペクトルデータ

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


      real(8)  :: q(jm/2*5)               ! 変換a用作業配列
      real(8)  :: ws(nn+1)                ! 変換用作業配列

      call ljtgzs(nm,nn,jm,l_y,y_data,p,q,r,ws,0)

    end function l_y
y_GradLat_l( l_data ) result(y_GradLat_l)
Function :
y_GradLat_l(1:jm) :real(8)
: (out) スペクトルデータを勾配型緯度微分した格子点データ
l_data(0:nn) :real(8), intent(in)
: (in) 入力スペクトルデータ

スペクトルデータに勾配型緯度微分 ∂/∂φ を作用させて 格子点データに変換して返す(1 層用).

[Source]

    function y_GradLat_l(l_data)
      !
      ! スペクトルデータに勾配型緯度微分 ∂/∂φ を作用させて
      ! 格子点データに変換して返す(1 層用).
      !
      real(8)              :: y_GradLat_l(1:jm)
      !(out) スペクトルデータを勾配型緯度微分した格子点データ

      real(8), intent(in)  :: l_data(0:nn)
      !(in) 入力スペクトルデータ

      real(8)              :: work(0:nn+1)
      ! 作業スペクトルデータ

      real(8)  :: q(jm/2*5)                ! 変換用作業配列
      real(8)  :: ws(nn+2)                 ! 変換用作業配列

      call ljcszy(nn,l_data,work,c)
      call ljtszg(nm,nn+1,jm,work,y_GradLat_l,p,q,r,ws,1)

    end function y_GradLat_l
y_GradMu_l( l_data ) result(y_GradMu_l)
Function :
y_GradMu_l(1:jm) :real(8)
: (out) スペクトルデータを勾配型緯度微分した格子点データ
l_data(0:nm) :real(8), intent(in)
: (in) 入力スペクトルデータ

スペクトルデータに勾配型緯度微分 (1-μ^2)∂/∂μ (μ=sinφ) を作用させて格子点データに変換して返す.

[Source]

    function y_GradMu_l(l_data)
      !
      ! スペクトルデータに勾配型緯度微分 (1-μ^2)∂/∂μ  (μ=sinφ)
      ! を作用させて格子点データに変換して返す. 
      !
      real(8)              :: y_GradMu_l(1:jm)
      !(out) スペクトルデータを勾配型緯度微分した格子点データ

      real(8), intent(in)  :: l_data(0:nm)
      !(in) 入力スペクトルデータ

      y_GradMu_l = y_GradLat_l(l_data)*cos(y_Lat)

    end function y_GradMu_l
y_Lat
Variable :
y_Lat(:) :real(8), allocatable
: 緯度経度
y_Lat_Weight
Variable :
y_Lat_Weight(:) :real(8), allocatable
: 座標重み
y_l( l_data ) result(y_l)
Function :
y_l(1:jm) :real(8)
: (out) 格子点データ
l_data(0:nn) :real(8), intent(in)
: (in) スペクトルデータ

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

[Source]

    function y_l(l_data)
      !
      ! スペクトルデータから格子データへ変換する(1 層用).
      !
      real(8)               :: y_l(1:jm)
      !(out) 格子点データ

      real(8), intent(in)   :: l_data(0:nn)
      !(in) スペクトルデータ

      real(8)  :: q(jm/2*5)                ! 変換用作業配列
      real(8)  :: ws(nn+1)                 ! 変換用作業配列

      call ljtszg(nm,nn,jm,l_data,y_l,p,q,r,ws,0)

    end function y_l