精华内容
下载资源
问答
  • What's the recommended package for constrained non-linear optimization in python ?The specific problem I'm trying to solve is this:I have an unknown X (Nx1), I have M (Nx1) u vectors and M (NxN) s mat...

    What's the recommended package for constrained non-linear optimization in python ?

    The specific problem I'm trying to solve is this:

    I have an unknown X (Nx1), I have M (Nx1) u vectors and M (NxN) s matrices.

    max [5th percentile of (ui_T*X), i in 1 to M]

    st

    0<=X<=1 and

    [95th percentile of (X_T*si*X), i in 1 to M]<= constant

    When I started out the problem I only had one point estimate for u and s and I was able to solve the problem above with cvxpy.

    I realized that instead of one estimate for u and s, I had the entire distribution of values so I wanted to change my objective function so that I could use the entire distribution. The problem description above is my attempt to include that information in a meaningful way.

    cvxpy cannot be used to solve this, I've tried scipy.optimize.anneal, but I can't seem to set bounds on the unknown values. I've looked at pulp too but it doesnt allow nonlinear constraints.

    解决方案

    scipy has a spectacular package for constrained non-linear optimization.

    You can get started by reading the optimize doc, but here's an example with SLSQP:

    minimize(func, [-1.0,1.0], args=(-1.0,), jac=func_deriv, constraints=cons, method='SLSQP', options={'disp': True})

    展开全文
  • 惩罚函数法在模拟退火算法求解非线性约束优化问题中的应用,陈思源,,本文首先介绍了模拟退火算法和惩罚函数法的基本原理和方法,然后将其结合成求解非线性约束优化问题的算法。在Matlab语言环境下编制
  • python求解带约束目标优化问题(非线性规划,粒子群,遗传,差分进化)

    情况1:输出值可以是浮点数

    算例1

    书上的答案

    该算例是一个带约束的目标问题

    方法1 非线性规划 scipy.optimize.minimize
    非线性规划原理就不讲解啦

    针对算例1
    求取一个函数的最小值。函数的参数可以是多个,但函数值只能是标量。

    参数

    • fun : callable 目标函数
    • x0 : ndarry初始值
    • args : tuple, optional额外的参数,传给目标函数和它的导数。
    • method : str or callable, optional求解问题的算法名,下面选其一:Nelder-Mead, Powell, CG, BFGS, Newton-CG, L-BFGS-B, TNC, COBYLA, SLSQP,dogleg, trust-ncg默认是 BFGS, L-BFGS-B, SLSQP 之一,根据问题是否含有约束和界限自动选择。
    • jac : bool or callable, optional目标函数的梯度矩阵。只适用于 CG, BFGS, Newton-CG, L-BFGS-B, TNC, SLSQP, dogleg, trust-ncg。如果jac是一个 Boolean 且为 True,则 fun 被认为是梯度与目标函数一起返回。如果是False,则梯度会被自动地计算。jac也可以是一个函数,返回目标函数的梯度,且参数必须与fun相同。
    • hess, hessp : callable, optional目标函数的二阶导矩阵,或者二阶导矩阵乘以一个随机向量p。只适用于Newton-CG,dogleg, trust-ncg。hess和hessp只需要给出一个即可。如果提供了hess,则hessp会被忽略。如果两者都没有提供,则二阶导矩阵会被自动计算
    • bounds : sequence, optionalbounds 是参数的界限,只适用于L-BFGS-B, TNC 和 SLSQP,每个参数对应一个 (min, max),表示参数的上下限。如果只有一边界限,则另一边置为None。当约束是针对 x x x xx x xxxx 中的单个元素的上下限时,就可以用 bounds 参数来设置。
    • constraints : dict or sequence of dict, optional
      约束定义,只适用于 COBYLA 和 SLSQP。每个约束定义为一个词典,键值对包括:
      fun : callable。定义了约束函数。
      type : str。约束类型: eq’ 表示等式约束(fun等于0),ineq 表示不等式约束(fun大于等于0)。COBYLA只支持不等式约束。
      jac : callable, optional。fun 的梯度矩阵,只适用于SLSQP
      args : sequence, optional。传递给fun和jac的额外参数。
    • tol : float, optional迭代终止的允许误差。
    • options : dict, optional求解器的选项字典。所有的算法都接受以下的通用选项:maxiter : int。迭代的最大次数。disp : bool。如果是True则打印出收敛信息。
    • callback : callable, optional
      每次迭代之后调用的函数,参数为xk,表示当前的参数向量。

    返回值
    res:优化结果。

    优化结果是OptimizeResult对象,重要属性如下:

    fun 是最优值
    x 是最优解
    success 表示求解器是否成功退出。
    message 描述了求解器退出的原因

    from scipy.optimize import minimize
    import numpy as np
    #目标函数即min(FG1+FG2+FG3)
    def fun(x):
    
        return (4+0.3*x[0]+0.0007*x[0]*x[0]+3+0.32*x[1]+0.0004*x[1]*x[1]+3.5+0.3*x[2]+0.00045*x[2]*x[2])
    
    
    def con():
        # 约束条件 分为eq 和ineq
        #eq表示 函数结果等于0 ; ineq 表示 表达式大于等于0
    
    
        cons = ({'type': 'eq', 'fun': lambda x: x[0]+x[1]+x[2]-700})
    
                #{'type': 'ineq', 'fun': lambda x: -x[2] + x2max}#如果有不等式约束
       #cons= ([con1, con2, con3, con4, con5, con6, con7, con8]) #如果有多个约束,则最后返回结果是这个
        
    
        return cons
    
    #上下限约束
    b1=(100,200)#
    b2=(120,250)#
    b3=(150,300)#
    bnds = (b1,b2,b3)#边界约束
    if __name__ == "__main__":
        cons=con()#约束
    
    
        # 设置x初始猜测值
        x0 = np.array((150, 250, 20))
        res = minimize(fun, x0, method='SLSQP', constraints=cons,bounds=bnds)
        print('代价',res.fun)
        print(res.success)
        print('解',res.x)
    
    

    方法2 粒子群pyswarm.pso
    粒子群PSO优化算法学习笔记 及其python实现(附讲解如何使用python语言sko.PSO工具包)

    pyswarm是一个支持带约束的粒子群优化包

    sko.PSO中的pso仅支持带上下限的约束,不支持等式和不等式约束。

    参数详解
    pso(func, lb, ub, ieqcons=[], f_ieqcons=None, args=(), kwargs={},
    swarmsize=100, omega=0.5, phip=0.5, phig=0.5, maxiter=100,
    minstep=1e-8, minfunc=1e-8, debug=False)

    • func : function要最小化的函数
    • lb : array设计变量的下界
    • ub : array设计变量的上界
    • ieqcons : list在一个成功优化的问题中,长度为n的函数列表,使ieqconsj >= 0.0(默认值:[])
    • f_ieqcons : function返回一个一维数组,其中每个元素都必须大于或等于成功优化的问题中的0.0。如果指定了f_ieqcons,则忽略ieqcons(默认值:None)
    • args : tuple传递给目标函数和约束函数的附加参数(默认:空元组)
    • kwargs : dict传递给目标函数和约束函数的附加关键字参数(默认值:空字典)
    • swarmsize : int集群中的粒子数(默认值:100)
    • omega : scalar粒子速度缩放因子(默认值:0.5)
    • phip : scalar缩放因子搜索远离粒子最知名的位置(默认:0.5)
    • phig : scalar从蜂群最知名的位置搜索的缩放因子(默认值:0.5)
    • maxiter : int 最大迭代次数
    • minstep : scalar在searc终止前,群最佳位置的最小步长(默认值:1e-8)
    • minfunc : scalar在searchterminates之前群的最佳目标值的最小变化(默认值:1e-8)
    • debug : boolean如果为True,则每次迭代都会显示进度语句(默认值:False)

    返回值

    • g : array群体最知名的位置(最优设计)
    • f : scalarg点的目标值
    from pyswarm import pso
    
    def object_func(x):
        return (4+0.3*x[0]+0.0007*x[0]*x[0]+3+0.32*x[1]+0.0004*x[1]*x[1]+3.5+0.3*x[2]+0.00045*x[2]*x[2])
    
    #不等式约束
    
    def cons1(x):
        return [x[0]+x[1]+x[2]-700]
    
    lb = [100, 120, 150]#
    ub = [200, 250, 300]
    
    xopt, fopt = pso(object_func,lb,ub,ieqcons=[cons1], maxiter=100,swarmsize=1000)
    print(xopt)
    print(fopt)
    
    

    粒子群求得的值不一定是最优解

    例子2

    from pyswarm import pso
    
    def object_func(x):
        x1 = x[0]
        x2 = x[1]
        return x1**4 - 2*x2*x1**2 + x2**2 + x1**2 - 2*x1 + 5 # 求解最小值
    
    def cons(x):
        x1 = x[0]
        x2 = x[1]
        con1 =  -(x1 + 0.25)**2 + 0.75*x2 # con1>=0
        con2 = -x1+0.2 # con2>=0
        return [-(x1 + 0.25)**2 + 0.75*x2,-x1+0.2]
    
    lb = [-3, -1] # -3<x1<2, -1<x2<6
    ub = [2, 6]
    
    xopt, fopt = pso(object_func, lb, ub, f_ieqcons=cons)
    print(xopt)
    print(fopt)
    
    

    方法3 遗传算法

    遗传算法原理及其python实现

    使用工具from sko.GA import GA
    参数详解
    func, n_dim,size_pop=50, max_iter=200,prob_mut=0.001,
    lb=-1, ub=1,constraint_eq=tuple(), constraint_ueq=tuple(),precision=1e-7

    func : function .The func you want to do optimal
    n_dim : int . number of variables of func
    lb : array_like The lower bound of every variables of func
    ub : array_like The upper bound of every vaiiables of func
    constraint_eq : list .equal constraint
    constraint_ueq : list . unequal constraint
    precision : array_like. The precision of every vaiiables of func
    size_pop : int .Size of population
    max_iter : int. Max of iter
    prob_mut : float .between 0 and Probability of mutation
    Attributes
    ----------------------
    Lind : array_like
    The num of genes of every variable of func(segments)
    generation_best_X : array_like. Size is max_iter.
    Best X of every generation
    generation_best_ranking : array_like. Size if max_iter.
    Best ranking of every generation

    
    def object_func(x):
        return (4+0.3*x[0]+0.0007*x[0]*x[0]+3+0.32*x[1]+0.0004*x[1]*x[1]+3.5+0.3*x[2]+0.00045*x[2]*x[2])
    
    #等式约束
    
    def cons1(x):
        return [x[0]+x[1]+x[2]-700]
    cons=cons1
    
    #导入包
    from sko.GA import GA
    ga = GA(func=object_func, n_dim=3, size_pop=200, max_iter=800, lb=[100, 120, 150], ub=[200, 250, 300],constraint_eq=[cons])
    best_x, best_y = ga.run()
    print('best_x:', best_x, '\n', 'best_y:', best_y)
    

    同理很难找到最优解

    方法4 差分进化算法

    '''
    min f(x1, x2, x3) = x1^2 + x2^2 + x3^2
    s.t.
        x1*x2 >= 1
        x1*x2 <= 5
        x2 + x3 = 1
        0 <= x1, x2, x3 <= 5
    '''
    
    
    def obj_func(p):
        x1, x2, x3 = p
        return x1 ** 2 + x2 ** 2 + x3 ** 2
    
    
    constraint_eq = [
        lambda x: 1 - x[1] - x[2]
    ]
    
    constraint_ueq = [
        lambda x: 1 - x[0] * x[1],
        lambda x: x[0] * x[1] - 5
    ]
    from sko.DE import DE
    
    de = DE(func=obj_func, n_dim=3, size_pop=50, max_iter=800, lb=[0, 0, 0], ub=[5, 5, 5],
            constraint_eq=constraint_eq, constraint_ueq=constraint_ueq)
    
    best_x, best_y = de.run()
    print('best_x:', best_x, '\n', 'best_y:', best_y)
    

    情况2:指定输出值为整数

    方法1:遗传算法做整数规划**
    在多维优化时,想让哪个变量限制为整数,就设定 precision 为 1即可
    例如,我想让我的自定义函数object_func 的第一个变量限制为整数,那么就设定 precision 的第一个数为1,例子如下

    
    def object_func(x):
        return (4+0.3*x[0]+0.0007*x[0]*x[0]+3+0.32*x[1]+0.0004*x[1]*x[1]+3.5+0.3*x[2]+0.00045*x[2]*x[2])
    
    #等式约束
    
    def cons1(x):
        return [x[0]+x[1]+x[2]-700]
    cons=cons1
    
    #导入包
    from sko.GA import GA
    ga = GA(func=object_func, n_dim=3, size_pop=200, max_iter=800, lb=[100, 120, 150], ub=[200, 250, 300],constraint_eq=[cons], precision=[1, 1e-7, 1e-7])
    best_x, best_y = ga.run()
    print('best_x:', best_x, '\n', 'best_y:', best_y)
    
    

    在这里插入图片描述
    作者:电气-余登武。写作不容易,请点个赞再走。在这里插入图片描述

    展开全文
  • 软件库:scipy.optimize, numpy, CVXPY,Gekko 软件:octave 5.1,matlab 本文将介绍三种计算非线性约束优化的方法: (1)scipy.optimize.minimize (2)cvxpy (3)Octave 5.1 sqp函数 (4)matlab ga函数 博主其他...

    软件库:scipy.optimize, numpy, CVXPY,Gekko
    软件:octave 5.1,matlab

    本文将介绍三种计算非线性约束优化的方法:
    (1)scipy.optimize.minimize
    (2)cvxpy
    (3)Octave 5.1 sqp函数
    (4)matlab ga函数

    博主其他相关博文:
    CVXPY/Scipy/Octave线性规划问题求解


    2020-08-13 更新
    使用cvxpy库的时候,对于有些优化问题需要注意转换一下形式,例如下面二次规划问题(https://www.inverseproblem.co.nz/OPTI/index.php/Probs/QP):
    在这里插入图片描述
    需要先将优化问题转成二次规划的一般形式:
    minf(x)=12xHx+fx\text{min} f(\bold{x})=\frac{1}{2} \bold{x}'\bold{H}\bold{x}+\bold{f}'\bold{x}
    其中x\bold{x}列向量,可以算出H=[1112]\bold{H}=\begin{bmatrix} 1 & -1 \\ -1 & 2\end{bmatrix}f=[26]\bold{f}=\begin{bmatrix} -2 \\ -6 \end{bmatrix},相应的cvxpy计算程序如下(参考官方例子https://www.cvxpy.org/examples/basic/quadratic_program.html):

    import numpy as np
    import cvxpy as cvx
    x = cvx.Variable(2)
    H = np.array([[1,-1],
                 [-1,2]])
    f = np.array([-2,-6])
    A = np.array([[1,1],
         [-1,2],
         [2,1]])
    b = np.array([2,2,3])
    
    # obj = cvx.Minimize(1/2*x[0]**2+x[1]**2-x[0]*x[1]-2*x[0]-6*x[1])
    obj = cvx.Minimize((1/2)*cvx.quad_form(x,H)+f.T@x)
    cons = [A@x<=b, x[0]>=0]
    
    prob = cvx.Problem(obj, cons)
    prob.solve()
    
    print(prob.status, prob.value, x.value)
    

    在上面例子中,如果采用obj = cvx.Minimize(1/2*x[0]**2+x[1]**2-x[0]*x[1]-2*x[0]-6*x[1])计算会提示求解错误,需要处理为obj = cvx.Minimize((1/2)*cvx.quad_form(x,H)+f.T@x)形式,具体为什么前面的形式无法求解,笔者也还需进一步学习。


    2019-10-19补充Gekko库非线性约束优化
    使用Gekko优化计算库求解博文中的一个非线性约束优化例子(无法用cvxpy求解,报错提示优化问题不遵循DCP规则),程序代码如下:

    from gekko import GEKKO
    import time
    
    start = time.time_ns()
    m = GEKKO() # Initialize gekko
    # Use IPOPT solver (default)
    m.options.SOLVER = 3
    # Change to parallel linear solver
    # m.solver_options = ['linear_solver ma97']
    # Initialize variables
    x1 = m.Var(value=1,lb=0,ub=5)
    x2 = m.Var(value=2,lb=0,ub=5)
    
    # Equations
    m.Equation(x1**2+x2**2<=6)
    m.Equation(x1*x2==2)
    
    m.Obj(x1**3-x2**2+x1*x2+2*x1**2) # Objective
    
    m.options.IMODE = 3 # Steady state optimization
    
    
    m.solve(disp=False) # Solve
    print('Results')
    print('x1: ' + str(x1.value))
    print('x2: ' + str(x2.value))
    end = time.time_ns()
    
    print('Objective: ' + str(m.options.objfcnval))
    print('Time used: {:.2f} s'.format((end-start)*1e-9))
    '''
    Results
    x1: [0.87403204792]
    x2: [2.2882456138]
    Objective: -1.040502879
    Time used: 24.37 s
    '''
    

    相比于scipy.optimize.minimize和octave,使用gekko求解该问题耗时很久,估计是解法器没有选择到最优,因为刚接触gekko,该问题后面再研究。
    Gekko默认使用公共服务器资源进行计算,网络传输会影响程序执行速度,在初始化Gekko时可设置为本地计算求解,只要设置m = Gekko(remote=False)即可,执行结果为:
    在这里插入图片描述


    MATLAB ga函数

    2019-05-10更新
    matlab有自带的遗传算法函数ga,相关语法如下图:
    在这里插入图片描述

    octave有个ga工具箱,同样提供了ga函数,程序和前面的类似,更改一下函数定义的位置即可。

    对于下面octave中求解Himmelblau函数:
    f(x,y)=(x2+y11)2+(x+y27)2f(x,y)=(x^{2} +y-11)^{2}+(x+y^2-7)^2的最优解,在matlab中程序如下:

    %% MATLAB ga求解最优化问题
    
    lb=[-5,-5];
    ub=[5,5]; 
    
    %% initial guess
    x0 = [-3,-2];
    %% MATLAB ga函数,详细信息请查看帮助文档
    x  = ga (@phi,2,[],[],[],[],lb,ub)
    
    %% object function,matlab中函数定义要放电一个脚本的最后面,这和octave脚本函数位置要求不同
    function obj = phi (x)
      obj = (x(1)^2+x(2)-11)^2+(x(1)+x(2)^2-7)^2;
    end
    
    

    结果为:

    x = 1×2    
        3.0000    2.0000
    
    

    从结构可以看出像ga这种进化算法只得到了一个最优解。因为matlab中的ga函数还可以处理线性约束和非线性约束问题,那么我们继续求解下文中的两个例子。

    Example1

    注意A,b系数矩阵要转换为Ax<=b的格式

    % MATLAB ga求解最优化问题,注意A,b系数矩阵要转换为Ax<=b的格式
    A = [-1,2;1,2;1,-2];
    b = [2;6;2];
    lb=[0,0];
    ub=[]; 
    option = optimoptions('ga','PlotFcn',@gaplotbestf)
    %% MATLAB ga函数
    
    x  = ga (@phi,2,A,b,[],[],lb,ub,[],option)
    
    %% object function
    
    function obj = phi (x)
      obj = (x(1)-1)^2+(x(2)-2.5)^2;
    end
    
    

    结果为:

    x = 1×2    
        1.4224    1.7117
    
    

    可见其结果离真实解还有一定偏差,为此我们可以控制ga的option来提供精度,将option修改为如下:

    option = optimoptions('ga','ConstraintTolerance',1e-6,'PlotFcn',@gaplotbestf)
    

    运行结果为:
    在这里插入图片描述

    x = 1×2    
        1.3941    1.6970
    
    

    Example 2

    option = optimoptions('ga','ConstraintTolerance',1e-6,'PlotFcn',@gaplotbestf)
    %% MATLAB ga函数
    x  = ga (@phi,2,[],[],[],[],[],[],@cons,option)
    
    %% object function
    function obj = phi (x)
      obj = x(1)^3-x(2)^3+x(1)*x(2)+2*x(1)^2;
    end
    
    %% nonlinear constraints,非线性约束条件
    function [c,ceq] = cons(x)
        c = x(1)^2+x(2)^2-6; % Compute nonlinear inequalities at x,不等式约束,c(x)<=0
        ceq = x(1)*x(2)-2;       % Compute nonlinear equalities at x,等式约束,ceq(x)=0
    end
    

    结果如下:
    在这里插入图片描述
    Optimization terminated: stall generations limit exceeded
    but constraints are not satisfied.

    x =

    1.0255    1.9512
    

    可见计算结果并不满足约束条件,也和真实解(x1 = 0.8740320488968545
    x2 = 2.2882456112715115)相差比较大。

    ga作为一种随机性的进化算法,ga本身的参数(如crossover、mutation等)都会对结果造成很大影响,这里就不继续探讨了。

    scipy.optimize.minimize

    相关函数:
    scipy.optimize.minimize(fun, x0, args=(), method=None, jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None)

    func:求解的目标函数,Python callable function
    x0:给定的初始值
    args:目标函数中的参数
    method:求解算法
    bounds:边界约束
    constraints:条件约束

    其他详细内容请查看官方帮助文档。

    这部分内容为scipy官方文档翻译过来,但官方文档对带约束的最优化问题介绍的比较少,其引用了一本书籍的内容,下面记录之。

    (1)根据文档所知,参考书籍为Nocedal, J, and S J Wright. 2006. Numerical Optimization. Springer New York.
    使用例子来自example16.3(p462),文档有误,问题描述如下:
    在这里插入图片描述

    代码如下:

    import numpy as np
    from scipy.optimize import minimize
    
    # 目标函数
    def objective(x):
        return (x[0] - 1)**2 + (x[1] - 2.5)**2
    
    # 约束条件
    def constraint1(x):
        return x[0] - 2 * x[1] + 2  #不等约束
    
    def constraint2(x):
        
        return -x[0] - 2 * x[1] + 6 #不等约束
    
    def constraint3(x):
            
        return -x[0] + 2 * x[1] + 2 #不等约束
    
    # 初始猜想
    n = 2
    x0 = np.zeros(n)
    x0[0] = 2
    x0[1] = 0
    
    
    # show initial objective
    print('Initial SSE Objective: ' + str(objective(x0)))
    
    # 边界约束
    b = (0.0,None)
    bnds = (b, b) # 注意是两个变量都要有边界约束
    
    con1 = {'type': 'ineq', 'fun': constraint1} 
    con2 = {'type': 'ineq', 'fun': constraint2}
    con3 = {'type': 'ineq', 'fun': constraint3} 
    cons = ([con1,con2,con3]) # 3个约束条件
    
    # 优化计算
    solution = minimize(objective,x0,method='SLSQP',\
                        bounds=bnds,constraints=cons)
    x = solution.x
    
    # show final objective
    print('Final SSE Objective: ' + str(objective(x)))
    
    # print solution
    print('Solution')
    print('x1 = ' + str(x[0]))
    print('x2 = ' + str(x[1]))
    
    

    结果如下:

    Initial SSE Objective: 7.25
    Final SSE Objective: 0.8000000011920985
    Solution
    x1 = 1.4000000033834283
    x2 = 1.7000000009466527
    

    就是当x1=1.4,x2=1.7时在约束条件下目标函数取得最小值。

    再举一个例子加深印象,参照https://jingyan.baidu.com/album/6c67b1d69508b52787bb1edf.html?picindex=2,求解下面优化问题:
    在这里插入图片描述

    使用python代码如下:

    import numpy as np
    from scipy.optimize import minimize
    
    def objective(x):
        return x[0]**3-x[1]**3+x[0]*x[1]+2*x[0]**2
    
    def constraint1(x):
        return -x[0]**2 - x[1]**2 + 6.0
    
    def constraint2(x):
        
        return 2.0-x[0]*x[1]
    
    # initial guesses
    n = 2
    x0 = np.zeros(n)
    x0[0] = 1
    x0[1] = 2
    
    
    # show initial objective
    print('Initial SSE Objective: ' + str(objective(x0)))
    
    # optimize
    b = (0.0,3.0)
    bnds = (b, b)
    con1 = {'type': 'ineq', 'fun': constraint1} 
    con2 = {'type': 'eq', 'fun': constraint2}
    cons = ([con1,con2])
    solution = minimize(objective,x0,method='SLSQP',\
                        bounds=bnds,constraints=cons)
    x = solution.x
    
    # show final objective
    print('Final SSE Objective: ' + str(objective(x)))
    
    # print solution
    print('Solution')
    print('x1 = ' + str(x[0]))
    print('x2 = ' + str(x[1]))
    
    

    结果如下:

    Initial SSE Objective: -3.0
    Final SSE Objective: -7.785844454002188
    Solution
    x1 = 0.8740320488968545
    x2 = 2.2882456112715115
    

    和原作者在matlab求解结果一致

    特别注意: 使用约束的时候,要转化为形如ax+by>=c的形式,像前面第二个例子中的不等约束为x2+y2<=6x^2+y^2<=6,则要转化为x2y2+6>=0-x^2-y^2+6>=0的形式,然后再编写约束函数的代码:

    def constraint1(x):
        return -x[0]**2 - x[1]**2 + 6.0
    

    接着定义约束类型:

    con1 = {'type': 'ineq', 'fun': constraint1} 
    

    cvxpy

    使用cvxpy求解前面第一个问题,起程序如下:

    import cvxpy as cp
    import numpy as np
    
    '''
    minimize f(x)=(x-1)^2+(y-2.5)^2
    s.t.
    x-2y+2>=0
    -x-2y+6>=0
    -x+2y+2>=0
    x>=0
    y>=0
    '''
    A = np.array([[1, -2],
                [-1, -2],
                [-1, 2]])
    
    b = np.array([-2, -6, -2])
    
    x = cp.Variable(2)
    prob=cp.Problem(cp.Minimize( (x[0]-1)**2 + (x[1]-2.5)**2),
                [A@x>=b,x[0]>=0,x[1]>=0])
    prob.solve()
    
    # show final objective
    print('Final SSE Objective: ' + str(prob.value))
    
    # print solution
    print('Solution of x and y: ' + str(x.value))
    
    

    运行结果为:

    Final SSE Objective: 0.8000000000000005
    Solution of x and y: [1.4 1.7]
    
    

    对于第二个例子笔者用cvxpy暂时解不出来,程序报错:

    cvxpy.error.DCPError: Problem does not follow DCP rules. Specifically:
    The objective is not DCP. Its following subexpressions are not:
    

    因为求解不成功就不放出程序了,以后研究透了再更新。


    Octave5.1

    相关函数:

    […] = sqp (x0, phi, g, h, lb, ub, maxiter, tol)
    

    该函数的官方文档介绍如下:

    Minimize an objective function using sequential quadratic programming (SQP). 
    Solve the nonlinear program 
    min phi (x)
     x 
    subject to 
    g(x)  = 0
    h(x) >= 0
    lb <= x <= ub 
    using a sequential quadratic programming method. 
    

    参数说明(重点)
    (1)x0:初始猜想值,向量形式,如果是n个变量优化问题,x0也为n维向量;
    (2)phi:目标函数,可以使用函数句柄
    (3)g:等式约束条件,如果没有约束用空矩阵[]代替
    (4)h:不等式约束条件,如果没有约束用空矩阵[]代替
    (5)lb:变量左边界约束,同样注意和变量个数有关
    (6)ub:变量右边界约束,和变量个数有关

    先举个简单的例子,目标函数为Himmelblau函数,关于了解更多优化计算中测试算法使用的函数,可以参考网站http://infinity77.net/global_optimization/index.html下面的Test Function部分内容,也可在维基百科中搜索,这里就不过多介绍了。

    Himmelblau函数:
    f(x,y)=(x2+y11)2+(x+y27)2f(x,y)=(x^{2} +y-11)^{2}+(x+y^2-7)^2
    在这里插入图片描述
    在这里插入图片描述
    使用octave求解代码如下:

    %%OCTAVE脚本里面有函数时,需要先有其他语句开始,不然就会当做是函数文件
    1;
    %% object function
    function obj = phi (x)
      obj = (x(1)^2+x(2)-11)^2+(x(1)+x(2)^2-7)^2;
    endfunction
    lb=[-5,-5];
    ub=[-5,5]; 
    
    %% initial guess
    x0 = [-3,-2];
    %% 没有约束条件,用[]代替
    [x, obj, info, iter, nf, lambda] = sqp (x0, @phi, [],[],lb,ub)
    

    结果如下:

    >> nonlinear_prob
    
    x =
    
      -3.7793
      -3.2832
    
    obj =    1.5554e-14
    info =  104
    iter =  8
    nf =  17
    lambda =
    
       0
       0
       0
       0
    

    可见在初始值(-3,-2)情况下找到了最优解(-3.7793,-3.2832)

    下面求解前面带约束的例子:

    1;
    %% equality constraints
    function r = g (x)
      r = [ x(1)*x(2)-2];
    endfunction
    
    function r = h (x)
      r = [ 6-x(1)^2-x(2)^2];
    endfunction
    
    %% object function
    function obj = phi (x)
      obj = x(1)^3-x(2)^3+x(1)*x(2)+2*x(1)^2;
    endfunction
    lb=[0,0];
    ub=[3,3]; 
    
    %% initial guess
    x0 = [1,1];
    
    [x, obj, info, iter, nf, lambda] = sqp (x0, @phi, @g, @h, lb, ub)
    
    

    结果如下:

    
    x =
    
       0.87403
       2.28825
    
    obj = -7.7858
    info =  104
    iter =  6
    nf =  8
    lambda =
    
       7.03150
       4.58428
       0.00000
       0.00000
       0.00000
       0.00000
    

    计算得出的x解和scipy.optimize求出的一致,这点以后有时间再验证一下。

    ————
    2019-04-17补充
    优路径最短问题可以归结为最优化问题,即求解一个点到其他点距离的最小值。

    minimize

    minimizef(x,y)=((xdataxi)2+(ydatayi)2) minimize f(x,y)=\sum((x-datax_i)^2+(y-datay_i)^2)
    这里为了方便计算,使用了距离的平方,没有开根号计算。

    根据前面的目标函数,我们写成符合scipy.optimzie.minimize的函数输入格式,最后完整的程序如下:

    import numpy as np
    from scipy.optimize import minimize
    import matplotlib.pyplot as plt
    
    datax = 10*np.random.randn(15)
    datay = 10*np.random.rand(15)
    
    def f(x):
        return sum((x[0]-datax[i])**2+(x[1]-datay[i])**2 for i in range(len(datax)))
    
    # 初值
    x0 = [1,1]
    # 计算结果
    opt = minimize(f,x0=x0)
    print(opt)
    
    plt.figure()
    plt.scatter(datax,datay,label='Original points')
    plt.scatter(opt.x[0], opt.x[1], s=240, marker='h',label='Center point')
    plt.legend()
    # plt.show()
    for i in range(len(datax)):
        plt.plot([opt.x[0],datax[i]],[opt.x[1],datay[i]])
    plt.show()
    

    结果如下:

          fun: 762.1598697487602
     hess_inv: array([[ 3.33359580e-02, -8.50862577e-06],
           [-8.50862577e-06,  3.33328322e-02]])
          jac: array([ 0.00000000e+00, -7.62939453e-06])
      message: 'Optimization terminated successfully.'
         nfev: 28
          nit: 5
         njev: 7
       status: 0
      success: True
            x: array([0.89233268, 5.09306968])
    

    在这里插入图片描述

    CVXPY

    import cvxpy as cp
    import numpy as np
    import matplotlib.pyplot as plt
    datax = [2.5,6.7,3,9.3,9.2,9,4,7.2,0.3,2.8,10.3,2.4,5.7,8.1,3.5,6.89,1.34]
    datay = [5.7,6.2,5.0,7.9,2.7,9.3,8,2,9,4,6.8,1.9,3.8,7.2,5.71,4.88,2.3]
    
    def f(x):
        return sum((x[0]-datax[i])**2+(x[1]-datay[i])**2 for i in range(len(datax)))
    
    x = cp.Variable(2)
    prob=cp.Problem(cp.Minimize(sum((x[0]-datax[i])**2+(x[1]-datay[i])**2 for i in range(len(datax)))))
    prob.solve()
    
    # show final objective
    print('Final SSE Objective: ' + str(prob.value))
    
    # print solution
    print('Solution of x and y: ' + str(x.value))
    
    plt.figure()
    plt.scatter(datax,datay,label='Original points')
    plt.scatter(x.value[0], x.value[1], s=240, marker='h',label='Center point')
    plt.legend()
    # plt.show()
    for i in range(len(datax)):
        plt.plot([x.value[0],datax[i]],[x.value[1],datay[i]])
    plt.show()
    
    

    对于求解线性规划问题,请参考文首给出的链接。

    参考链接:
    [1]https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html#rdd2e1855725e-5
    [2]https://apmonitor.com/pdc/index.php/Main/NonlinearProgramming
    [3]https://jingyan.baidu.com/article/6c67b1d69508b52787bb1edf.html

    展开全文
  • python求解非线性优化问题

    千次阅读 2020-05-04 18:12:02
    一元 例题: 代码: 二元 例题: 代码:

    一元

    例题:
    在这里插入图片描述
    代码:
    在这里插入图片描述

    二元

    例题:

    在这里插入图片描述
    代码:
    在这里插入图片描述
    二元问题也可转化为一元的方法来解决
    例二是简便后的结果,例一方法便于理解
    欢迎讨论、提问和交流

    展开全文
  • 本文提纲一维搜索/单变量优化问题无约束多元优化问题非线性最小二乘问题约束优化问题非线性规划问题的目标函数或约束条件是非线性的。本文使用SciPy的optimize模块来求解非线性规划问题。目标函数和约束条件是否...
  • 今天小编就为大家分享一篇python 非线性规划方式(scipy.optimize.minimize),具有很好的参考价值,希望对大家有所帮助。一起跟随小编过来看看吧
  • python 非线性规划(scipy.optimize.minimize)

    万次阅读 多人点赞 2018-08-09 13:48:34
    背景:现在项目上有一个用python 实现非线性规划的需求。非线性规划可以简单分两种,目标函数为凸函数 or 非凸函数。 凸函数的 非线性规划,比如fun=x^2+y^2+x*y,有很多常用的python库来完成,网上也有很多资料,...
  • 非线性规划问题分为两种问题:无约束约束。 不受约束优化问题是形式上的问题 min⁡f(x),x=[x1,x2,⋯ ,xn]T∈Rn  (1)\min \boldsymbol{f}\left( \boldsymbol{x} \right) ,\boldsymbol{x}=\left[ \boldsymbol{x...
  • x0:非线性优化问题的求解并不是用解析的方法进行的,而是从初始值开始进行迭代,x0就是猜测的初始值,在很大程度上决定了算法的效果。另外要注意,如果是多元函数,那么就要将x0视为一个向量,每一个分量都要单独...
  • 各种规划一直神秘的不敢去触碰...规划论包括线性规划、非线性规划和动态规划线性规划问题概念理解:当目标函数与约束条件都是线形的,则称为线性规划。线性规划问题是在一组线性约束条件的限制下,求一线性目标函数...
  • 经典遗传算法(SGA)解非线性优化问题的原理及其python代码实现
  • 系统的真实模型都是非线性,在实际模型构建或使用过程中,...如果模型简单,姑且可以手动求解,模型一旦复杂些(例如大型的非线性模型,约束函数为微分方程or差分方程形式),非常容易出错。虽然也有一些库能够用...
  • 我有这样的优化:Maximize (x1*a1 + x2*a2+...+xn * an)such that :any j, (xj*aj)/sum(all xi*ai) = 0.15any three j, (xj*aj)/sum(all xi*ai) >= 0.07any j 1
  • 1、遗传算法(GA)介绍 ...在求解较为复杂的组合优化问题时,相对一些常规的优化算法,通常能够较快地获得较好的优化结果。遗传算法已被人们广泛地应用于组合优化、机器学习、信号处理、自适应控制和人工生
  • 模拟退火算法(Simulated Annealing,SA)是一种随机寻优算法,其出发点是基于物理中固体物质的退火过程与一般组合优化问题之间的相似性。 模拟退火算法从某一较高初温出发,伴随温度参数的不断下降,结合概率突跳特性...
  • 今天小编就为大家分享一篇使用Python求解带约束的最优化问题详解,具有很好的参考价值,希望对大家有所帮助。一起跟随小编过来看看吧
  • 这个函数是用来求解双重线性约束问题的,满足公式如下 其中 minimize 部分中约束条件的第一和第三条可以推导出:Gx <= h,而这种形式也是我们在日常应用中会遇到的最常见的形式,有了这种不等关系,我们就可以将...
  • 如前位答者所说,Python中的SciPy库可以...提纲一维搜索/单变量优化问题无约束多元优化问题非线性最小二乘问题约束优化问题非线性规划问题的目标函数或约束条件是非线性的。这里使用SciPy的optimize模块来求解非...
  • 1.多峰非线性优化问题 非线性优化问题 2. 经典遗传算法 2.1 遗传算法概述
  • NLopt是一个用于非线性局部和全局优化的库,用于具有和不具有梯度信息的函数。 它被设计为简单,统一的界面,并包装了多个免费/开源非线性优化库。 可以从Github上的页面下载最新版本,并且托管在readthedocs上。 ...
  • 本篇针对非线性规(优)划(化)问题...两个工具都可以接受非线性模型形式,并且有与各种优化求解器的接口。 首先声明,这篇文章是转载的,先上链接 https://blog.csdn.net/u013468614/article/details/103624851 另...
  • 一般使用默认 # constraints:约束条件,针对fun中为参数的部分进行约束限制 # 1.计算 1/x+x 的最小值 from scipy.optimize import minimize import numpy as np #计算 1/x+x 的最小值 def fun(args): a=args v=...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 6,455
精华内容 2,582
关键字:

python非线性约束优化

python 订阅