!c Description:
!c  スプライン補間式
!c
!c Current Code Owner: 
!c   sugiyama@gfd-dennou.org
!c
!c Histry: 
!c   Version    Date          Comment
!c   -------    ----------    --------
!c   1.0        2003-02-05    作成
!c
!c Copyright (C) SUGIYAMA Ko-ichiro, 2003, All rights reserved

module interface_spl
  interface

     subroutine spl(cef, t0, t1, sw, var, var_int, var_dvt_int)
       use prm_type
       real(8), intent(in)            :: t0
       real(8), intent(in)            :: t1
       real(8), intent(out), optional :: var
       real(8), intent(out), optional :: var_int
       real(8), intent(out), optional :: var_dvt_int
       logical, intent(in)            :: sw
       type(type_spl), intent(in)     :: cef
     end subroutine spl

  end interface
end module interface_spl



subroutine spl(cef, t0, t1, sw, var, var_int, var_dvt_int)
  use prm_type, only: type_spl

  real(8), intent(in)            :: t0
  real(8), intent(in)            :: t1
  real(8), intent(out), optional :: var
  real(8), intent(out), optional :: var_int
  real(8), intent(out), optional :: var_dvt_int
  logical, intent(in)            :: sw
  type(type_spl), intent(in)     :: cef


  !スプライン補間
  if (present(var)) then 
     var = cef%dat &
          & + cef%C * (T1 - T0)  &
          & + cef%D * ((T1 - T0) ** 2.0D0)
     
     !外挿しない場合には式に項が追加される. 
     if (sw) then
        var = var + cef%E * ((T1 - T0) ** 3.0D0)
     end if
  end if


  !スプライン補間の積分
  if(present(var_int)) then 
     var_int = &
          & + (cef%dat * (T1 - T0)) &
          & + cef%C &
          & * (((T1**2.0D0 - T0**2.0D0) / 2.0D0) - cef%temp * (T1 - T0) ) &
          & + cef%D &
          & * ( ((T1**3.0D0 - T0**3.0D0) / 3.0D0) &
          & - (cef%temp * (T1**2.0D0 - T0**2.0D0)) &
          & + (cef%temp**2.0D0 * (T1 - T0)) ) 
     
     !外挿しない場合には式に項が追加される. 
     if (sw) then 
        var_int = var_int& 
             & + cef%E  &
             & * (((T1**4.0D0 - T0**4.0D0) / 4.0D0) &
             & - (cef%temp * (T1**3.0D0 - T0**3.0D0)) &
             & + ((cef%temp**2.0D0 * (T1**2.0D0 - T0**2.0D0)) *1.5D0) &
             & - (cef%temp**3.0D0 * (T1 - T0)))
     end if
  end if


  !変数/積分変数の積分
  if (present(var_dvt_int)) then 
     var_dvt_int = &
          & + cef%dat * dlog(T1/T0) &
          & + cef%C * ((T1 - T0) - cef%temp * dlog(T1/T0)) &
          & + cef%D &
          & * ((T1**2.0D0 - T0**2.0D0)/2.0D0 &
          & - 2.0D0 * cef%temp * (T1 - T0) &
          & + cef%temp**2.0D0 * dlog(T1/T0)) 
     
     !外挿しない場合には補間式に項が追加される. 
     if (sw) then
        var_dvt_int = var_dvt_int &
             & + cef%E &
             & * ((T1**3.0D0 - T0**3.0D0)/3.0D0 &
             & - (3.0D0 * cef%temp &
             & * (T1**2.0D0 - T0**2.0D0))/2.0D0 &
             & + 3.0D0 * cef%temp**2.0D0 * (T1 - T0) & 
             & - cef%temp**3.0D0 * dlog(T1/T0))  
     end if
  end if

end subroutine spl





