cnmn02 Subroutine

    subroutine cnmn02(me, ncalc, slope, dftdf1, df, s)

        !!  Routine to determine conjugate direction vector or direction
        !!  of steepest descent for unconstrained function minimization.
        !!
        !!  BY G. N. VANDERPLAATS, APRIL, 1972.
        !!
        !!  Conjugate direction is found by fletcher-reeves algorithm.

        class(conmin_class), intent(inout) :: me
        integer, intent(inout)     :: ncalc !! NCALC = CALCULATION CONTROL.
                                            !!
                                            !!  * NCALC = 0,     S = STEEPEST DESCENT.
                                            !!  * NCALC = 1,     S = CONJUGATE DIRECTION.
        real(wp), intent(out)     :: slope
        real(wp), intent(inout)   :: dftdf1
        real(wp), intent(in)      :: df(:)
        real(wp), intent(inout)   :: s(:)

        integer  :: i
        real(wp) :: beta, dfi, dftdf, s1, s2, si
        logical  :: fletcher_reeves !! if the fletcher-reeves conjugate direction was computed

        dftdf = sum(df(1:me%ndv)**2) ! CALCULATE NORM OF GRADIENT VECTOR

        ! FIND DIRECTION S
        fletcher_reeves = .false.
        if (ncalc == 1) then
            if (dftdf1 >= me%small) then
                ! FIND FLETCHER-REEVES CONJUGATE DIRECTION
                beta = dftdf/dftdf1
                slope = 0.0_wp
                do i = 1, me%ndv
                    dfi = df(i)
                    si = beta*s(i) - dfi
                    slope = slope + si*dfi
                    s(i) = si
                end do
                fletcher_reeves = .true.
            end if
        end if
        if (.not. fletcher_reeves) then
            ncalc = 0
            ! CALCULATE DIRECTION OF STEEPEST DESCENT
            s(1:me%ndv) = -df(1:me%ndv)
            slope = -dftdf
        end if
        ! NORMALIZE S TO MAX ABS VALUE OF UNITY
        s1 = 0.0_wp
        do i = 1, me%ndv
            s2 = abs(s(i))
            if (s2 > s1) s1 = s2
        end do
        if (s1 < me%small) s1 = me%small
        s1 = 1.0_wp/s1
        dftdf1 = dftdf*s1
        s(1:me%ndv) = s1*s(1:me%ndv)
        slope = s1*slope

    end subroutine cnmn02