cnmn08 Subroutine

    subroutine cnmn08(ndb, ner, c, ms1, b, n3, n4)

        !!  Routine to solve special linear problem for imposing s-transpose
        !!  times s<=1 bounds in the modified method of feasible directions.
        !!
        !!  FORM OF L. P. IS BX=C WHERE 1ST NDB COMPONENTS OF X CONTAIN VECTOR
        !!  U AND LAST NDB COMPONENTS CONTAIN VECTOR V.
        !!  CONSTRAINTS ARE U>=0, V>=0, AND U-TRANSPOSE TIMES V = 0.
        !!
        !!  BY G. N. VANDERPLAATS, APRIL, 1972.
        !!
        !!### Reference
        !!  * "[Structural optimization by methods of feasible directions](https://www.sciencedirect.com/science/article/abs/pii/0045794973900552)",
        !!    G. N. Vanderplaats and F. Moses, Journal of computers
        !!    and structures, vol 3, pp 739-755, 1973.

        !  ------------------------------------------------------------------
        !  CHOOSE INITIAL BASIC VARIABLES AS V, AND INITIALIZE VECTOR MS1
        !  ------------------------------------------------------------------

        integer, intent(in)        :: ndb, n3, n4
        integer, intent(out)       :: ner !! NER = ERROR FLAG.  IF NER/=0 ON RETURN, PROCESS HAS NOT
                                          !! CONVERGED IN 5*NDB ITERATIONS.
        integer, intent(out)       :: ms1(:) !! VECTOR MS1 IDENTIFIES THE SET OF BASIC VARIABLES.
        real(wp), intent(inout)   :: c(n4), b(n3, n3)

        integer   :: i, ichk, iter1, j, jj, kk, m2, nmax
        real(wp)  :: bb, bb1, bi, c1, cb, cbmax, cbmin, eps

        ner = 1
        m2 = 2*ndb
        ! CALCULATE CBMIN AND EPS AND INITIALIZE MS1.
        eps = -1.0e+10_wp
        cbmin = 0.0_wp
        do i = 1, ndb
            bi = b(i, i)
            cbmax = 0.
            if (bi < -1.0e-6_wp) cbmax = c(i)/bi
            if (bi > eps) eps = bi
            if (cbmax > cbmin) cbmin = cbmax
            ms1(i) = 0
        end do
        eps = 0.0001_wp*eps
!       IF (EPS<-1.0E-10) EPS=-1.0E-10
!
!  E-10 CHANGED TO E-03 ON 1/12/81
!
        if (eps < -1.0e-03_wp) eps = -1.0e-03_wp
        if (eps > -0.0001_wp) eps = -0.0001_wp
        cbmin = cbmin*1.0e-6_wp
!       IF (CBMIN<1.0e-10_wp) CBMIN=1.0e-10_wp
!
!  E-10 CHANGED TO E-05 ON 1/12/81
!
        if (cbmin < 1.0e-05_wp) cbmin = 1.0e-05_wp
        iter1 = 0
        nmax = 5*ndb
        main: do
            ! ------------------------------------------------------------------
            ! **********             BEGIN NEW ITERATION              **********
            ! ------------------------------------------------------------------
            iter1 = iter1 + 1
            if (iter1 > nmax) return
            ! FIND MAX. C(I)/B(I,I) FOR I=1,NDB.
            cbmax = 0.9_wp*cbmin
            ichk = 0
            do i = 1, ndb
                c1 = c(i)
                bi = b(i, i)
                ! IF (BI>EPS .OR. C1>0.) GO TO 30
                if (bi <= eps .and. c1 <= -1.0e-05_wp) then
                    ! 0. CHANGED TO -1.0E-05 ON 1/12/81
                    cb = c1/bi
                    if (cb > cbmax) then
                        ichk = i
                        cbmax = cb
                    end if
                end if
            end do
            if (cbmax >= cbmin) then
                if (ichk /= 0) then
                    ! UPDATE VECTOR MS1.
                    jj = ichk
                    if (ms1(jj) == 0) jj = ichk + ndb
                    kk = jj + ndb
                    if (kk > m2) kk = jj - ndb
                    ms1(kk) = ichk
                    ms1(jj) = 0
                    ! ------------------------------------------------------------------
                    !                 PIVOT OF B(ICHK,ICHK)
                    ! ------------------------------------------------------------------
                    bb = 1.0_wp/b(ichk, ichk)
                    do j = 1, ndb
                        b(ichk, j) = bb*b(ichk, j)
                    end do
                    c(ichk) = cbmax
                    b(ichk, ichk) = bb
                    ! ELIMINATE COEFICIENTS ON VARIABLE ENTERING BASIS AND STORE
                    ! COEFICIENTS ON VARIABLE LEAVING BASIS IN THEIR PLACE.
                    do i = 1, ndb
                        if (i /= ichk) then
                            bb1 = b(i, ichk)
                            b(i, ichk) = 0.0_wp
                            do j = 1, ndb
                                b(i, j) = b(i, j) - bb1*b(ichk, j)
                            end do
                            c(i) = c(i) - bb1*cbmax
                        end if
                    end do
                    cycle main
                end if
            end if
            exit main
        end do main
        ner = 0
        ! ------------------------------------------------------------------
        ! STORE ONLY COMPONENTS OF U-VECTOR IN 'C'.  USE B(I,1) FOR
        ! TEMPORARY STORAGE
        ! ------------------------------------------------------------------
        b(1:ndb, 1) = c(1:ndb)
        do i = 1, ndb
            c(i) = 0.0_wp
            j = ms1(i)
            if (j > 0) c(i) = b(j, 1)
            if (c(i) < 0.0_wp) c(i) = 0.0_wp
        end do

    end subroutine cnmn08