average_vector Function

    function average_vector(this)
        ! Returns average of a vector
        ! General routine, does not require global data

        real(dp) :: this(:)
        integer :: length
        integer :: i
        real(dp) :: s ! sum of array elements
        real(dp) :: average_vector

        s = 0
        length = size(this, 1)
        do i = 1, length
            s = s + this(i)
        end do
        if (length == 0) then
            average_vector = 0
        else
            average_vector = s / length
        end if

    end function average_vector