【问题标题】:Fast CVX solvers in MatlabMatlab 中的快速 CVX 求解器
【发布时间】:2017-04-20 10:17:56
【问题描述】:

我想知道 Matlab 中最快的凸优化器是什么,或者有什么方法可以加快当前求解器的速度?我正在使用 CVX,但要花很长时间才能解决我遇到的优化问题。 我的优化是解决

minimize norm(Ax-b, 2)
subject to
x >= 0
and x d <= delta

其中 A 和 b 的大小非常大。

有什么方法可以通过最小二乘求解器来解决这个问题,然后将其转移到约束版本以使其更快?

【问题讨论】:

  • norm(Ax,b) 对我来说看起来很奇怪。你的意思是norm(Ax-b,2)
  • x.d 是什么意思?
  • d 是另一个向量。理想情况下,我希望 x 和 d 是正交的,这由 delta 的值控制。

标签: matlab optimization least-squares convex-optimization cvx


【解决方案1】:

我不确定x.d &lt;= delta 是什么意思,但我假设它应该是x &lt;= delta

您可以使用投影梯度法或加速投影梯度法(这只是对投影梯度法的轻微修改,“神奇地”收敛得更快)来解决这个问题。这是一些显示如何最小化 .5|| 的 python 代码Ax - b ||^2 受限于 0 manuscript 中关于近端算法的文章。

import numpy as np
import matplotlib.pyplot as plt

def fista(gradf,proxg,evalf,evalg,x0,params):
    # This code does FISTA with line search

    maxIter = params['maxIter']
    t = params['stepSize'] # Initial step size
    showTrigger = params['showTrigger']

    increaseFactor = 1.25
    decreaseFactor = .5

    costs = np.zeros((maxIter,1))

    xkm1 = np.copy(x0)
    vkm1 = np.copy(x0)

    for k in np.arange(1,maxIter+1,dtype = np.double):

        costs[k-1] = evalf(xkm1) + evalg(xkm1)
        if k % showTrigger == 0:
            print "Iteration: " + str(k) + "    cost: " + str(costs[k-1])

        t = increaseFactor*t

        acceptFlag = False
        while acceptFlag == False:
            if k == 1:
                theta = 1
            else:
                a = tkm1
                b = t*(thetakm1**2)
                c = -t*(thetakm1**2)
                theta = (-b + np.sqrt(b**2 - 4*a*c))/(2*a)

            y = (1 - theta)*xkm1 + theta*vkm1
            (gradf_y,fy) = gradf(y)
            x = proxg(y - t*gradf_y,t)
            fx = evalf(x)
            if fx <= fy + np.vdot(gradf_y,x - y) + (.5/t)*np.sum((x - y)**2):
                acceptFlag = True
            else:
                t = decreaseFactor*t

        tkm1 = t
        thetakm1 = theta
        vkm1 = xkm1 + (1/theta)*(x - xkm1)
        xkm1 = x

    return (xkm1,costs)


if __name__ == '__main__':

    delta = 5.0
    numRows = 300
    numCols = 50
    A = np.random.randn(numRows,numCols)
    ATrans = np.transpose(A)
    xTrue = delta*np.random.rand(numCols,1)
    b = np.dot(A,xTrue)
    noise = .1*np.random.randn(numRows,1)
    b = b + noise

    def evalf(x):
        AxMinusb = np.dot(A, x) - b
        val = .5 * np.sum(AxMinusb ** 2)
        return val

    def gradf(x):
        AxMinusb = np.dot(A, x) - b
        grad = np.dot(ATrans, AxMinusb)
        val = .5 * np.sum(AxMinusb ** 2)
        return (grad, val)

    def evalg(x):
        return 0.0

    def proxg(x,t):
        return np.maximum(np.minimum(x,delta),0.0)

    x0 = np.zeros((numCols,1))
    params = {'maxIter': 500, 'stepSize': 1.0, 'showTrigger': 5}
    (x,costs) = fista(gradf,proxg,evalf,evalg,x0,params)

    plt.figure()
    plt.plot(x)
    plt.plot(xTrue)

    plt.figure()
    plt.semilogy(costs)

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2017-03-02
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多