| Class | gauss_quad |
| In: |
radiation/gauss_quad.f90
|
Note that Japanese and English are described in parallel.
| !$ ! 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) |
| 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) |
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
| Constant : | |||
| version = ’$Name: dcpam5-20101015 $’ // ’$Id: gauss_quad.f90,v 1.1 2010-09-12 02:11:35 yot Exp $’ : | character(*), parameter
|