使用 NLOPT/Gurobi 解决混合约束优化答案

【问题标题】：Using NLOPT/Gurobi for solving mixed constraint optimization使用 NLOPT/Gurobi 解决混合约束优化
【发布时间】：2018-01-23 15:58:45
【问题描述】：

我目前正在做一个项目，我想使用R和NLOPT包（或Gurobi）来解决以下优化问题：

求最小值 ||y-y_h||_L^2 使得 x = Ay_h, y >= 0，其中 x, y 是给定大小为 16*1 的向量，A = 16* 24矩阵也给出了。

我的尝试：

R 代码

nrow=16;
ncol = 24;
lambda = matrix(sample.int(100, size = ncol*nrow, replace = T),nrow,ncol);
lambda = lambda - diag(lambda)*diag(x=1, nrow, ncol);
y = rpois(ncol,lambda) + rtruncnorm(ncol,0,1,mean = 0, sd = 1); 

x = matrix (0, nrow, 1);
x_A1 = y[1]+y[2]+y[3];
x_A2 = y[4]+y[7]+y[3];
x_B1 = y[4]+y[5]+y[6];
x_B2 = y[11]+y[1];
x_C1 = y[7]+y[8]+y[9];
x_C2 = y[2]+y[5]+y[12];
x_D1 = y[10]+y[11]+y[12];
x_D2 = y[3]+y[6]+y[9];
x_E1 = y[13]+y[14]+y[15];
x_E2 = y[18]+y[19]+y[23];
x_F1 = y[20]+y[21]+y[19];
x_F2 = y[22]+y[16]+y[13];
x_G1 = y[23]+y[22]+y[24];
x_G2 = y[14]+y[17]+y[20];
x_H1 = y[16]+y[17]+y[18];
x_H2 = y[15]+y[21]+y[24];

d <- c(x_A1, x_A2,x_B1, x_B2,x_C1, x_C2,x_D1, x_D2,x_E1, 
       x_E2,x_F1, x_F2,x_G1, x_G2,x_H1, x_H2)
x <- matrix(d, nrow, byrow=TRUE)

A = matrix(c(1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, #x_A^1 
             0,0,0,1,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0, #x_A^2 
             0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, #x_B^1
             1,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0, #x_B^2 
             0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, #x_C^1 
             0,1,0,0,1,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0, #x_C^2 
             0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0, #x_D^1 
             0,0,1,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, #x_D^2 
             0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0, #x_E^1 
             0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,1,0, #x_E^2
             0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0, #x_F^1 
             0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,0,0,0,1,0,0, #x_F^2 
             0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,1,0,0,0,0, #x_G^2 
             0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1, #x_G^1 
             0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0, #x_H^1 
             0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,1,0,0,1), #x_H^2
           nrow, ncol, byrow= TRUE)

尝试了两个代码来解决问题：min ||y - y_h||_L^2 其中 x= Ay_h, y>=0 其中 x,y,A 都在上面给出。

# f(x) = ||yhat-y||_L2

eval_f <- function( yhat ) {
  return( list( "objective" = norm((mean(yhat-y))^2, type = "2")))
}

# inequality constraint
eval_g_ineq <- function( yhat ) {
  constr <- c(0 - yhat)
  return( list( "constraints"=constr ))
}

# equalities constraint
eval_g_eq <- function( yhat ) {
  constr <- c( x-A%*%yhat )
  return( list( "constraints"=constr ))
}

x0 <- y

#lower bound of control variable
lb <- c(matrix (0, ncol, 1))

local_opts <- list( "algorithm" = "NLOPT_LD_MMA",
                    "xtol_rel" = 1.0e-7 )
opts <- list( "algorithm" = "NLOPT_LD_AUGLAG",
              "xtol_rel" = 1.0e-7,
              "maxeval" = 1000,
              "local_opts" = local_opts )
res <- nloptr( x0=x0,
               eval_f=eval_f, 
               eval_grad_f = NULL,
               lb=lb,
               eval_g_ineq = eval_g_ineq, 
               eval_g_eq=eval_g_eq,
               opts=opts)
print(res)

古罗比代码：

**#model <- list()
#model$B <- A
#model$obj <- norm((y-yhat)^2, type = "2")
#model$modelsense <- "min"
#model$rhs <- c(x,0)
#model$sense <- c('=', '>=')
#model$vtype <- 'C'
#result <- gurobi(model, params)
#print('Solution:')
#print(result$objval)
#print(result$yhat)**

我的问题：首先，当我运行上面的 R 代码时，它一直给我这样的信息： is.nloptr(ret) 中的错误：目标梯度中的元素数量错误另外：警告信息：在 is.na(f0$gradient) 中： is.na() 应用于“NULL”类型的非（列表或向量）

我试图避免计算梯度，因为我没有关于 y 的密度函数的任何信息。谁能帮我解决上面的错误？

对于 Gurobi 代码，我收到以下消息：错误：is(model$A, "matrix") ||是（模型$A，“稀疏矩阵”）|| is(model$A, .... is not TRUE

但是我的矩阵A输入正确，那么这个错误是什么意思？

【问题讨论】：

请不要在此处发布代码图片。包含一个reproducible example 和示例输入数据，我们可以复制/粘贴到 R 中运行。
@MrFlick：我修好了。非常感谢，因为我是新来的。但是得知 LATEX 在这个论坛上不工作很烦人；p 你能帮我解决这个问题吗，因为我已经花了大约 7 个小时在它上面，仍然无法看到如何避免使用渐变，但仍然使“NLOPT”代码工作：（
这里没有人想帮我修复上面的代码？？；）

标签： r nonlinear-optimization

【解决方案1】：

我几天前才开始使用 nloptr。这个问题已经是老问题了，但我仍然会回答。当您使用带有 'NLOPT_LD_AUGLAG' 算法的 'nloptr' 时，'LD' 代表本地并使用渐变......所以你需要选择中间有 'LN' 的其他东西。例如，“NLOPT_LN_COBYLA”应该可以在没有渐变的情况下正常工作。其实你可以看看 nloptr 包手册。

【讨论】：