Results 1 to 9 of 9

Thread: Linear Programming DENSE_LP() producing wrong solution?

  1. #1
    Junior Member
    Join Date
    Mar 2006
    Posts
    11

    Exclamation Linear Programming DENSE_LP() producing wrong solution?

    The Linear Programming subroutine DENSE_LP() seems to be producing an incorrect solution. Its probably a programming problem I'm doing with the usual FORTRAN column-row issue as I know IMSL subrountines are thoroughly tested. I have been using IMSL libraries for 10+ years.
    The full program is shown below. The solution does not agree with the correct solution that Mathematic-a produces. You should know that both Octave and "R" also fail at getting to a solution, so this small system might be require something "special"? The solution should be {0, 0, 0, 0, 40, 0, 0, 40, 40, 0, 0, 0, 0, 0, 0, 0}, but I'm getting zeros....

    ================================================== ============
    PROGRAM linreg
    INCLUDE 'link_fnl_static.h'
    !DEC$ OBJCOMMENT LIB:"libiomp5md.lib"
    ! INCLUDE 'link_fnl_shared.h' ! (To use the dynamic link library form of the libraries)
    USE UMACH_INT
    USE WRRRN_INT
    USE DENSE_LP_INT
    IMPLICIT NONE
    INTEGER NOUT, M, NVAR, I, J
    PARAMETER (M=23, NVAR=16) ! for S8 23 rows (eqns) and 16 columns (column is an x)
    REAL*8 A(M, NVAR), B(M), C(NVAR), XSOL(NVAR), DSOL(M), BL(M), BU(M), OBJ, XLB(M), XUB(M)
    INTEGER IRTYPE(M)

    INTEGER ITREF, ITERS, IERR

    DATA A/1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., &
    0., 0., -1., 0., 0., 1., 0., 1., -1., -1., 1., 0., -1., 0., -1., 0., &
    -2., 0., 1., 0., 0., -1., 0., -1., 1., 1., -1., -2., 1., 0., 0., 0., &
    -1., 0., 0., 0., 1., 0., 1., -1., 0., 0., 0., -1., 0., 0., 0., 0., &
    0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., &
    -1., -1., 0., 0., 0., 0., 0., 0., 0., 0., 0., -1., 0., 0., 0., -1., &
    0., 0., 2., 0., 0., 0., 0., 0., 0., 1., -1., -1., 0., 1., 0., 0., &
    0., 0., -1., 0., 0., 0., 0., 0., 0., -1., 1., 1., 0., 0., 0., 0., &
    -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1./

    DATA C/-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1./
    DATA B/0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0. ,0.,0.,0.,0.,0.,-120./

    ! solution is {0, 0, 0, 0, 40, 0, 0, 40, 40, 0, 0, 0, 0, 0, 0, 0}

    BL=B
    BU=-1.0d30
    IRTYPE=3 ! type cof constraints
    XLB=1.0d30
    XUB=-1.0d30
    ITERS=0; IERR=0
    print *,'m,nvar=',m,nvar

    print *,'1,B='
    print 3,(b(j), j=1,m)
    print *,'2.C='
    print 1,(c(j), j=1,nvar)
    print *,'BOUNDS'
    print *,'3. bl='
    print 3,(bl(j), j=1,m)
    print *,'4. bu='
    print 4,(bu(j), j=1,m)
    print *,'>Before matrix a (bottom matrix of S8):'
    do i=1,m
    print 1,(a(i,j), j=1,nvar)
    end do
    ITREF=2
    CALL UMACH(2, NOUT)
    ! Solve the LP problem
    CALL DENSE_LP (A, BL, BU, C, IRTYPE, OBJ, XSOL, DSOL) ! pagr 1346 Math
    !CALL DENSE_LP (A, BL, BU, C, IRTYPE, OBJ, XSOL, DSOL,m, nvar,m,XLB, XUB,ITREF, ITERS, IERR)

    99999 FORMAT (' Objective', F9.4)
    WRITE(NOUT, 99999) OBJ
    CALL WRRRN('Solution', XSOL, 1, NVAR, 1)

    print *,'*********************************************** *******************************'

    ! XSOL(NVAR), DSOL(M)
    print *,'XSOL'
    print 3,(XSOL(j), j=1,NVAR)
    print *,'DSOL'
    print 3,(DSOL(j), j=1,m)

    print *,'iterations to solution=',ITERS
    print *,'IERR=',IERR

    print *,'2 B='
    print 3,(b(j), j=1,m)
    print *,'2 C='
    print 1,(c(j), j=1,nvar)
    ! print *,'**************************************'
    ! print *,'MATMUL=',MATMUL(DSOL,A)
    1 format(1x, 60(F5.1))
    2 format(1x, 60(I3))
    3 format(1x, 60(F6.1))
    4 format(1x, 60(E8.1))
    END
    ================================================== =========================
    Compiled with full debugging and "-check all -g -traceback" with Intel compilers.
    Output is:
    m,nvar= 23 16
    1,B=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    2.C=
    -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0
    BOUNDS
    3. bl=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    4. bu=
    -0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31
    >Before matrix a (bottom matrix of S8):
    1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 2.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 -2.0 -1.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0
    0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 -1.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 1.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 -1.0 0.0 0.0 -1.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 -1.0
    0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 -1.0
    0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 1.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 1.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 -2.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 1.0 1.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 -1.0 -1.0
    1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 -1.0 0.0 0.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 1.0 0.0 0.0 -1.0 0.0 -1.0
    0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0
    0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 0.0 -1.0
    Objective 0.0000

    Solution
    1 2 3 4 5 6 7 8
    0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000

    9 10 11 12 13 14 15 16
    0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
    ************************************************** ****************************
    XSOL
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    DSOL
    1.4 1.0 1.0 1.4 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.4 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0
    iterations to solution= 0
    IERR= 0
    2 B=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    Last edited by jianni; 06-30-2020 at 12:50 AM.

  2. #2
    Senior Member mecej4's Avatar
    Join Date
    Dec 2009
    Posts
    139
    I suspect that there is at least one error, from looking at your DATA statement for A.

    In Fortran 2-D arrays are column-major, so the first line of that DATA statement causes the 16 values to be entered into rows 1-16 or column-1 of A. From the 16 values in the next line of the source 7 will by used to complete the remaining part of column-1, and the remaining 9 will be placed in rows 1-9 of column 2. And so forth.

    However, since you did not give a mathematical description of the LP problem, I am only guessing. In the DATA statement you have formatted the numbers in standard matrix format. In the output you have displayed, after the label "Before matrix a (bottom matrix of S8):", the transpose of that matrix. Which one is the one desired, or is the actual A something else?

  3. #3
    Junior Member
    Join Date
    Mar 2006
    Posts
    11
    Unfortunately I get the same incorrect result when the DATA lines are reformatted in row-major like this:

    DATA (A(1,n),n=1,NVAR)/1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(2,n),n=1,NVAR)/0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(3,n),n=1,NVAR)/0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(4,n),n=1,NVAR)/ 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(5,n),n=1,NVAR)/ 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(6,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(7,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(8,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(9,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(10,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0./
    DATA (A(11,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0./
    DATA (A(12,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0./
    DATA (A(13,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0./
    DATA (A(14,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0./
    DATA (A(15,n),n=1,NVAR)/0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0./
    DATA (A(16,n),n=1,NVAR)/0., 0., -1., 0., 0., 1., 0., 1., -1., -1., 1., 0., -1., 0., -1., 0./
    DATA (A(17,n),n=1,NVAR)/ -2., 0., 1., 0., 0., -1., 0., -1., 1., 1., -1., -2., 1., 0., 0., 0./
    DATA (A(18,n),n=1,NVAR)/-1., 0., 0., 0., 1., 0., 1., -1., 0., 0., 0., -1., 0., 0., 0., 0./
    DATA (A(19,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1./
    DATA (A(20,n),n=1,NVAR)/-1., -1., 0., 0., 0., 0., 0., 0., 0., 0., 0., -1., 0., 0., 0., -1./
    DATA (A(21,n),n=1,NVAR)/ 0., 0., 2., 0., 0., 0., 0., 0., 0., 1., -1., -1., 0., 1., 0., 0./
    DATA (A(22,n),n=1,NVAR)/ 0., 0., -1., 0., 0., 0., 0., 0., 0., -1., 1., 1., 0., 0., 0., 0./
    DATA (A(23,n),n=1,NVAR)/ -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1./

    OUTPUT
    m,nvar= 23 16
    1,B=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    2.C=
    -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0
    BOUNDS
    3. bl=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    4. bu=
    -0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31-0.1E+31
    >Before matrix a (bottom matrix of S8):
    1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0
    0.0 0.0 -1.0 0.0 0.0 1.0 0.0 1.0 -1.0 -1.0 1.0 0.0 -1.0 0.0 -1.0 0.0
    -2.0 0.0 1.0 0.0 0.0 -1.0 0.0 -1.0 1.0 1.0 -1.0 -2.0 1.0 0.0 0.0 0.0
    -1.0 0.0 0.0 0.0 1.0 0.0 1.0 -1.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0
    -1.0 -1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 -1.0
    0.0 0.0 2.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 -1.0 -1.0 0.0 1.0 0.0 0.0
    0.0 0.0 -1.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 1.0 1.0 0.0 0.0 0.0 0.0
    -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0
    Objective 0.0000

    Solution
    1 2 3 4 5 6 7 8
    0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000

    9 10 11 12 13 14 15 16
    0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
    ************************************************** ****************************
    XSOL
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    DSOL
    0.7 1.0 0.7 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 0.7 1.0 1.0 1.0 1.0 1.4 1.0 1.0 1.0 1.4 1.0 1.0
    iterations to solution= 0
    IERR= 0
    2 B=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    2 C=
    -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0

  4. #4
    Junior Member
    Join Date
    Mar 2006
    Posts
    11
    Mathomoitica solves it correct:

    c = ( - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 )


    Out[3]= {{- 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1}}

    In[14]:= Dimensions[c]

    Out[14]= {1, 16}

    In[2]:= b = ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 - 120 )

    Out[2]= {{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, - 120}}



    In[16]:= MatrixForm[%15]

    Out[16]//MatrixForm=
    1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
    0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
    0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
    0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
    0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
    0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0
    0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
    0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0
    0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
    0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
    0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0
    0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
    0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
    0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
    0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
    0 0 - 1 0 0 1 0 1 - 1 - 1 1 0 - 1 0 - 1 0
    - 2 0 1 0 0 - 1 0 - 1 1 1 - 1 - 2 1 0 0 0
    - 1 0 0 0 1 0 1 - 1 0 0 0 - 1 0 0 0 0
    0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1
    - 1 - 1 0 0 0 0 0 0 0 0 0 - 1 0 0 0 - 1
    0 0 2 0 0 0 0 0 0 1 - 1 - 1 0 1 0 0
    0 0 - 1 0 0 0 0 0 0 - 1 1 1 0 0 0 0
    - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1

    In[17]:= LinearProgramming[c[[1]], A, b[[1]]]


    Out[17]= {0, 0, 0, 0, 40, 0, 0, 40, 40, 0, 0, 0, 0, 0, 0, 0}

  5. #5
    Senior Member mecej4's Avatar
    Join Date
    Dec 2009
    Posts
    139
    I'd like to see a standard mathematical problem statement before responding further. Different software packages have different conventions, and the data may need conversion/reformatting to suit.

  6. #6
    Junior Member
    Join Date
    Mar 2006
    Posts
    11
    Please see the free supporting information PDF at
    https://pubs.acs.org/doi/suppl/10.10...771_si_001.pdf
    The LP problem is on Page S8. That's what I'm trying to solve. Please let me knowif you need any further details.

  7. #7
    Senior Member mecej4's Avatar
    Join Date
    Dec 2009
    Posts
    139
    The main point is that these rather special LP problems may not give you a unique solution.

    A source of errors in comparing IMSL results and other results, e.g., from Mathematica, is that the mathematical conventions are not consistent.

    The Mathematica problem (for the case that you used, with LinearProgramming(c,A,b) is:

    Find a vector x that minimizes the quantity c^T.x subject to the constraints A.x>=b and x>=0.

    In contrast, IMSL DENSE_LP minimizes c^T.x subject to and b_l <= A.x <= b_u and x_l <= x <= x_u

    After you fix your code to accommodate these different conventions and run with IMSL DENSE_LP, please add the optional argument IERR, and note the value returned for it.

    You will probably find that IERR = 1, which signifies that "Multiple solutions giving essentially the same solution exist". For the purposes of the ACS paper that you referred to, that does not matter, since they only wish to find whether c^T.x is equal to zero. In such circumstances, it should not matter that different LP solvers give different solutions as long as c^T.x as computed from those solutions comes out the same .

  8. #8
    Junior Member
    Join Date
    Mar 2006
    Posts
    11
    I'm getting one of the "solutions" with all zeros. Program and output below:
    PROGRAM linreg
    INCLUDE 'link_fnl_static.h'
    !DEC$ OBJCOMMENT LIB:"libiomp5md.lib"
    ! INCLUDE 'link_fnl_shared.h' ! (To use the dynamic link library form of the libraries)
    USE UMACH_INT
    USE WRRRN_INT
    USE DENSE_LP_INT
    IMPLICIT NONE
    INTEGER NOUT, M, NVAR, I, J, n

    PARAMETER (M=23, NVAR=16) ! for S8 23 rows (eqns) and 16 columns (column is an x)
    REAL*8 A(M, NVAR), B(M), C(NVAR), XSOL(NVAR), DSOL(M), BL(M), BU(M), OBJ, XLB(M), XUB(M)
    INTEGER IRTYPE(M)

    INTEGER ITREF, ITERS, IERR

    DATA (A(1,n),n=1,NVAR)/1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(2,n),n=1,NVAR)/0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(3,n),n=1,NVAR)/0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(4,n),n=1,NVAR)/ 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(5,n),n=1,NVAR)/ 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(6,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(7,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(8,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(9,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0., 0./
    DATA (A(10,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0., 0./
    DATA (A(11,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0., 0./
    DATA (A(12,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0., 0./
    DATA (A(13,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 0./
    DATA (A(14,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0./
    DATA (A(15,n),n=1,NVAR)/0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0./
    DATA (A(16,n),n=1,NVAR)/0., 0., -1., 0., 0., 1., 0., 1., -1., -1., 1., 0., -1., 0., -1., 0./
    DATA (A(17,n),n=1,NVAR)/ -2., 0., 1., 0., 0., -1., 0., -1., 1., 1., -1., -2., 1., 0., 0., 0./
    DATA (A(18,n),n=1,NVAR)/-1., 0., 0., 0., 1., 0., 1., -1., 0., 0., 0., -1., 0., 0., 0., 0./
    DATA (A(19,n),n=1,NVAR)/ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1./
    DATA (A(20,n),n=1,NVAR)/-1., -1., 0., 0., 0., 0., 0., 0., 0., 0., 0., -1., 0., 0., 0., -1./
    DATA (A(21,n),n=1,NVAR)/ 0., 0., 2., 0., 0., 0., 0., 0., 0., 1., -1., -1., 0., 1., 0., 0./
    DATA (A(22,n),n=1,NVAR)/ 0., 0., -1., 0., 0., 0., 0., 0., 0., -1., 1., 1., 0., 0., 0., 0./
    DATA (A(23,n),n=1,NVAR)/ -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1., -1./


    DATA C/-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1.,-1./
    DATA B/0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0. ,0.,0.,0.,0.,0.,-120./

    BL=0.d0 ! lower constraints
    ! BL=B
    BU=10000.0d0
    IRTYPE=2 ! type cof constraints
    XLB=0.0d0
    XUB=100.0d0
    ITERS=0; IERR=0
    print *,'m,nvar=',m,nvar

    print *,'B='
    print 3,(b(j), j=1,m)
    print *,'C='
    print 1,(c(j), j=1,nvar)
    print *,'BOUNDS:'
    print *,' bl='
    print 3,(bl(j), j=1,m)
    print *,' bu='
    print 4,(bu(j), j=1,m)
    print *,'XLB='
    print 4,(XLB(j), j=1,nvar)
    print *,'Xub='
    print 4,(XUB(j), j=1,nvar)

    print *,'>Before matrix a (bottom matrix of S8):'
    do i=1,m
    print 1,(a(i,j), j=1,nvar)
    end do
    ITREF=2 ! refinement 1 iterative, 2=extended iterative, 0=none
    CALL UMACH(2, NOUT)
    ! Solve the LP problem
    !CALL DENSE_LP (A, BL, BU, C, IRTYPE, OBJ, XSOL, DSOL) ! pagr 1346 Math
    CALL DENSE_LP (A, BL, BU, C, IRTYPE, OBJ, XSOL, DSOL,m, nvar,m,XLB, XUB,ITREF, ITERS, IERR)

    99999 FORMAT (' Objective', F9.4)
    WRITE(NOUT, 99999) OBJ
    CALL WRRRN('Solution', XSOL, 1, NVAR, 1)

    print *,'*********************************************** *******************************'

    print *,'XSOL'
    print 3,(XSOL(j), j=1,NVAR)
    print *,'DSOL'
    print 3,(DSOL(j), j=1,m)

    print *,'iterations to solution=',ITERS
    print *,'IERR=',IERR

    1 format(1x, 60(F5.1))
    2 format(1x, 60(I3))
    3 format(1x, 60(F6.1))
    4 format(1x, 60(E8.1))
    END


    ================================================== ===============================================
    OUTPUT:
    m,nvar= 23 16
    B=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0-120.0
    C=
    -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0
    BOUNDS:
    bl=
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    bu=
    0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05 0.1E+05
    XLB=
    0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00 0.0E+00
    Xub=
    0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03 0.1E+03
    >Before matrix a (bottom matrix of S8):
    1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0
    0.0 0.0 -1.0 0.0 0.0 1.0 0.0 1.0 -1.0 -1.0 1.0 0.0 -1.0 0.0 -1.0 0.0
    -2.0 0.0 1.0 0.0 0.0 -1.0 0.0 -1.0 1.0 1.0 -1.0 -2.0 1.0 0.0 0.0 0.0
    -1.0 0.0 0.0 0.0 1.0 0.0 1.0 -1.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 0.0
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0
    -1.0 -1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 0.0 0.0 0.0 -1.0
    0.0 0.0 2.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 -1.0 -1.0 0.0 1.0 0.0 0.0
    0.0 0.0 -1.0 0.0 0.0 0.0 0.0 0.0 0.0 -1.0 1.0 1.0 0.0 0.0 0.0 0.0
    -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0 -1.0
    Objective 0.0000

    Solution
    1 2 3 4 5 6 7 8
    0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000

    9 10 11 12 13 14 15 16
    0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
    ************************************************** ****************************
    XSOL
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
    DSOL
    0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0
    iterations to solution= 3
    IERR= 1

  9. #9
    Senior Member mecej4's Avatar
    Join Date
    Dec 2009
    Posts
    139
    I note that you did obtain IERR = 1 from DENSE_LP. The implication is that the solution is not unique, and the solution found was x_sol = 0.

    Since there are other solutions with some components of x_sol not equal to 0, such as the Mathematica solver's solution, the method of the ACS paper appears to have a weakness. That is, depending on the LP solver used, the solution may or may not have some components of x_sol non-zero. The IMSL solver DENSE_LP would lead to the conclusion "no illegal loops", whereas use of the Mathematica solver would lead to a different conclusion.

    I am afraid that I do not know the subject matter of the ACS paper. Perhaps you can communicate these findings to the authors and request a response.

Tags for this Thread

Posting Permissions

  • You may not post new threads
  • You may not post replies
  • You may not post attachments
  • You may not edit your posts
  •