R version 2.9.0 (2009-04-17) Copyright (C) 2009 The R Foundation for Statistical Computing ISBN 3-900051-07-0 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > # smplx.r simplex algorithm > # J F Monahan, July 2009 > # > # first the code > # > # phase 1 / phase 2 algorithm > p1p2 <- function(A,b,c,m,n) { + p <- n+(1:m) # indices of starting basic feasible solution + print(dim(A)) + Anew <- cbind(A*sign(b),diag(rep(1,m))) # append identity + print(dim(Anew)) + cnew <- c(rep(0,n),rep(1,m)) # new cost vector + feas <- simplex(p,Anew,abs(b),cnew,m,n+m) # solve lp problem + print("******donewithphase1****************") + print(feas) + print("phase1resultsabove") + if( (feas[[3]] == 0)*(feas[[4]] == 0 ) ) { # found feasible soln + p <- feas[[1]] + p1p2 <- simplex(p,A,b,c,m,n) } + else { # no feasible soln + print("nofeasiblesolution") + p1p2 <- list(feas[[1]],feas[[2]],3,0) } + } # end of function p1p2 > # > # simplex algorithm > simplex <- function(p,A,b,c,m,n) { # dumb simplex algorithm + print(dim(A)) + print(p) + B <- A[,p] # matrix of starting solution + BiAb <- solve(B,cbind(A,b)) + # start loop + for( i in (1:100)) { # begin loop + print(p) + pn <- (1:n)[-p] + print(pn) + cp <- c[p] + cn <- c[pn] + x <- rep(0,n) + x[p] <- BiAb[,n+1] + print(BiAb[,n+1]) + print(sum(cp*BiAb[,n+1])) # current obj fn value + # find feasible directions + zjmcj <- cp %*% BiAb[,pn] - cn + print(zjmcj) + flist <- (1:(n-m))*(zjmcj > 1.e-10) + print(flist) + #p print((BiAb[,flist]>1.e-10)) + #p print(apply((BiAb[,flist]> 1.e-10)+0,2,sum)) + if( sum(flist) == 0 ) { print("done") + code <- 0 + break } + ibig <- which.max(zjmcj) # ibig is in, what is out? + print(c(ibig,max(zjmcj))) + ibig <- pn[ibig] + print(ibig) + if( max(BiAb[,ibig]) < 1.e-10 ) { + print("unbound") + code <- 2 + break } + v <- BiAb[,n+1]/BiAb[,ibig] + v[ (BiAb[,ibig]<1.e-10) ] <- 1.e30 + lmin <- which.min( v ) # lmin is out + print(c(lmin,v[lmin])) + v <- BiAb[,ibig]/BiAb[lmin,ibig] + v[lmin] <- (BiAb[lmin,ibig]-1)/BiAb[lmin,ibig] + BiAb <- BiAb - outer(v,BiAb[lmin,]) + print(BiAb) + p[lmin] <- ibig + code <- 1 # too many iterations if leave with this + } # end of loop + x <- rep(0,n) + x[p] <- BiAb[,n+1] + simplex <- list(p,x,code,sum(c*x)) } # done with simplex > # call phase 1/ phase2 algorithm > # > # first problem from Wagner, p103, 119-120 > m <- 3 > n <- 7 > alldat <- scan("smplx11.dat") > Ab <- matrix(alldat,m,n+1,byrow=TRUE) > A <- Ab[,1:n] > b <- Ab[,n+1] > c <- alldat[m*(n+1)+(1:n)] > A [,1] [,2] [,3] [,4] [,5] [,6] [,7] [1,] 1 1 1 1 1 0 0 [2,] 7 5 3 2 0 1 0 [3,] 3 5 10 15 0 0 1 > b [1] 15 120 100 > c [1] -4 -5 -9 -11 0 0 0 > prob1 <- p1p2(A,b,c,m,n) [1] 3 7 [1] 3 10 [1] 3 10 [1] 8 9 10 [1] 8 9 10 [1] 1 2 3 4 5 6 7 [1] 15 120 100 [1] 235 [1,] 11 11 14 18 1 1 1 [1,] 1 2 3 4 5 6 7 [1] 4 18 [1] 4 [1] 3.000000 6.666667 b [1,] 0.8 0.6666667 0.3333333 0 1 0 -0.06666667 1 0 -0.06666667 8.333333 [2,] 6.6 4.3333333 1.6666667 0 0 1 -0.13333333 0 1 -0.13333333 106.666667 [3,] 0.2 0.3333333 0.6666667 1 0 0 0.06666667 0 0 0.06666667 6.666667 [1] 8 9 4 [1] 1 2 3 5 6 7 10 [1] 8.333333 106.666667 6.666667 [1] 115 [1,] 7.4 5 2 1 1 -0.2 -1.2 [1,] 1 2 3 4 5 0 0 [1] 1.0 7.4 [1] 1 [1] 1.00000 10.41667 [1,] 1.000000e+00 0.8333333 0.4166667 0 1.25 0 -0.08333333 1.25 0 [2,] 8.881784e-16 -1.1666667 -1.0833333 0 -8.25 1 0.41666667 -8.25 1 [3,] 0.000000e+00 0.1666667 0.5833333 1 -0.25 0 0.08333333 -0.25 0 b [1,] -0.08333333 10.416667 [2,] 0.41666667 37.916667 [3,] 0.08333333 4.583333 [1] 1 9 4 [1] 2 3 5 6 7 8 10 [1] 10.416667 37.916667 4.583333 [1] 37.91667 [1,] -1.166667 -1.083333 -8.25 1 0.4166667 -9.25 -0.5833333 [1,] 0 0 0 4 5 0 0 [1] 4 1 [1] 6 [1] 2.00000 37.91667 [1,] 1.000000e+00 0.8333333 0.4166667 0 1.25 0 -0.08333333 1.25 0 [2,] 8.881784e-16 -1.1666667 -1.0833333 0 -8.25 1 0.41666667 -8.25 1 [3,] 0.000000e+00 0.1666667 0.5833333 1 -0.25 0 0.08333333 -0.25 0 b [1,] -0.08333333 10.416667 [2,] 0.41666667 37.916667 [3,] 0.08333333 4.583333 [1] 1 6 4 [1] 2 3 5 7 8 9 10 [1] 10.416667 37.916667 4.583333 [1] 0 [1,] 0 0 0 0 -1 -1 -1 [1,] 0 0 0 0 0 0 0 [1] "done" [1] "******donewithphase1****************" [[1]] [1] 1 6 4 [[2]] [1] 10.416667 0.000000 0.000000 4.583333 0.000000 37.916667 0.000000 [8] 0.000000 0.000000 0.000000 [[3]] [1] 0 [[4]] [1] 0 [1] "phase1resultsabove" [1] 3 7 [1] 1 6 4 [1] 1 6 4 [1] 2 3 5 7 [1] 10.416667 37.916667 4.583333 [1] -92.08333 [1,] -0.1666667 0.9166667 -2.25 -0.5833333 [1,] 0 2 0 0 [1] 2.0000000 0.9166667 [1] 3 [1] 3.000000 7.857143 b [1,] 1 0.7142857 0 -0.7142857 1.4285714 0 -0.1428571 7.142857 [2,] 0 -0.8571429 0 1.8571429 -8.7142857 1 0.5714286 46.428571 [3,] 0 0.2857143 1 1.7142857 -0.4285714 0 0.1428571 7.857143 [1] 1 6 3 [1] 2 4 5 7 [1] 7.142857 46.428571 7.857143 [1] -99.28571 [1,] -0.4285714 -1.571429 -1.857143 -0.7142857 [1,] 0 0 0 0 [1] "done" > prob1 [[1]] [1] 1 6 3 [[2]] [1] 7.142857 0.000000 7.857143 0.000000 0.000000 46.428571 0.000000 [[3]] [1] 0 [[4]] [1] -99.28571 > # > # second problem from VJ Bowman > m <- 6 > n <- 10 > alldat <- scan("smplxy.dat") > Ab <- matrix(alldat,m,n+1,byrow=TRUE) > A <- Ab[,1:n] > b <- Ab[,n+1] > c <- alldat[m*(n+1)+(1:n)] > A [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] 1 4 5 2 1 0 0 0 0 0 [2,] 2 3 0 0 0 1 0 0 0 0 [3,] 5 1 0 0 0 0 1 0 0 0 [4,] 0 0 1 0 0 0 0 1 0 0 [5,] 0 0 0 1 0 0 0 0 1 0 [6,] 0 0 -3 -4 0 0 0 0 0 1 > b [1] 7 6 5 4 3 -12 > c [1] -1 -8 -5 -6 0 0 0 0 0 0 > prob2 <- p1p2(A,b,c,m,n) [1] 6 10 [1] 6 16 [1] 6 16 [1] 11 12 13 14 15 16 [1] 11 12 13 14 15 16 [1] 1 2 3 4 5 6 7 8 9 10 [1] 7 6 5 4 3 12 [1] 37 [1,] 8 8 9 7 1 1 1 1 1 -1 [1,] 1 2 3 4 5 6 7 8 9 0 [1] 3 9 [1] 3 [1] 1.0 1.4 b [1,] 0.2 0.8 1 0.4 0.2 0 0 0 0 0 0.2 0 0 0 0 0 1.4 [2,] 2.0 3.0 0 0.0 0.0 1 0 0 0 0 0.0 1 0 0 0 0 6.0 [3,] 5.0 1.0 0 0.0 0.0 0 1 0 0 0 0.0 0 1 0 0 0 5.0 [4,] -0.2 -0.8 0 -0.4 -0.2 0 0 1 0 0 -0.2 0 0 1 0 0 2.6 [5,] 0.0 0.0 0 1.0 0.0 0 0 0 1 0 0.0 0 0 0 1 0 3.0 [6,] -0.6 -2.4 0 2.8 -0.6 0 0 0 0 -1 -0.6 0 0 0 0 1 7.8 [1] 3 12 13 14 15 16 [1] 1 2 4 5 6 7 8 9 10 11 [1] 1.4 6.0 5.0 2.6 3.0 7.8 [1] 24.4 [1,] 6.2 0.8 3.4 -0.8 1 1 1 1 -1 -1.8 [1,] 1 2 3 0 5 6 7 8 0 0 [1] 1.0 6.2 [1] 1 [1] 3 1 b [1,] 0 0.76 1 0.4 0.2 0 -0.04 0 0 0 0.2 0 -0.04 0 0 0 1.2 [2,] 0 2.60 0 0.0 0.0 1 -0.40 0 0 0 0.0 1 -0.40 0 0 0 4.0 [3,] 1 0.20 0 0.0 0.0 0 0.20 0 0 0 0.0 0 0.20 0 0 0 1.0 [4,] 0 -0.76 0 -0.4 -0.2 0 0.04 1 0 0 -0.2 0 0.04 1 0 0 2.8 [5,] 0 0.00 0 1.0 0.0 0 0.00 0 1 0 0.0 0 0.00 0 1 0 3.0 [6,] 0 -2.28 0 2.8 -0.6 0 0.12 0 0 -1 -0.6 0 0.12 0 0 1 8.4 [1] 3 12 1 14 15 16 [1] 2 4 5 6 7 8 9 10 11 13 [1] 1.2 4.0 1.0 2.8 3.0 8.4 [1] 18.2 [1,] -0.44 3.4 -0.8 1 -0.24 1 1 -1 -1.8 -1.24 [1,] 0 2 0 4 0 6 7 0 0 0 [1] 2.0 3.4 [1] 4 [1] 1 3 b [1,] 0 1.9 2.5 1 0.5 0 -0.1 0 0 0 0.5 0 -0.1 0 0 0 3.000000e+00 [2,] 0 2.6 0.0 0 0.0 1 -0.4 0 0 0 0.0 1 -0.4 0 0 0 4.000000e+00 [3,] 1 0.2 0.0 0 0.0 0 0.2 0 0 0 0.0 0 0.2 0 0 0 1.000000e+00 [4,] 0 0.0 1.0 0 0.0 0 0.0 1 0 0 0.0 0 0.0 1 0 0 4.000000e+00 [5,] 0 -1.9 -2.5 0 -0.5 0 0.1 0 1 0 -0.5 0 0.1 0 1 0 8.881784e-16 [6,] 0 -7.6 -7.0 0 -2.0 0 0.4 0 0 -1 -2.0 0 0.4 0 0 1 3.552714e-15 [1] 4 12 1 14 15 16 [1] 2 3 5 6 7 8 9 10 11 13 [1] 3.000000e+00 4.000000e+00 1.000000e+00 4.000000e+00 8.881784e-16 [6] 3.552714e-15 [1] 8 [1,] -6.9 -8.5 -2.5 1 0.1 1 1 -1 -3.5 -0.9 [1,] 0 0 0 4 5 6 7 0 0 0 [1] 4 1 [1] 6 [1] 2 4 b [1,] 0 1.9 2.5 1 0.5 0 -0.1 0 0 0 0.5 0 -0.1 0 0 0 3.000000e+00 [2,] 0 2.6 0.0 0 0.0 1 -0.4 0 0 0 0.0 1 -0.4 0 0 0 4.000000e+00 [3,] 1 0.2 0.0 0 0.0 0 0.2 0 0 0 0.0 0 0.2 0 0 0 1.000000e+00 [4,] 0 0.0 1.0 0 0.0 0 0.0 1 0 0 0.0 0 0.0 1 0 0 4.000000e+00 [5,] 0 -1.9 -2.5 0 -0.5 0 0.1 0 1 0 -0.5 0 0.1 0 1 0 8.881784e-16 [6,] 0 -7.6 -7.0 0 -2.0 0 0.4 0 0 -1 -2.0 0 0.4 0 0 1 3.552714e-15 [1] 4 6 1 14 15 16 [1] 2 3 5 7 8 9 10 11 12 13 [1] 3.000000e+00 4.000000e+00 1.000000e+00 4.000000e+00 8.881784e-16 [6] 3.552714e-15 [1] 4 [1,] -9.5 -8.5 -2.5 0.5 1 1 -1 -3.5 -1 -0.5 [1,] 0 0 0 4 5 6 0 0 0 0 [1] 5 1 [1] 8 [1] 4 4 b [1,] 0 1.9 2.5 1 0.5 0 -0.1 0 0 0 0.5 0 -0.1 0 0 0 3.000000e+00 [2,] 0 2.6 0.0 0 0.0 1 -0.4 0 0 0 0.0 1 -0.4 0 0 0 4.000000e+00 [3,] 1 0.2 0.0 0 0.0 0 0.2 0 0 0 0.0 0 0.2 0 0 0 1.000000e+00 [4,] 0 0.0 1.0 0 0.0 0 0.0 1 0 0 0.0 0 0.0 1 0 0 4.000000e+00 [5,] 0 -1.9 -2.5 0 -0.5 0 0.1 0 1 0 -0.5 0 0.1 0 1 0 8.881784e-16 [6,] 0 -7.6 -7.0 0 -2.0 0 0.4 0 0 -1 -2.0 0 0.4 0 0 1 3.552714e-15 [1] 4 6 1 8 15 16 [1] 2 3 5 7 9 10 11 12 13 14 [1] 3.000000e+00 4.000000e+00 1.000000e+00 4.000000e+00 8.881784e-16 [6] 3.552714e-15 [1] 4.440892e-15 [1,] -9.5 -9.5 -2.5 0.5 1 -1 -3.5 -1 -0.5 -1 [1,] 0 0 0 4 5 0 0 0 0 0 [1] 5 1 [1] 9 [1] 5.000000e+00 8.881784e-16 b [1,] 0 1.9 2.5 1 0.5 0 -0.1 0 0 0 0.5 0 -0.1 0 0 0 3.000000e+00 [2,] 0 2.6 0.0 0 0.0 1 -0.4 0 0 0 0.0 1 -0.4 0 0 0 4.000000e+00 [3,] 1 0.2 0.0 0 0.0 0 0.2 0 0 0 0.0 0 0.2 0 0 0 1.000000e+00 [4,] 0 0.0 1.0 0 0.0 0 0.0 1 0 0 0.0 0 0.0 1 0 0 4.000000e+00 [5,] 0 -1.9 -2.5 0 -0.5 0 0.1 0 1 0 -0.5 0 0.1 0 1 0 8.881784e-16 [6,] 0 -7.6 -7.0 0 -2.0 0 0.4 0 0 -1 -2.0 0 0.4 0 0 1 3.552714e-15 [1] 4 6 1 8 9 16 [1] 2 3 5 7 10 11 12 13 14 15 [1] 3.000000e+00 4.000000e+00 1.000000e+00 4.000000e+00 8.881784e-16 [6] 3.552714e-15 [1] 3.552714e-15 [1,] -7.6 -7 -2 0.4 -1 -3 -1 -0.6 -1 -1 [1,] 0 0 0 4 0 0 0 0 0 0 [1] 4.0 0.4 [1] 7 [1] 5.000000e+00 8.881784e-15 [1,] 0 0.000000e+00 4.440892e-16 1 0.000000e+00 0 0 0 1 0 0.000000e+00 0 [2,] 0 -5.000000e+00 -1.000000e+01 0 -2.000000e+00 1 0 0 4 0 -2.000000e+00 1 [3,] 1 4.000000e+00 5.000000e+00 0 1.000000e+00 0 0 0 -2 0 1.000000e+00 0 [4,] 0 0.000000e+00 1.000000e+00 0 0.000000e+00 0 0 1 0 0 0.000000e+00 0 [5,] 0 -1.900000e+01 -2.500000e+01 0 -5.000000e+00 0 1 0 10 0 -5.000000e+00 0 [6,] 0 -8.881784e-16 3.000000e+00 0 -2.220446e-16 0 0 0 -4 -1 -2.220446e-16 0 b [1,] 0 0 1 0 3.000000e+00 [2,] 0 0 4 0 4.000000e+00 [3,] 0 0 -2 0 1.000000e+00 [4,] 0 1 0 0 4.000000e+00 [5,] 1 0 10 0 8.881784e-15 [6,] 0 0 -4 1 3.944305e-31 [1] 4 6 1 8 7 16 [1] 2 3 5 9 10 11 12 13 14 15 [1] 3.000000e+00 4.000000e+00 1.000000e+00 4.000000e+00 8.881784e-15 [6] 3.944305e-31 [1] 3.944305e-31 [1,] -8.881784e-16 3 -2.220446e-16 -4 -1 -1 -1 -1 -1 -5 [1,] 0 2 0 0 0 0 0 0 0 0 [1] 2 3 [1] 3 [1] 6.000000e+00 1.314768e-31 [1,] 0 1.314768e-31 0 1 3.286920e-32 0 0 0 1.000000 1.480297e-16 [2,] 0 -5.000000e+00 0 0 -2.000000e+00 1 0 0 -9.333333 -3.333333e+00 [3,] 1 4.000000e+00 0 0 1.000000e+00 0 0 0 4.666667 1.666667e+00 [4,] 0 2.960595e-16 0 0 7.401487e-17 0 0 1 1.333333 3.333333e-01 [5,] 0 -1.900000e+01 0 0 -5.000000e+00 0 1 0 -23.333333 -8.333333e+00 [6,] 0 -2.960595e-16 1 0 -7.401487e-17 0 0 0 -1.333333 -3.333333e-01 b [1,] 3.286920e-32 0 0 0 1.000000 -1.480297e-16 3.000000e+00 [2,] -2.000000e+00 1 0 0 -9.333333 3.333333e+00 4.000000e+00 [3,] 1.000000e+00 0 0 0 4.666667 -1.666667e+00 1.000000e+00 [4,] 7.401487e-17 0 0 1 1.333333 -3.333333e-01 4.000000e+00 [5,] -5.000000e+00 0 1 0 -23.333333 8.333333e+00 8.881784e-15 [6,] -7.401487e-17 0 0 0 -1.333333 3.333333e-01 1.314768e-31 [1] 4 6 1 8 7 3 [1] 2 5 9 10 11 12 13 14 15 16 [1] 3.000000e+00 4.000000e+00 1.000000e+00 4.000000e+00 8.881784e-15 [6] 1.314768e-31 [1] 0 [1,] 0 0 0 0 -1 -1 -1 -1 -1 -1 [1,] 0 0 0 0 0 0 0 0 0 0 [1] "done" [1] "******donewithphase1****************" [[1]] [1] 4 6 1 8 7 3 [[2]] [1] 1.000000e+00 0.000000e+00 1.314768e-31 3.000000e+00 0.000000e+00 [6] 4.000000e+00 8.881784e-15 4.000000e+00 0.000000e+00 0.000000e+00 [11] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 [16] 0.000000e+00 [[3]] [1] 0 [[4]] [1] 0 [1] "phase1resultsabove" [1] 6 10 [1] 4 6 1 8 7 3 [1] 4 6 1 8 7 3 [1] 2 5 9 10 [1] 3 4 1 4 0 0 [1] -19 [1,] 4 -1 -4 0 [1,] 1 0 0 0 [1] 1 4 [1] 2 [1] 3.00 0.25 b [1,] 0.00 0 0 1 0.00 0 0 0 1.000000 0.0000000 3.00 [2,] 1.25 0 0 0 -0.75 1 0 0 -3.500000 -1.2500000 5.25 [3,] 0.25 1 0 0 0.25 0 0 0 1.166667 0.4166667 0.25 [4,] 0.00 0 0 0 0.00 0 0 1 1.333333 0.3333333 4.00 [5,] 4.75 0 0 0 -0.25 0 1 0 -1.166667 -0.4166667 4.75 [6,] 0.00 0 1 0 0.00 0 0 0 -1.333333 -0.3333333 0.00 [1] 4 6 2 8 7 3 [1] 1 5 9 10 [1] 3.00 5.25 0.25 4.00 4.75 0.00 [1] -20 [1,] -1 -2 -8.666667 -1.666667 [1,] 0 0 0 0 [1] "done" > prob2 [[1]] [1] 4 6 2 8 7 3 [[2]] [1] 0.00 0.25 0.00 3.00 0.00 5.25 4.75 4.00 0.00 0.00 [[3]] [1] 0 [[4]] [1] -20 > # > # third problem from Wagner -- unbounded > m <- 2 > n <- 4 > A <- matrix(c(1,-1,-1,1,1,0,0,1),m,n) > b <- c(1,1) > c <- c(-1,-1,0,0) > A [,1] [,2] [,3] [,4] [1,] 1 -1 1 0 [2,] -1 1 0 1 > b [1] 1 1 > c [1] -1 -1 0 0 > prob3 <- p1p2(A,b,c,m,n) [1] 2 4 [1] 2 6 [1] 2 6 [1] 5 6 [1] 5 6 [1] 1 2 3 4 [1] 1 1 [1] 2 [1,] 0 0 1 1 [1,] 0 0 3 4 [1] 3 1 [1] 3 [1] 1 1 b [1,] 1 -1 1 0 1 0 1 [2,] -1 1 0 1 0 1 1 [1] 3 6 [1] 1 2 4 5 [1] 1 1 [1] 1 [1,] -1 1 1 -1 [1,] 0 2 3 0 [1] 2 1 [1] 2 [1] 2 1 b [1,] 0 0 1 1 1 1 2 [2,] -1 1 0 1 0 1 1 [1] 3 2 [1] 1 4 5 6 [1] 2 1 [1] 0 [1,] 0 0 -1 -1 [1,] 0 0 0 0 [1] "done" [1] "******donewithphase1****************" [[1]] [1] 3 2 [[2]] [1] 0 1 2 0 0 0 [[3]] [1] 0 [[4]] [1] 0 [1] "phase1resultsabove" [1] 2 4 [1] 3 2 [1] 3 2 [1] 1 4 [1] 2 1 [1] -1 [1,] 2 -1 [1,] 1 0 [1] 1 2 [1] 1 [1] "unbound" > prob3 [[1]] [1] 3 2 [[2]] [1] 0 1 2 0 [[3]] [1] 2 [[4]] [1] -1 > # > # fourth problem from Wagner -- no solution > # same A, c, but negate b > b <- -b > A [,1] [,2] [,3] [,4] [1,] 1 -1 1 0 [2,] -1 1 0 1 > b [1] -1 -1 > c [1] -1 -1 0 0 > prob4 <- p1p2(A,b,c,m,n) [1] 2 4 [1] 2 6 [1] 2 6 [1] 5 6 [1] 5 6 [1] 1 2 3 4 [1] 1 1 [1] 2 [1,] 0 0 -1 -1 [1,] 0 0 0 0 [1] "done" [1] "******donewithphase1****************" [[1]] [1] 5 6 [[2]] [1] 0 0 0 0 1 1 [[3]] [1] 0 [[4]] [1] 2 [1] "phase1resultsabove" [1] "nofeasiblesolution" > prob4 [[1]] [1] 5 6 [[2]] [1] 0 0 0 0 1 1 [[3]] [1] 3 [[4]] [1] 0 > # > # fifth problem is L1 regression > m <- 6 # number of obs > n <- 16 # 2*nobs + 2*(cols of X) > X <- cbind(rep(1,m),(1:m)) # design matrix > A <- cbind(diag(rep(1,m)),diag(rep(-1,m)),X,-X) > b <- c(4.5,5.5,6.5,8,10,12) # y's > c <- c(rep(1,2*m),0,0,0,0) > A [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [1,] 1 0 0 0 0 0 -1 0 0 0 0 0 1 1 [2,] 0 1 0 0 0 0 0 -1 0 0 0 0 1 2 [3,] 0 0 1 0 0 0 0 0 -1 0 0 0 1 3 [4,] 0 0 0 1 0 0 0 0 0 -1 0 0 1 4 [5,] 0 0 0 0 1 0 0 0 0 0 -1 0 1 5 [6,] 0 0 0 0 0 1 0 0 0 0 0 -1 1 6 [,15] [,16] [1,] -1 -1 [2,] -1 -2 [3,] -1 -3 [4,] -1 -4 [5,] -1 -5 [6,] -1 -6 > b [1] 4.5 5.5 6.5 8.0 10.0 12.0 > c [1] 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 > p <- (1:m) # we have feasible soln, b positive > prob5 <- simplex(p,A,b,c,m,n) [1] 6 16 [1] 1 2 3 4 5 6 [1] 1 2 3 4 5 6 [1] 7 8 9 10 11 12 13 14 15 16 [1] 4.5 5.5 6.5 8.0 10.0 12.0 [1] 46.5 [1,] -2 -2 -2 -2 -2 -2 6 21 -6 -21 [1,] 0 0 0 0 0 0 7 8 0 0 [1] 8 21 [1] 14 [1] 4 2 b [1,] 1 0 0 -0.25 0 0 -1 0 0 0.25 0 0 0.75 0 -0.75 0 2.5 [2,] 0 1 0 -0.50 0 0 0 -1 0 0.50 0 0 0.50 0 -0.50 0 1.5 [3,] 0 0 1 -0.75 0 0 0 0 -1 0.75 0 0 0.25 0 -0.25 0 0.5 [4,] 0 0 0 0.25 0 0 0 0 0 -0.25 0 0 0.25 1 -0.25 -1 2.0 [5,] 0 0 0 -1.25 1 0 0 0 0 1.25 -1 0 -0.25 0 0.25 0 0.0 [6,] 0 0 0 -1.50 0 1 0 0 0 1.50 0 -1 -0.50 0 0.50 0 0.0 [1] 1 2 3 14 5 6 [1] 4 7 8 9 10 11 12 13 15 16 [1] 2.5 1.5 0.5 2.0 0.0 0.0 [1] 4.5 [1,] -5.25 -2 -2 -2 3.25 -2 -2 0.75 -0.75 0 [1,] 0 0 0 0 5 0 0 8 0 0 [1] 5.00 3.25 [1] 10 [1] 5 0 b [1,] 1 0 0 0 -0.2 0 -1 0 0 0 0.2 0 0.8 0 -0.8 0 2.5 [2,] 0 1 0 0 -0.4 0 0 -1 0 0 0.4 0 0.6 0 -0.6 0 1.5 [3,] 0 0 1 0 -0.6 0 0 0 -1 0 0.6 0 0.4 0 -0.4 0 0.5 [4,] 0 0 0 0 0.2 0 0 0 0 0 -0.2 0 0.2 1 -0.2 -1 2.0 [5,] 0 0 0 -1 0.8 0 0 0 0 1 -0.8 0 -0.2 0 0.2 0 0.0 [6,] 0 0 0 0 -1.2 1 0 0 0 0 1.2 -1 -0.2 0 0.2 0 0.0 [1] 1 2 3 14 10 6 [1] 4 5 7 8 9 11 12 13 15 16 [1] 2.5 1.5 0.5 2.0 0.0 0.0 [1] 4.5 [1,] -2 -2.6 -2 -2 -2 0.6 -2 1.4 -1.4 0 [1,] 0 0 0 0 0 6 0 8 0 0 [1] 8.0 1.4 [1] 13 [1] 3.00 1.25 b [1,] 1 0 -2.0 0 1.0 0 -1 0 2.0 0 -1.0 0 0 0 0 0 1.50 [2,] 0 1 -1.5 0 0.5 0 0 -1 1.5 0 -0.5 0 0 0 0 0 0.75 [3,] 0 0 2.5 0 -1.5 0 0 0 -2.5 0 1.5 0 1 0 -1 0 1.25 [4,] 0 0 -0.5 0 0.5 0 0 0 0.5 0 -0.5 0 0 1 0 -1 1.75 [5,] 0 0 0.5 -1 0.5 0 0 0 -0.5 1 -0.5 0 0 0 0 0 0.25 [6,] 0 0 0.5 0 -1.5 1 0 0 -0.5 0 1.5 -1 0 0 0 0 0.25 [1] 1 2 13 14 10 6 [1] 3 4 5 7 8 9 11 12 15 16 [1] 1.50 0.75 1.25 1.75 0.25 0.25 [1] 2.75 [1,] -3.5 -2 -0.5 -2 -2 1.5 -1.5 -2 0 0 [1,] 0 0 0 0 0 6 0 0 0 0 [1] 6.0 1.5 [1] 9 [1] 2.0 0.5 [1,] 1 -1.3333333 0 0 0.3333333 0 -1 1.3333333 0 0 -0.3333333 0 0 0 0 0 [2,] 0 0.6666667 -1 0 0.3333333 0 0 -0.6666667 1 0 -0.3333333 0 0 0 0 0 [3,] 0 1.6666667 0 0 -0.6666667 0 0 -1.6666667 0 0 0.6666667 0 1 0 -1 0 [4,] 0 -0.3333333 0 0 0.3333333 0 0 0.3333333 0 0 -0.3333333 0 0 1 0 -1 [5,] 0 0.3333333 0 -1 0.6666667 0 0 -0.3333333 0 1 -0.6666667 0 0 0 0 0 [6,] 0 0.3333333 0 0 -1.3333333 1 0 -0.3333333 0 0 1.3333333 -1 0 0 0 0 b [1,] 0.5 [2,] 0.5 [3,] 2.5 [4,] 1.5 [5,] 0.5 [6,] 0.5 [1] 1 9 13 14 10 6 [1] 2 3 4 5 7 8 11 12 15 16 [1] 0.5 0.5 2.5 1.5 0.5 0.5 [1] 2 [1,] -1 -2 -2 -1 -2 -1 -1 -2 0 0 [1,] 0 0 0 0 0 0 0 0 0 0 [1] "done" > prob5 [[1]] [1] 1 9 13 14 10 6 [[2]] [1] 0.5 0.0 0.0 0.0 0.0 0.5 0.0 0.0 0.5 0.5 0.0 0.0 2.5 1.5 0.0 0.0 [[3]] [1] 0 [[4]] [1] 2 > # done > rm(list=ls()) > q()