cnmn05 Subroutine

    subroutine cnmn05(me, g, df, a, s, b, c, slope, phi, isc, ic, ms1, nvc, n1, n2, n3, n4, n5)

        !!  Routine to solve direction finding problem in modified method of
        !!  feasible directions.
        !!
        !!  BY G. N. VANDERPLAATS, MAY, 1972.
        !!
        !!  NORM OF S VECTOR USED HERE IS S-TRANSPOSE TIMES S<=1.
        !!  IF NVC = 0 FIND DIRECTION BY ZOUTENDIJK'S METHOD.  OTHERWISE
        !!  FIND MODIFIED DIRECTION.

        class(conmin_class), intent(inout) :: me
        integer, intent(in)       :: n1, n2, n3, n4, n5
        real(wp), intent(inout)   :: df(:), g(n2), a(n1, n3), s(:), c(n4), b(n3, n3)
        real(wp), intent(inout)   :: slope, phi
        integer, intent(inout)    :: isc(n2), ic(n3), ms1(n5), nvc

        real(wp)  :: a1, c1, ct1, ct2, cta, ctam, ctb, ctbm, ctc, ctd, gg, &
                     s1, sg, thmax, tht
        integer   :: i, j, j1, k, nac1, nci, ncj, ndb, ndv1, ndv2, ner

        ! ------------------------------------------------------------------
        ! ***  NORMALIZE GRADIENTS, CALCULATE THETA'S AND DETERMINE NVC  ***
        ! ------------------------------------------------------------------
        ndv1 = me%ndv + 1
        ndv2 = me%ndv + 2
        nac1 = me%nac + 1
        nvc = 0
        thmax = 0.0_wp
        cta = abs(me%ct)
        ct1 = 1.0_wp/cta
        ctam = abs(me%ctmin)
        ctb = abs(me%ctl)
        ct2 = 1.0_wp/ctb
        ctbm = abs(me%ctlmin)
        a1 = 1.0_wp
        do i = 1, me%nac
            ! CALCULATE THETA
            nci = ic(i)
            ncj = 1
            if (nci <= me%ncon) ncj = isc(nci)
            c1 = g(nci)
            ctd = ct1
            ctc = ctam
            if (ncj > 0) then
                ctc = ctbm
                ctd = ct2
            end if
            if (c1 > ctc) nvc = nvc + 1
            tht = 0.
            gg = 1.+ctd*c1
            if (ncj == 0 .or. c1 > ctc) tht = me%theta*gg*gg
            if (tht > 50.0_wp) tht = 50.0_wp
            if (tht > thmax) thmax = tht
            a(ndv1, i) = tht
            ! ------------------------------------------------------------------
            !                NORMALIZE GRADIENTS OF CONSTRAINTS
            ! ------------------------------------------------------------------
            a(ndv2, i) = 1.0_wp
            if (nci <= me%ncon) then
                a1 = 0.
                do j = 1, me%ndv
                    a1 = a1 + a(j, i)**2
                end do
                if (a1 < me%small) a1 = me%small
                a1 = sqrt(a1)
                a(ndv2, i) = a1
                a1 = 1.0_wp/a1
                do j = 1, me%ndv
                    a(j, i) = a1*a(j, i)
                end do
            end if
        end do
        ! ------------------------------------------------------------------
        ! CHECK FOR ZERO GRADIENT.  PROGRAM CHANGE-FEB, 1981, GV.
        ! ------------------------------------------------------------------
        i = 0
        main : do
            i = i + 1
            do
                if (a(ndv2, i) > 1.0e-6_wp) exit
                ! ZERO GRADIENT IS FOUND.  WRITE ERROR MESSAGE.
                if (me%iprint >= 2) write (me%iunit, 5000) ic(i)
                ! REDUCE NAC BY ONE.
                me%nac = me%nac - 1
                ! SHIFT COLUMNS OF A AND ROWS OF IC IF I<=NAC.
                if (i > me%nac) exit main
                ! SHIFT.
                do j = i, me%nac
                    j1 = j + 1
                    ic(j) = ic(j1)
                    do k = 1, ndv2
                        a(k, j) = a(k, j1)
                    end do
                end do
                if (i > me%nac) exit
            end do
            if (i >= me%nac) exit
        end do main
        if (me%nac <= 0) return
        nac1 = me%nac + 1
        ! DETERMINE IF CONSTRAINTS ARE VIOLATED.
        nvc = 0
        do i = 1, me%nac
            nci = ic(i)
            ncj = 1
            if (nci <= me%ncon) ncj = isc(nci)
            ctc = ctam
            if (ncj > 0) ctc = ctbm
            if (g(nci) > ctc) nvc = nvc + 1
        end do
        ! ------------------------------------------------------------------
        ! NORMALIZE GRADIENT OF OBJECTIVE FUNCTION AND STORE IN NAC+1
        ! COLUMN OF A
        ! ------------------------------------------------------------------
        a1 = 0.
        do i = 1, me%ndv
            a1 = a1 + df(i)**2
        end do
        if (a1 < me%small) a1 = me%small
        a1 = sqrt(a1)
        a1 = 1.0_wp/a1
        do i = 1, me%ndv
            a(i, nac1) = a1*df(i)
        end do
        ! BUILD C VECTOR.
        if (nvc <= 0) then
            ! ------------------------------------------------------------------
            !             BUILD C FOR CLASSICAL METHOD
            ! ------------------------------------------------------------------
            ndb = nac1
            a(ndv1, ndb) = 1.0_wp
            do i = 1, ndb
                c(i) = -a(ndv1, i)
            end do
        else
            ! ------------------------------------------------------------------
            !               BUILD C FOR MODIFIED METHOD
            ! ------------------------------------------------------------------
            ndb = me%nac
            a(ndv1, nac1) = -phi
            ! ------------------------------------------------------------------
            !       SCALE THETA'S SO THAT MAXIMUM THETA IS UNITY
            ! ------------------------------------------------------------------
            if (thmax > 0.00001_wp) thmax = 1./thmax
            do i = 1, ndb
                a(ndv1, i) = a(ndv1, i)*thmax
            end do
            do i = 1, ndb
                c(i) = 0.0_wp
                do j = 1, ndv1
                    c(i) = c(i) + a(j, i)*a(j, nac1)
                end do
            end do
        end if
        ! ------------------------------------------------------------------
        !                  BUILD B MATRIX
        ! ------------------------------------------------------------------
        do i = 1, ndb
            do j = 1, ndb
                b(i, j) = 0.0_wp
                do k = 1, ndv1
                    b(i, j) = b(i, j) - a(k, i)*a(k, j)
                end do
            end do
        end do
        ! ------------------------------------------------------------------
        !                SOLVE SPECIAL L. P. PROBLEM
        ! ------------------------------------------------------------------
        call cnmn08(ndb, ner, c, ms1, b, n3, n4)
        if (me%iprint > 1 .and. ner > 0) write (me%iunit, 5200)
        ! CALCULATE RESULTING DIRECTION VECTOR, S.
        slope = 0.0_wp
        ! ------------------------------------------------------------------
        !              USABLE-FEASIBLE DIRECTION
        ! ------------------------------------------------------------------
        do i = 1, me%ndv
            s1 = 0.0_wp
            if (nvc > 0) s1 = -a(i, nac1)
            do j = 1, ndb
                s1 = s1 - a(i, j)*c(j)
            end do
            slope = slope + s1*df(i)
            s(i) = s1
        end do
        s(ndv1) = 1.0_wp
        if (nvc > 0) s(ndv1) = -a(ndv1, nac1)
        do j = 1, ndb
            s(ndv1) = s(ndv1) - a(ndv1, j)*c(j)
        end do
        ! ------------------------------------------------------------------
        ! CHECK TO INSURE THE S-VECTOR IS FEASIBLE.
        ! PROGRAM MOD-FEB, 1981, GV.
        ! ------------------------------------------------------------------
        do j = 1, me%nac
            ! S DOT DEL(G).
            sg = dot_product(s(1:me%ndv), a(1:me%ndv, j))
            ! IF(SG>0.) GO TO 176
            ! THIS CHANGE MADE ON 4/8/81 FOR G. VANDERPLAATS
            if (sg > 1.0e-04_wp) then
                ! S-VECTOR IS NOT FEASIBLE DUE TO SOME NUMERICAL PROBLEM.
                if (me%iprint >= 2) write (me%iunit, 5100)
                s(ndv1) = 0.0_wp
                nvc = 0
                return
            end if
            ! FEASIBLE FOR THIS CONSTRAINT.  CONTINUE.
        end do

        ! ------------------------------------------------------------------
        !              NORMALIZE S TO MAX ABS OF UNITY
        ! ------------------------------------------------------------------
        s1 = 0.0_wp
        do i = 1, me%ndv
            a1 = abs(s(i))
            if (a1 > s1) s1 = a1
        end do
        ! IF (S1<1.0E-10) RETURN
        ! E-10 CHANGED TO E-04 ON 1/12/81

        if (s1 < 1.0e-04_wp) return
        s1 = 1.0_wp/s1
        do i = 1, me%ndv
            s(i) = s1*s(i)
        end do
        slope = s1*slope
        s(ndv1) = s1*s(ndv1)
        return

! ------------------------------------------------------------------
!                       FORMATS
! ------------------------------------------------------------------

5000    format(t6, '** CONSTRAINT', i5, ' HAS ZERO GRADIENT'/t6, &
               'DELETED FROM ACTIVE SET')
5100    format(t6, '** CALCULATED S-VECTOR IS NOT FEASIBLE'/t6, &
               'BETA IS SET TO ZERO')
5200    format(//t6, '* * DIRECTION FINDING PROCESS DID NOT CONVERGE'/t6, &
                '* * S-VECTOR MAY NOT BE VALID')

    end subroutine cnmn05