rk_variable_step_class Derived Type

type, public :: rk_variable_step_class

Main integration class for variable step size Runge-Kutta methods


Inherits

type~~rk_variable_step_class~~InheritsGraph type~rk_variable_step_class rk_variable_step_class type~stepsize_class stepsize_class type~rk_variable_step_class->type~stepsize_class stepsize_method

Inherited by

type~~rk_variable_step_class~~InheritedByGraph type~rk_variable_step_class rk_variable_step_class type~rkf108_class rkf108_class type~rkf108_class->type~rk_variable_step_class type~rkf1210_class rkf1210_class type~rkf1210_class->type~rk_variable_step_class type~rkf1412_class rkf1412_class type~rkf1412_class->type~rk_variable_step_class type~rkf78_class rkf78_class type~rkf78_class->type~rk_variable_step_class type~rkf89_class rkf89_class type~rkf89_class->type~rk_variable_step_class type~rkv89_class rkv89_class type~rkv89_class->type~rk_variable_step_class

Components

Type Visibility Attributes Name Initial
integer, private :: n = 0

user specified number of variables

procedure(deriv_func), private, pointer :: f => null()

user-specified derivative function

procedure(report_func), private, pointer :: report => null()

user-specified report function

procedure(event_func), private, pointer :: g => null()

event function (stop when this is zero)

class(stepsize_class), private, allocatable :: stepsize_method

the method for varying the step size

real(kind=wp), private, dimension(:), allocatable :: rtol

relative tolerance (size(n))

real(kind=wp), private, dimension(:), allocatable :: atol

absolute tolerance (size(n))

integer, private :: p = 0

order of the method

integer, private :: hinit_method = 1

if automatically computing the inital step size, which method to use. 1 = hstart, 2 = hinit.

integer, private :: num_rejected_steps = 0

number of rejected steps


Type-Bound Procedures

procedure, public :: initialize

initialize the class (set n,f, and report)

  • private subroutine initialize(me, n, f, rtol, atol, stepsize_method, hinit_method, report, g)

    Initialize the rk_variable_step_class.

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(inout) :: me
    integer, intent(in) :: n

    number of equations

    procedure(deriv_func) :: f

    derivative function

    real(kind=wp), intent(in), dimension(:) :: rtol

    relative tolerance (if size=1, then same tol used for all equations)

    real(kind=wp), intent(in), dimension(:) :: atol

    absolute tolerance (if size=1, then same tol used for all equations)

    class(stepsize_class), intent(in) :: stepsize_method

    method for varying the step size

    integer, intent(in), optional :: hinit_method

    which method to use for automatic initial step size computation. 1 = use hstart, 2 = use hinit.

    procedure(report_func), optional :: report

    for reporting the steps

    procedure(event_func), optional :: g

    for stopping at an event

procedure, public :: destroy

destructor

procedure, public, non_overridable :: integrate

main integration routine

  • private subroutine integrate(me, t0, x0, h, tf, xf, ierr)

    Main integration routine for the rk_variable_step_class.

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(inout) :: me
    real(kind=wp), intent(in) :: t0

    initial time

    real(kind=wp), intent(in), dimension(:) :: x0

    initial state

    real(kind=wp), intent(in) :: h

    initial abs(time step)

    real(kind=wp), intent(in) :: tf

    final time

    real(kind=wp), intent(out), dimension(:) :: xf

    final state

    integer, intent(out), optional :: ierr

    0 = no errors, <0 = error. if not present, an error will stop program.

procedure, public, non_overridable :: integrate_to_event

integration with event finding

  • private subroutine integrate_to_event(me, t0, x0, h, tmax, tol, tf, xf, gf, ierr)

    Event-finding integration routine for the rk_variable_step_class. Integrates until g(t,x)=0, or until t=tf (whichever happens first).

    Read more…

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(inout) :: me
    real(kind=wp), intent(in) :: t0

    initial time

    real(kind=wp), intent(in), dimension(me%n) :: x0

    initial state

    real(kind=wp), intent(in) :: h

    abs(time step)

    real(kind=wp), intent(in) :: tmax

    max final time if event not located

    real(kind=wp), intent(in) :: tol

    function tolerance for root finding

    real(kind=wp), intent(out) :: tf

    actual final time reached

    real(kind=wp), intent(out), dimension(me%n) :: xf

    final state (at tf)

    real(kind=wp), intent(out) :: gf

    g value at tf

    integer, intent(out), optional :: ierr

    0 = no errors, <0 = error. if not present, an error will stop program.

procedure(step_func), private, deferred :: step

the step routine for the rk method

  • subroutine step_func(me, t, x, h, xf, terr) Prototype

    rk step function

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(inout) :: me
    real(kind=wp), intent(in) :: t

    initial time

    real(kind=wp), intent(in), dimension(me%n) :: x

    initial state vector

    real(kind=wp), intent(in) :: h

    time step

    real(kind=wp), intent(out), dimension(me%n) :: xf

    final state vector

    real(kind=wp), intent(out), dimension(me%n) :: terr

    truncation error estimate

procedure(order_func), private, deferred :: order

returns p, the order of the method

  • pure function order_func(me) result(p) Prototype

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(in) :: me

    Return Value integer

    order of the method

procedure, private :: hstart

for automatically computing the initial step size [this is from DDEABM]

  • private subroutine hstart(me, a, b, y, yprime, etol, h)

    Computes a starting step size to be used in solving initial value problems in ordinary differential equations.

    Read more…

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(inout) :: me
    real(kind=wp), intent(in) :: a

    the initial point of integration.

    real(kind=wp), intent(in) :: b

    a value of the independent variable used to define the direction of integration. a reasonable choice is to set b to the first point at which a solution is desired. you can also use b, if necessary, to restrict the length of the first integration step because the algorithm will not compute a starting step length which is bigger than abs(b-a), unless b has been chosen too close to a. (it is presumed that hstart has been called with b different from a on the machine being used. also see the discussion about the parameter small.)

    real(kind=wp), intent(in), dimension(me%n) :: y

    the vector of initial values of the neq solution components at the initial point a.

    real(kind=wp), intent(in), dimension(me%n) :: yprime

    the vector of derivatives of the neq solution components at the initial point a. (defined by the differential equations in subroutine me%f)

    real(kind=wp), intent(in), dimension(me%n) :: etol

    the vector of error tolerances corresponding to the neq solution components. it is assumed that all elements are positive. following the first integration step, the tolerances are expected to be used by the integrator in an error test which roughly requires that abs(local error) <= etol for each vector component.

    real(kind=wp), intent(out) :: h

    appropriate starting step size to be attempted by the differential equation method.

procedure, private :: hinit

for automatically computing the initial step size [this is from DOP853]

  • private function hinit(me, x, y, posneg, f0, hmax, atol, rtol)

    computation of an initial step size guess

    Read more…

    Arguments

    Type IntentOptional Attributes Name
    class(rk_variable_step_class), intent(inout) :: me
    real(kind=wp), intent(in) :: x
    real(kind=wp), intent(in), dimension(:) :: y

    dimension(n)

    real(kind=wp), intent(in) :: posneg

    posneg = sign(1.0_wp,xend-x)

    real(kind=wp), intent(in), dimension(:) :: f0

    dimension(n)

    real(kind=wp), intent(in) :: hmax
    real(kind=wp), intent(in), dimension(:) :: atol
    real(kind=wp), intent(in), dimension(:) :: rtol

    Return Value real(kind=wp)

Source Code

    type,abstract,public :: rk_variable_step_class

        !! Main integration class for variable step size Runge-Kutta methods

        private

        integer :: n = 0  !! user specified number of variables
        procedure(deriv_func),pointer  :: f      => null()  !! user-specified derivative function
        procedure(report_func),pointer :: report => null()  !! user-specified report function
        procedure(event_func),pointer  :: g      => null()  !! event function (stop when this is zero)

        class(stepsize_class),allocatable :: stepsize_method  !! the method for varying the step size

        real(wp),dimension(:),allocatable :: rtol  !! relative tolerance (`size(n)`)
        real(wp),dimension(:),allocatable :: atol  !! absolute tolerance (`size(n)`)

        integer :: p = 0 !! order of the method

        integer :: hinit_method = 1 !! if automatically computing the inital step size, which
                                    !! method to use. 1 = `hstart`, 2 = `hinit`.

        integer :: num_rejected_steps = 0 !! number of rejected steps

        contains

        private

        procedure,public                 :: initialize         !! initialize the class (set n,f, and report)
        procedure,public                 :: destroy            !! destructor
        procedure,non_overridable,public :: integrate          !! main integration routine
        procedure,non_overridable,public :: integrate_to_event !! integration with event finding
        procedure(step_func),deferred    :: step               !! the step routine for the rk method
        procedure(order_func),deferred   :: order              !! returns `p`, the order of the method
        procedure :: hstart  !! for automatically computing the initial step size [this is from DDEABM]
        procedure :: hinit   !! for automatically computing the initial step size [this is from DOP853]

    end type rk_variable_step_class