代码之家 › 专栏 › 技术社区 › coolhand

应用lpSolveAPI

lpsolve mathematical-optimization r

coolhand · 技术社区 · 2 年前

我有许多变量表示我可以更换的设备。当被替换时,它们改善了影响度量, [I] 。每一项都有相应的年度成本节约 [S] 以及重置成本 [C] .

n <- 1000 # variable count

# impact
# negative for use in minimization function
I <- -rnorm(n, mean=20000, sd=8000)
# cost savings
s <- rnorm(n, mean=2500, sd=1000)
# replacement cost
c <- rnorm(n, mean=15000, sd=5000)

我想选择要替换的组件,以在预算范围内全面最大化影响,同时确保整个项目(作为一个整体)达到一个简单的回报目标。

payback_goal <- 3
budget <- 1000000

这个问题由以下方程式描述。

我正在努力设置这个 lpSolveAPI 具体来说,我不知道如何合并方程式3。

library(lpSolveAPI)
m <- 2 # number of constraints, disregarding binary constraint set by type
my.lp <- make.lp(m, n)
set.row(my.lp, 1, c)
# i don't think this is the correct way to set up the calculation
set.row(my.lp, 2, c/s) # cost divided by savings

# create obj fn coefficients
set.objfn(my.lp, I)
set.type(my.lp, 1:n, "binary")

eq <- c("<=", "<=")
set.constr.type(my.lp, eq)
set.rhs(my.lp, c(budget, payback_goal))

solve(my.lp)
soln1 <- get.variables(my.lp)

# total cost
s1_cost <- sum(soln1 * c)
# total impact
s1_impact <- -get.objective(my.lp)
# total simple payback
s1_pb <- sum(soln1*c) / sum(soln1*s)

我尝试了一种有点棘手的解决方法,即假设总成本将相当接近预算,这就得出了第三个方程式

但这只是一个近似值,我想知道如何在R代码中更精确地实现它。

0 回复 | 直到 2 年前

Erwin Kalvelagen 2 年前

您可以使用符号将(3)线性化如下:

 sum([X]*[C]) - payback*sum([X]*[S]) <= 0

这可以改写为:;

 sum([X]*([C]-payback*[S])) <= 0

只要满足以下条件,这就是有效的 sum([X]*[S])>0 .

ThomasIsCoding 2 年前

正如所指出的那样 @Erwin Kalvelagen ,您可以重新制定约束使其呈线性。

但是,你应该意识到,这种转变是有条件的 sum([X]*[S]) ,这意味着

如果 sum([X]*[S])>0 ,我们应该有 sum([X]*([C]-payback*[S])) <= 0
如果 sum([X]*[S])<0 ,我们应该有 sum([X]*([C]-payback*[S])) >= 0

然后,结合上述两种情况(比较并取最小值),我们可以通过以下方式排除影响 总和([X]*[S]) .

代码示例(in `CVXR` )

set.seed(0)
n <- 1000 # variable count

# impact
# negative for use in minimization function
I <- -rnorm(n, mean = 20000, sd = 8000)
# cost savings
s <- rnorm(n, mean = 2500, sd = 1000)
# replacement cost
c <- rnorm(n, mean = 15000, sd = 5000)


library(CVXR)
X <- Variable(1, n, boolean = TRUE)
payback_goal <- 3
budget <- 1000000

obj <- Minimize(X %*% I)
# case 1: conditional on X %*% s > 0
constr1 <- list(
  X %*% c <= budget,
  X %*% (c - payback_goal * s) <= 0,
  X %*% s >= 0
)
problem1 <- Problem(obj, constr1)
sol1 <- solve(problem1)

# case 2: conditional on X %*% s < 0
constr2 <- list(
  X %*% c <= budget,
  X %*% (c - payback_goal * s) >= 0,
  X %*% s <= 0
)
problem2 <- Problem(obj, constr2)
sol2 <- solve(problem2)

# compare two cases and take the minimum between them as the solution
if (round(sol1$getValue(X)) %*% I > round(sol2$getValue(X)) %*% I) {
  sol <- sol2
} else {
  sol <- sol1
}

你会看到的

> round(sol1$getValue(X)) %*% I
         [,1]
[1,] -3064996

> round(sol2$getValue(X)) %*% I
          [,1]
[1,] -287548.7

这样 sol 被分配了 sol1

应用lpSolveAPI

代码示例(in CVXR )

代码示例(in `CVXR` )