Programación Lineal

Ejemplo

TOYCO arma tres juguetes: trenes, camiones y coches, con tres operaciones. Los límites diarios de tiempo disponible para las tres operaciones son 430, 460 y 420 minutos, respectivamente, y las utilidades por tren, camión y coche de juguete son $3, $2 y $5, respectivamente. Los tiempos de armado por tren, en las tres operaciones son 1, 3 y 1 minutos, respectivamente. Los tiempos respectivos por camión y por coche son (2, 0, 4) y (1, 2, 0) minutos (un tiempo de cero indica que no se usa la operación).

Usando Paquete lpSolve

library('lpSolve')

Función Objetivo

\[ max \space z = 3x_1 + 2x_2 + 5x_3 \] Sujeto a restricciones

\[ \begin{aligned} x_1 + 2x_2 + x_3 \leq 430 \\\\ 3x_1 + 2x_3 \leq 460 \\\\ x_1 + 4x_2 \leq 420 \\\\ x_1,x_2 \geq 0 \end{aligned} \]

Paso 1: Definiendo los límtes

obj.fun <- c(3,2,5)
restricciones <- matrix(c(1, 2, 1, 
                          3, 0, 2,
                          1, 4, 0),
                          nrow=3, byrow = TRUE)
direccion <- c("<=", "<=", "<=")
rhs <- c(430, 460, 420)

Paso 2: Resolviendo el modelo

solucion <- lp("max", obj.fun, restricciones, direccion, rhs, compute.sens = TRUE )
solucion$status
## [1] 0

status: 0 = Exito, 2 = Solución No Factible

Paso 3: Obtención de solución

solucion$objval
## [1] 1350
solucion$solution
## [1]   0 100 230
optimo <- solucion$solution
names(optimo) <- c("x1","x2","x3")
print(optimo)
##  x1  x2  x3 
##   0 100 230
print(paste("Ganancia Total: ", solucion$objval, sep=""))
## [1] "Ganancia Total: 1350"
solucion$constraints
##               [,1] [,2] [,3]
##                  1    3    1
##                  2    0    4
##                  1    2    0
## const.dir.num    1    1    1
## const.rhs      430  460  420
solucion$duals
## [1]  1  2  0 -4  0  0

Usando Paquete lpSolveAPI

library('lpSolveAPI')

4 restricciones y 2 variables de decisión

lprec <- make.lp(nrow = 3, ncol = 3)

Definiendo el tipo de problema a resolver

lp.control(lprec, sense="max")
## $anti.degen
## [1] "fixedvars" "stalling" 
## 
## $basis.crash
## [1] "none"
## 
## $bb.depthlimit
## [1] -50
## 
## $bb.floorfirst
## [1] "automatic"
## 
## $bb.rule
## [1] "pseudononint" "greedy"       "dynamic"      "rcostfixing" 
## 
## $break.at.first
## [1] FALSE
## 
## $break.at.value
## [1] 1e+30
## 
## $epsilon
##       epsb       epsd      epsel     epsint epsperturb   epspivot 
##      1e-10      1e-09      1e-12      1e-07      1e-05      2e-07 
## 
## $improve
## [1] "dualfeas" "thetagap"
## 
## $infinite
## [1] 1e+30
## 
## $maxpivot
## [1] 250
## 
## $mip.gap
## absolute relative 
##    1e-11    1e-11 
## 
## $negrange
## [1] -1e+06
## 
## $obj.in.basis
## [1] TRUE
## 
## $pivoting
## [1] "devex"    "adaptive"
## 
## $presolve
## [1] "none"
## 
## $scalelimit
## [1] 5
## 
## $scaling
## [1] "geometric"   "equilibrate" "integers"   
## 
## $sense
## [1] "maximize"
## 
## $simplextype
## [1] "dual"   "primal"
## 
## $timeout
## [1] 0
## 
## $verbose
## [1] "neutral"

Definiendo el tipo de variables de decisión

Opciones:

integer : Entero binary : Binario (0,1) real : Real

set.type(lprec, 1:3, type=c("real"))

Definición de los coeficientes de la función objetivo

set.objfn(lprec, obj.fun)

Adición de restricción

add.constraint(lprec, restricciones[1, ], "<=", rhs[1])
add.constraint(lprec, restricciones[2, ], "<=", rhs[2])
add.constraint(lprec, restricciones[3, ], "<=", rhs[3])

Visualización de matrix del modelo

lprec
## Model name: 
##             C1    C2    C3           
## Maximize     3     2     5           
## R1           0     0     0  free    0
## R2           0     0     0  free    0
## R3           0     0     0  free    0
## R4           1     2     1    <=  430
## R5           3     0     2    <=  460
## R6           1     4     0    <=  420
## Kind       Std   Std   Std           
## Type      Real  Real  Real           
## Upper      Inf   Inf   Inf           
## Lower        0     0     0

Solución del problema

solve(lprec)
## [1] 0

Valores de las variables de decisión

get.variables(lprec)
## [1]   0 100 230

Valor de la función objetivo

get.objective(lprec)
## [1] 1350

Tenga en cuenta que los límites predeterminados en la variable de decisión son c(0, 0, 0) y c(Inf, Inf, Inf)

get.bounds(lprec)
## $lower
## [1] 0 0 0
## 
## $upper
## [1] Inf Inf Inf

Programación Cuadrática

La programación cuadrática (QP) funciona cuando la función objetivo es una función cuadrática, es decir, contiene hasta dos productos. Aquí las funciones de restricción siguen siendo una combinación lineal de variables.

Podemos expresar el problema en forma matricial.

Minimizar Objetivo: \[\frac{1}{2} X^TDX - d^TX\]

donde \(X\) es el vector \([x_1, \ldots, x_n]^T\), \(D\) es la matriz de pesos de cada par \(x_ix_j\) y \(d\) son los pesos de cada \(x_i\). \(\frac{1}{2}\) proviene del hecho de que \(D\) es simétrico y, por lo tanto, cada \(x_ix_j\) se cuenta dos veces.

con restricciones:

\[A^TX[=| \geq] ~ b\]

donde los primeros operadores \(k\) son igualdad, los otros son \(\geq\) y \(b\) los valores a los que las restricciones deberían ser iguales.

Por ejemplo, una función objetivo cuadrática: \[f(x_1, x_2, x_3) = 2x_1^2 - x_1x_2 + 2 x_2^2 + x_2x_3 + 2x_3^2 - 5x_2 + 3x_3\]

Sujeto a restricciones:

\[ \begin{aligned} - 4x_1 & + -3x_2 & & = -8 \\ 2x_1 & + x_2 & & = 2 \\ & -2x_2 & + x_3 & \geq 0 \end{aligned} \]

library(quadprog)

Matriz de coeficientes de la función objetivo

D <- matrix(c( 2,-1, 0,
               -1, 2,-1,
                0,-1, 2),3,3)

Vector de Coeficientes de la función objetivo lineales

dvec <- c(0,-5,3)

Matriz de restricciones y lado derecho

Amat <- matrix(c(-4,-3, 0,
                        2, 1, 0,
                        0,-2, 1),3,3)
bvec <- c(-8,2,0)

Las dos primeras restricciones son igualdades

n.eqs <- 2

Solución

sol <- solve.QP(D,dvec,Amat,bvec=bvec,meq=2)
sol$solution
## [1] -1  4  8
solucion.q <- sol$solution
names(solucion.q) <- c("x1","x2","x3")
print(solucion.q)
## x1 x2 x3 
## -1  4  8
sol$value
## [1] 49
print(paste("Valor Mínimo: ", sol$val, sep=""))
## [1] "Valor Mínimo: 49"