Matlab到Julia优化：JuMP中的函数@SetNLObjective

Question

我将尝试在Julia中重写Matlab的fmincon优化函数。

以下是Matlab代码：

function [x,fval] = example3()

  x0 = [0; 0; 0; 0; 0; 0; 0; 0];  
  A = [];                         
  b = [];
  Ae = [1000 1000 1000 1000 -1000 -1000 -1000 -1000];   
  be = [100];                     
  lb = [0; 0; 0; 0; 0; 0; 0; 0];  
  ub = [1; 1; 1; 1; 1; 1; 1; 1];  
  noncon = [];                  

  options = optimset('fmincon');
  options.Algorithm = 'interior-point';

  [x,fval] = fmincon(@objfcn,x0,A,b,Ae,be,lb,ub,@noncon,options);

end

function f = objfcn(x) 

  % user inputs
  Cr = [ 0.0064 0.00408 0.00192 0; 
       0.00408 0.0289 0.0204 0.0119;
       0.00192 0.0204 0.0576 0.0336; 
       0 0.0119 0.0336 0.1225 ];

  w0 = [ 0.3; 0.3; 0.2; 0.1 ];
  Er = [0.05; 0.1; 0.12; 0.18];

  % calculate objective function
  w = w0+x(1:4)-x(5:8);
  Er_p = w'*Er;
  Sr_p = sqrt(w'*Cr*w);

  % f = objective function
  f = -Er_p/Sr_p;

end

以下是我的Julia代码：

using JuMP
using Ipopt

m = Model(solver=IpoptSolver())

# INPUT DATA
w0 = [ 0.3; 0.3; 0.2; 0.1 ]
Er = [0.05; 0.1; 0.12; 0.18]
Cr = [ 0.0064 0.00408 0.00192 0; 
    0.00408 0.0289 0.0204 0.0119;
    0.00192 0.0204 0.0576 0.0336; 
    0 0.0119 0.0336 0.1225 ]

# VARIABLES
@defVar(m, 0 <= x[i=1:8] <= 1, start = 0.0)
@defNLExpr(w, w0+x[1:4]-x[5:8])
@defNLExpr(Er_p, w'*Er)
@defNLExpr(Sr_p, w'*Cr*w)
@defNLExpr(f, Er_p/Sr_p)

# OBJECTIVE
@setNLObjective(m, Min, f)

# CONSTRAINTS 
@addConstraint(m, 1000*x[1] + 1000*x[2] + 1000*x[3] + 1000*x[4] - 
1000*x[5] - 1000*x[6] - 1000*x[7] - 1000*x[8] == 100)

# SOLVE
status = solve(m)

# DISPLAY RESULTS
println("x = ", round(getValue(x),4))
println("f = ", round(getObjectiveValue(m),4))

朱莉娅优化功能在我明确定义@setNLObjective中的目标函数时有效，但这是不合适的，因为用户的输入可能会发生变化，从而导致不同的目标函数，您可以从创建目标函数的方式中看到这一点。

问题似乎在于JuMP在@setNLObjective参数中输入目标函数的限制：

所有表达式必须是简单的标量操作。您不能使用点、矩阵-向量乘积、向量切片等。将向量操作转换为显式的sum{}操作。

有没有办法解决这个问题？还是Julia中有其他软件包可以解决这个问题，需要注意的是我将没有雅可比矩阵或海森矩阵。

非常感谢。

- kulsuri

---

NLopt 软件包没有这个限制。您可以传递任何函数或匿名函数。 - Colin T Bowers

我已经看了一下NLopt，它似乎是我需要的，但缺乏示例正在阻碍我的进展。您能否给我一个使用NLopt的简单示例？甚至可以使用我的函数在NLopt中进行演示。谢谢。 - kulsuri

我大约一周前在SO上询问了关于NLopt的问题，其中包括一个简单问题的工作示例。这里是链接。确保您已经在您的机器上运行它，如果您仍然有问题，请在此回复我，我会尽力提供更多帮助。 - Colin T Bowers

我使用你的示例成功地让它工作了，这应该能帮助我用NLopt完成项目。我遇到了一个小问题，目标函数是通过矩阵乘法创建的，并将其传递给min_objective!时出现错误。我通过创建一个变量来解决了这个问题，该变量的值是1x1矩阵的第一个元素，并将其传递给目标函数。我将我的代码作为下面的答案发布。感谢你的帮助！ - kulsuri

很高兴听到你已经解决了它！ - Colin T Bowers

Answer 1

使用Julia和NLopt优化包的Matlab代码的工作示例。

using NLopt

function objective_function(x::Vector{Float64}, grad::Vector{Float64})

    w0 = [ 0.3; 0.3; 0.2; 0.1 ]
    Er = [0.05; 0.1; 0.12; 0.18]
    Cr = [ 0.0064 0.00408 0.00192 0; 
            0.00408 0.0289 0.0204 0.0119;
            0.00192 0.0204 0.0576 0.0336; 
            0 0.0119 0.0336 0.1225 ]

    w = w0 + x[1:4] - x[5:8]

    Er_p = w' * Er

    Sr_p = sqrt(w' * Cr * w)

    f = -Er_p/Sr_p

    obj_func_value = f[1]

    return(obj_func_value)
end

function constraint_function(x::Vector, grad::Vector)

    constraintValue = 1000*x[1] + 1000*x[2] + 1000*x[3] + 1000*x[4] -
1000*x[5] - 1000*x[6] - 1000*x[7] - 1000*x[8] - 100        

return constraintValue
end

opt1 = Opt(:LN_COBYLA, 8)

lower_bounds!(opt1, [0, 0, 0, 0, 0, 0, 0, 0])
upper_bounds!(opt1, [1, 1, 1, 1, 1, 1, 1, 1])

#ftol_rel!(opt1, 0.5)
#ftol_abs!(opt1, 0.5)

min_objective!(opt1, objective_function)
equality_constraint!(opt1,  constraint_function)

(fObjOpt, xOpt, flag) = optimize(opt1, [0, 0, 0, 0, 0, 0, 0, 0])

- kulsuri