ガウス重み, 分点の計算

Calculate Gauss node and Gaussian weight

Note that Japanese and English are described in parallel.

References

Procedures List

!$ ! RadiationFluxDennouAGCM :	放射フラックスの計算
!$ ! RadiationDTempDt :	放射フラックスによる温度変化の計算
!$ ! RadiationFluxOutput :	放射フラックスの出力
!$ ! RadiationFinalize :	終了処理 (モジュール内部の変数の割り付け解除)
!$ ! ———— :	————
!$ ! RadiationFluxDennouAGCM :	Calculate radiation flux
!$ ! RadiationDTempDt :	Calculate temperature tendency with radiation flux
!$ ! RadiationFluxOutput :	Output radiation fluxes
!$ ! RadiationFinalize :	Termination (deallocate variables in this module)

NAMELIST

!$ ! NAMELIST#radiation_DennouAGCM_nml

Methods

GauLeg module_name version

Included Modules

dc_types constants

Public Instance methods

GauLeg( x1, x2, n, x, w )

*Subroutine* :
`x1` :	`real(DP), intent(in )`
`x2` :	`real(DP), intent(in )`
`n` :	`integer , intent(in )`
`x(n)` :	`real(DP), intent(out)`
`w(n)` :	`real(DP), intent(out)`

[Source]

  subroutine GauLeg( x1, x2, n, x, w )

    !
    ! Physical constants settings
    !
    use constants, only: PI   ! $ \pi $ .
                              ! Circular constant

    real(DP), intent(in ) :: x1,x2
    integer , intent(in ) :: n
    real(DP), intent(out) :: x(n),w(n)


    real(DP):: eps
    real(DP):: z1, z, xm, xl, pp, p3, p2, p1
    integer :: i, j, m


    !
    ! changed at 2005/09/14
    !
!!$      eps = 1.0e-11
    eps = 1.0e-15

    m  = ( n+1 ) / 2
    xm = ( x2+x1 ) / 2.0d0
    xl = ( x2-x1 ) / 2.0d0

    do i = 1, m
      z = cos( pi*dble(i-0.25)/(dble(n+0.5)) )

100   continue
      p1 = 1.0d0
      p2 = 0.0d0
!!$         do j = 1, n
      do j = 1, n
        p3 = p2
        p2 = p1
        p1 = ( (2.0d0*dble(j)-1.0d0)*z*p2-(dble(j)-1.0d0)*p3 ) / dble(j)
      enddo

      pp = dble(n) * (z*p1-p2) / (z*z-1.0d0)
      z1 = z
      z  = z1 - p1 / pp

      !
      ! changed at 2005/09/14
      !
!!$         if( abs( z-z1 ) .gt. eps ) go to 100
!!$         if( abs( z-z1 ) / z1 .gt. eps ) go to 100
      !
      ! changed before 2008/08/03
      !
      if( abs( z-z1 ) / abs( z1 + 1.0d-200 ) .gt. eps ) go to 100

      x( i     ) = xm-xl*z
      x( n+1-i ) = xm+xl*z
      w( i     ) = 2.0d0*xl/((1.0d0-z*z)*pp*pp)
      w( n+1-i ) = w(i)
    end do


  end subroutine GauLeg

Private Instance methods

module_name

Constant :

module_name = ‘gauss_quad‘ :

character(*), parameter

:	モジュールの名称. Module name

version

Constant :

version = ’$Name: dcpam5-20110221-2 $’ // ’$Id: gauss_quad.f90,v 1.1 2010-09-12 02:11:35 yot Exp $’ :

character(*), parameter

:	モジュールのバージョン Module version