dag_toposort Subroutine

    subroutine dag_toposort(me,order,istat)

    class(dag),intent(inout) :: me
    integer(ip),dimension(:),allocatable,intent(out) :: order  !! the toposort order
    integer(ip),intent(out) :: istat !! Status flag:
                                     !!
                                     !! * 0 if no errors
                                     !! * -1 if circular dependency
                                     !!  (in this case, `order` will not be allocated)

    integer(ip) :: i,iorder

    if (me%n==0) return

    ! initialize internal variables, in case
    ! we have called this routine before.
    call me%init_internal_vars()

    allocate(order(me%n))
    iorder = 0  ! index in order array
    istat = 0   ! no errors so far
    do i=1,me%n
      if (.not. me%vertices(i)%marked) call dfs(me%vertices(i))
      if (istat==-1) exit
    end do

    if (istat==-1) deallocate(order)

    contains

    recursive subroutine dfs(v)

    !! depth-first graph traversal

    type(vertex),intent(inout) :: v
    integer(ip) :: j

    if (istat==-1) return

    if (v%checking) then
      ! error: circular dependency
      istat = -1
    else
      if (.not. v%marked) then
        v%checking = .true.
        if (allocated(v%edges)) then
          do j=1,size(v%edges)
            call dfs(me%vertices(v%edges(j)%ivertex))
            if (istat==-1) return
          end do
        end if
        v%checking = .false.
        v%marked = .true.
        iorder = iorder + 1
        order(iorder) = v%ivertex
      end if
    end if

    end subroutine dfs

    end subroutine dag_toposort