backtracking_linesearch Subroutine

private subroutine backtracking_linesearch(me, xold, p, x, f, fvec, fjac, fjac_sparse)

Backtracking line search.

See also

Note

Either fjac or fjac_sparse should be present.

Arguments

Type IntentOptional Attributes Name
class(nlesolver_type), intent(inout) :: me
real(kind=wp), intent(in), dimension(me%n) :: xold

previous value of x

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

search direction

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

new x

real(kind=wp), intent(inout) :: f

magnitude of fvec

real(kind=wp), intent(inout), dimension(me%m) :: fvec

function vector

real(kind=wp), intent(in), optional, dimension(:,:) :: fjac

jacobian matrix [dense]

real(kind=wp), intent(in), optional, dimension(:) :: fjac_sparse

jacobian matrix [sparse]


Source Code

    subroutine backtracking_linesearch(me,xold,p,x,f,fvec,fjac,fjac_sparse)

    implicit none

    class(nlesolver_type),intent(inout) :: me
    real(wp),dimension(me%n),intent(in) :: xold      !! previous value of `x`
    real(wp),dimension(me%n),intent(in) :: p         !! search direction
    real(wp),dimension(me%n),intent(out) :: x        !! new `x`
    real(wp),intent(inout) :: f                      !! magnitude of `fvec`
    real(wp),dimension(me%m),intent(inout) :: fvec   !! function vector
    real(wp),dimension(:,:),intent(in),optional :: fjac !! jacobian matrix [dense]
    real(wp),dimension(:),intent(in),optional :: fjac_sparse !! jacobian matrix [sparse]

    integer    :: i                  !! counter
    real(wp)   :: slope              !! local slope of the function of `alpha` along the search direction used for line search
    logical    :: min_alpha_reached  !! if the minimum step size is reached during the line search
    real(wp)   :: alpha              !! `alpha` for the line search
    real(wp)   :: ftmp               !! `f` value for linesearch
    real(wp)   :: t                  !! used for line search
    real(wp),dimension(:),allocatable :: gradf      !! line search objective function gradient vector
    real(wp),dimension(:),allocatable :: xtmp       !! `x` value for linesearch
    real(wp),dimension(:),allocatable :: fvectmp    !! `fvec` value for linesearch

    ! allocate arrays:
    allocate(gradf(me%n))
    allocate(xtmp(me%n))
    allocate(fvectmp(me%m))

    ! compute the gradient of the function to be minimized
    ! (which in this case is 1/2 the norm of fvec). Use the chain
    ! rule and the Jacobian matrix already computed.
    if (present(fjac)) then
        ! dense
        do i=1,me%n
            gradf(i) = dot_product(fvec,fjac(:,i))
        end do
    else
        ! sparse
        do i=1,me%n
            gradf(i) = dot_product(fvec,pack(fjac_sparse,mask=me%icol==i))
        end do
    end if
    slope = dot_product(p, gradf)
    t = -me%c * slope

    if (me%verbose) then
        write(me%iunit,'(1P,*(A,1X,E16.6))') '        slope    = ', slope
        write(me%iunit,'(1P,*(A,1X,E16.6))') '        t        = ', t
    end if

    ! perform the line search:
    min_alpha_reached = .false.
    alpha = me%alpha_max  ! start with the largest step
    do

        xtmp = xold + p * alpha
        call me%func(xtmp,fvectmp)
        ftmp = norm2(fvectmp)

        if (me%verbose) then
            write(me%iunit,'(1P,*(A,1X,E16.6))')          '        alpha    = ', alpha,    ' f       = ', ftmp
            if (f - ftmp >= alpha*t) then
                write(me%iunit,'(1P,2(A,1X,E16.6),1X,A)') '        f - ftmp = ', f - ftmp, ' alpha*t = ', alpha*t, ' [ACCEPTED]'
            else
                write(me%iunit,'(1P,*(A,1X,E16.6))')      '        f - ftmp = ', f - ftmp, ' alpha*t = ', alpha*t
            end if
        end if

        if (((f - ftmp) / 2.0_wp >= alpha*t) .or. min_alpha_reached) then
            if (me%verbose .and. min_alpha_reached) then
                write(me%iunit,'(A)') '        Minimum alpha reached'
            end if
            ! Armijo-Goldstein condition is satisfied
            ! (or the min step has been reached)
            x    = xtmp
            fvec = fvectmp
            f    = ftmp
            exit
        end if
        alpha = alpha * me%tau ! reduce step size

        if (alpha<=me%alpha_min) then
            alpha = me%alpha_min
            min_alpha_reached = .true. ! will stop on the next step
        end if

    end do

    end subroutine backtracking_linesearch