integrate Subroutine

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

Main integration routine for the rk_class.

Type Bound

rk_class

Arguments

Type IntentOptional Attributes Name
class(rk_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) :: tf

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

final state

Calls

proc~~integrate~2~~CallsGraph proc~integrate~2 rk_module::rk_class%integrate step step proc~integrate~2->step

Called by

proc~~integrate~2~~CalledByGraph proc~integrate~2 rk_module::rk_class%integrate proc~compute_halo_monodromy_matrix halo_orbit_module::compute_halo_monodromy_matrix proc~compute_halo_monodromy_matrix->proc~integrate~2 proc~halo_to_rv_diffcorr halo_orbit_module::halo_to_rv_diffcorr proc~halo_to_rv_diffcorr->proc~integrate~2 proc~rk_test rk_module::rk_test proc~rk_test->proc~integrate~2

Source Code

    subroutine integrate(me,t0,x0,h,tf,xf)

    implicit none

    class(rk_class),intent(inout)        :: me
    real(wp),intent(in)                  :: t0    !! initial time
    real(wp),dimension(me%n),intent(in)  :: x0    !! initial state
    real(wp),intent(in)                  :: h     !! abs(time step)
    real(wp),intent(in)                  :: tf    !! final time
    real(wp),dimension(me%n),intent(out) :: xf    !! final state

    real(wp) :: t,dt,t2
    real(wp),dimension(me%n) :: x
    logical :: last,export

    if (.not. associated(me%f)) error stop 'Error in integrate: f is not associated.'

    export = associated(me%report)

    if (export) call me%report(t0,x0)  !first point

    if (h==zero) then
        xf = x0
    else

        t = t0
        x = x0
        dt = sign(h,tf-t0)  !time step (correct sign)
        do
            t2 = t + dt
            last = ((dt>=zero .and. t2>=tf) .or. &  !adjust last time step
                    (dt<zero .and. t2<=tf))         !
            if (last) dt = tf-t                     !
            call me%step(t,x,dt,xf)
            if (last) exit
            if (export) call me%report(t2,xf)   !intermediate point
            x = xf
            t = t2
        end do

    end if

    if (export) call me%report(tf,xf)   !last point

    end subroutine integrate