ddeabm_wrapper Subroutine

private subroutine ddeabm_wrapper(me, t, y, tout, tstop, idid, integration_mode, tstep)

Wrapper routine for ddeabm.

If starting a new problem, must first call me%first_call().

Note

Currently not using the recommended tols if idid=-2.

Type Bound

ddeabm_class

Arguments

Type IntentOptional Attributes Name
class(ddeabm_class), intent(inout) :: me
real(kind=wp), intent(inout) :: t
real(kind=wp), intent(inout), dimension(me%neq) :: y
real(kind=wp), intent(in) :: tout

time at which a solution is desired.

real(kind=wp), intent(in), optional :: tstop

for some problems it may not be permissible to integrate past a point tstop because a discontinuity occurs there or the solution or its derivative is not defined beyond tstop. when you have declared a tstop point (see info(4)), you have told the code not to integrate past tstop. in this case any tout beyond tstop is invalid input. [not used if not present]

integer, intent(out) :: idid

indicator reporting what the code did. you must monitor this integer variable to decide what action to take next.

integer, intent(in), optional :: integration_mode

Step mode: 1 - normal integration from t to tout, no reporting [default]. 2 - normal integration from t to tout, report each step.

real(kind=wp), intent(in), optional :: tstep

Fixed time step to use for integration_mode=2. If not present, then default integrator steps are used. If integration_mode=1, then this is ignored.


Calls

proc~~ddeabm_wrapper~~CallsGraph proc~ddeabm_wrapper ddeabm_class%ddeabm_wrapper proc~ddeabm ddeabm_class%ddeabm proc~ddeabm_wrapper->proc~ddeabm proc~report_error report_error proc~ddeabm_wrapper->proc~report_error proc~ddeabm->proc~report_error proc~ddes ddeabm_class%ddes proc~ddeabm->proc~ddes proc~ddes->proc~report_error proc~dintp dintp proc~ddes->proc~dintp proc~dsteps ddeabm_class%dsteps proc~ddes->proc~dsteps proc~dhstrt ddeabm_class%dhstrt proc~dsteps->proc~dhstrt proc~dhvnrm dhvnrm proc~dhstrt->proc~dhvnrm

Source Code

   subroutine ddeabm_wrapper(me, t, y, tout, tstop, idid, integration_mode, tstep)

      implicit none

      class(ddeabm_class), intent(inout)        :: me
      real(wp), intent(inout)                   :: t
      real(wp), dimension(me%neq), intent(inout) :: y
      real(wp), intent(in)                      :: tout    !! time at which a solution is desired.
      real(wp), intent(in), optional             :: tstop   !! for some problems it may not be permissible to integrate
                                                        !! past a point `tstop` because a discontinuity occurs there
                                                        !! or the solution or its derivative is not defined beyond
                                                        !! `tstop`.  when you have declared a `tstop` point (see `info(4)`),
                                                        !! you have told the code not to integrate past `tstop`.
                                                        !! in this case any `tout` beyond `tstop` is invalid input.
                                                        !! [not used if not present]
      integer, intent(out)                      :: idid    !! indicator reporting what the code did.
                                                        !! you must monitor this integer variable to
                                                        !! decide what action to take next.
      integer, intent(in), optional :: integration_mode     !! Step mode:
                                                        !! *1* - normal integration from `t` to `tout`, no reporting [default].
                                                        !! *2* - normal integration from `t` to `tout`, report each step.
      real(wp), intent(in), optional :: tstep    !! Fixed time step to use for `integration_mode=2`.
                                             !! If not present, then default integrator steps are used.
                                             !! If `integration_mode=1`, then this is ignored.

      integer :: mode                 !! local copy of `integration_mode`
      logical :: fixed_output_step    !! if reporting output at a fixed step size
      real(wp) :: t2                  !! for fixed step size: an intermediate step
      logical :: last                 !! for fixed step size: the last step
      real(wp) :: dt                  !! for fixed step size: actual signed time step
      real(wp) :: direction           !! direction of integration for
                                    !! fixed step size: `+1: dt>=0`, `-1: dt<0`

      !optional input:
      if (present(integration_mode)) then
         mode = integration_mode
      else
         mode = 1  !default
      end if

      ! if we are reporting the output at a fixed step size:
      fixed_output_step = present(tstep) .and. mode == 2
      if (fixed_output_step) then
         direction = sign(1.0_wp, tout - t)
         dt = direction*abs(tstep)
         last = .false.
      end if

      !check for invalid inputs:
      if (mode /= 1 .and. mode /= 2) then
         call report_error('ddeabm_wrapper', 'integration_mode must be 1 or 2.', 0, 0)
         idid = -33
         return
      end if
      if ((mode == 2) .and. (.not. associated(me%report))) then
         call report_error('ddeabm_wrapper', 'REPORT procedure must be associated for integration_mode=2.', 0, 0)
         idid = -33
         return
      end if

      !set info array:

      !info(1) is set when ddeabm_new_problem is called

      !info(2)
      if (me%scalar_tols) then
         me%info(2) = 0
      else
         me%info(2) = 1
      end if

      !info(3)
      if (mode == 2 .and. .not. fixed_output_step) then  !requires the intermediate steps
         me%info(3) = 1
      else
         me%info(3) = 0
      end if

      !info(4)
      if (present(tstop)) then
         me%info(4) = 1
         me%tstop = tstop
      else
         me%info(4) = 0
         me%tstop = 0.0_wp !not used
      end if

      !make a copy of the tols, since the routine might change them:
      me%rtol_tmp = me%rtol
      me%atol_tmp = me%atol

      select case (mode)
      case (1)  !normal integration with no reporting

         !call the lower-level routine:
         call me%ddeabm(neq=me%neq, &
                        t=t, &
                        y=y, &
                        tout=tout, &
                        info=me%info, &
                        rtol=me%rtol_tmp, &
                        atol=me%atol_tmp, &
                        idid=idid)

         !if there was a user-triggered error:
         if (me%error) idid = -1000

      case (2) !normal integration, reporting the default intermediate points

         call me%report(t, y)  !initial point

         do

            if (fixed_output_step) then
               t2 = t + dt  ! take one step to t2 and return
               last = direction*(tout - t2) <= 0.0_wp  ! if last point
               if (last) t2 = tout  ! adjust last time step if necessary
            else
               t2 = tout  ! take a step in the direction of tout and return
            end if

            !take one step (note: info(3)=1 for this case):
            call me%ddeabm(neq=me%neq, &
                           t=t, &
                           y=y, &
                           tout=t2, &
                           info=me%info, &
                           rtol=me%rtol_tmp, &
                           atol=me%atol_tmp, &
                           idid=idid)

            !if there was a user-triggered error:
            if (me%error) then
               idid = -1000
               exit
            end if

            !check status (see ddeabm or idid codes):
            if (idid > 0) call me%report(t, y)  ! report a successful intermediate or final step

            if (fixed_output_step) then
               if (last) exit  ! exit if finished
               if (.not. (idid == 2 .or. idid == 3)) exit  ! exit on error
            else
               if (idid /= 1) exit      ! exit if finished, or if there was an error
            end if

         end do

      end select

   end subroutine ddeabm_wrapper