subroutine easysimplex(n,l,u,ml,meq,alin,acol,aval,annz,b,bdstype,
     +ivar,ctx,icase)

C     SCALAR ARGUMENTS
      character * 5 bdstype
      integer annz,icase,ivar,meq,ml,n
      double precision ctx

C     ARRAY ARGUMENTS
      integer alin(annz),acol(annz)
      double precision aval(annz),b(ml),l(n),u(n)

C     Input of this subroutine: 
C
C     Min c^t x 
C
C     subject to ai^t x  = bi, i = 1, ..., meq,
C                ai^t x >= bi, i = meq + 1, ..., ml,
C                l <= x <= u.
C
C     Intermediate result:
C
C     Min c^t y + c^t l
C
C     subject to ai^t y <= btmpi, i = 1, ..., m1,
C                ai^t y >= btmpi, i = m1 + 1, ..., m1 + m2
C                ai^t y  = btmpi, i = m1 + m2 + 1, ..., m1 + m2 + m3,
C                0 <= y <= utmp,
C
C     with y = x - l,
C          mltmp = m1 + m2 + m3, 
C          utmp  = u - l,
C          btmp  = b - A l,
C          and btmpi "redefined adequately" to have btmpi >= 0.

C     PARAMETERS
      include "dim.f"

C     LOCAL SCALARS
      integer i,j,m1,m2,m3,mltmp
      double precision ctl

C     LOCAL ARRAYS
      integer izrow(nmax),iposv(mlmax)
      double precision A(mlmax,nmax),btmp(mlmax),T(0:mlmax+1,0:nmax),
     +        utmp(nmax)

C     ==================================================================
C     Eliminate lower bound on x 
C     ==================================================================

C     Change of variables y = x - l.
C
C     New problem: 
C
C     Min c^t y + c^t l subject to A y >= btmp, 0 <= y <= utmp

C     btmp = b - A l

      do i = 1,ml
          btmp(i) = b(i)
      end do

      do i = 1,annz
          btmp(alin(i)) = btmp(alin(i)) - aval(i) * l(acol(i))
      end do

C     utmp = u - l

      do i = 1,n
          utmp(i) = u(i) - l(i)
C         Obviously, u(i) >= l(i). However, if they are very very 
C         similar, due to round-off errors, u(i) - l(i) might be
C         e very small negative number. As it will be used as the
C         RHS of the simples tableau, and a negative number in the
C         RHS cause a crash in the simplex subroutine being used,
C         we take the maximum with zero.
          utmp(i) = max( 0.0d0, utmp(i) )
      end do

C     ==================================================================
C     Transform A from sparse to dense format 
C     ==================================================================

      do j = 1,n
          do i = 1,ml
              A(i,j) = 0.0d0
          end do
      end do

      do i = 1,annz
         A(alin(i),acol(i)) = A(alin(i),acol(i)) + aval(i)
      end do

C     ==================================================================
C     SET TABLEAU 
C     ==================================================================

C     x(i) <= utmp(i)

      do j = 1,n
          do i = 1,n
              T(i,j) = 0.0d0
          end do
          T(j,j) = - 1.0d0
          T(j,0) = utmp(j)
      end do

      m1 = n

C     Inequality constraints ai^t x >= bi with negative b(i)

      do i = meq + 1,ml
          if ( btmp(i) .lt. 0.0d0 ) then
              m1 = m1 + 1
              do j = 1,n
                  T(m1,j) = A(i,j)
              end do
              T(m1,0) = - btmp(i)
          end if
      end do
 
C     Inequality constraints ai^t x >= bi with non-negative b(i)

      m2 = 0
      do i = meq + 1,ml
          if ( btmp(i) .ge. 0.0d0 ) then
              m2 = m2 + 1
              do j = 1,n
                  T(m1+m2,j) = - A(i,j)
              end do
              T(m1+m2,0) = btmp(i)
          end if
      end do
 
C     Equality constraints A x = b

      do i = 1,meq
          if ( btmp(i) .ge. 0.0d0 ) then
              do j = 1,n
                  T(m1+m2+i,j) = - A(i,j)
              end do
              T(m1+m2+i,0) =   btmp(i)
          else
              do j = 1,n
                  T(m1+m2+i,j) =   A(i,j)
              end do
              T(m1+m2+i,0) = - btmp(i)
          end if
      end do

      m3 = meq

C     Objective function value

      T(0,0) = 0.0d0

C     Objective function coefficients

      do j = 1,n
          T(0,j) = 0.0d0
      end do

      if ( bdstype .eq. 'LOWER' ) then
          T(0,ivar) = - 1.0d0
          ctl       = - l(ivar)
      else
          T(0,ivar) =   1.0d0
          ctl       =   l(ivar)
      end if

C     Number of constraints

      mltmp = m1 + m2 + m3

C     Call the Numerical Recipes Simplex subroutine (modified
C     to work in double precision)

      call simplx(T,mltmp,n,mlmax+2,nmax+1,m1,m2,m3,icase,izrow,iposv)

C     Objective function value: c^t y + c^t l

      ctx = T(0,0) + ctl

      end