ddatrp Subroutine

pure subroutine ddatrp(t, tout, yout, ydotout, neq, kold, phi, psi)

Uses

  • proc~~ddatrp~~UsesGraph proc~ddatrp ddatrp module~daskr_kinds daskr_kinds proc~ddatrp->module~daskr_kinds iso_fortran_env iso_fortran_env module~daskr_kinds->iso_fortran_env

The methods in subroutine DDSTP use polynomials to approximate the solution. This routine approximates the solution and its derivative at tout by evaluating one of these polynomials, and its derivative, there. Information defining the polynomial is passed from the caller.

Arguments

Type IntentOptional Attributes Name
real(kind=rk), intent(in) :: t

Current value of the independent variable.

real(kind=rk), intent(in) :: tout

Value of the independent variable at which the solution is desired.

real(kind=rk), intent(out) :: yout(neq)

Interpolated approximation to y at tout.

real(kind=rk), intent(out) :: ydotout(neq)

Interpolated approximation to ydot at tout.

integer, intent(in) :: neq

Problem size.

integer, intent(in) :: kold

Order used on last successful step.

real(kind=rk), intent(in) :: phi(neq,kold+1)

Array of scaled divided differences of y.

real(kind=rk), intent(in) :: psi(kold+1)

Array of past stepsize history.


Variables

Type Visibility Attributes Name Initial
integer, public :: j
real(kind=rk), public :: c
real(kind=rk), public :: d
real(kind=rk), public :: gama
real(kind=rk), public :: dt

Source Code

pure subroutine ddatrp(t, tout, yout, ydotout, neq, kold, phi, psi)
!! The methods in subroutine [[ddstp]] use polynomials to approximate the solution.
!! This routine approximates the solution and its derivative at `tout` by evaluating
!! one of these polynomials, and its derivative, there. Information defining the
!! polynomial is passed from the caller.

   use daskr_kinds, only: rk, zero, one
   implicit none

   real(rk), intent(in) :: t
      !! Current value of the independent variable.
   real(rk), intent(in) :: tout
      !! Value of the independent variable at which the solution is desired.
   real(rk), intent(out) :: yout(neq)
      !! Interpolated approximation to `y` at `tout`.
   real(rk), intent(out) :: ydotout(neq)
      !! Interpolated approximation to `ydot` at `tout`.
   integer, intent(in) :: neq
      !! Problem size.
   integer, intent(in) :: kold
      !! Order used on last successful step.
   real(rk), intent(in) :: phi(neq, kold + 1)
      !! Array of scaled divided differences of `y`.
   real(rk), intent(in) :: psi(kold + 1)
      !! Array of past stepsize history.

   integer :: j
   real(rk) :: c, d, gama, dt

   dt = tout - t
   yout = phi(:, 1)
   ydotout = zero
   c = one
   d = zero
   gama = dt/psi(1)

   do j = 2, kold + 1
      d = d*gama + c/psi(j - 1)
      c = c*gama
      gama = (dt + psi(j - 1))/psi(j)
      yout = yout + c*phi(:, j)
      ydotout = ydotout + d*phi(:, j)
   end do

end subroutine ddatrp