Быстрые решатели CVX в Matlab

Question

Мне интересно, что является самым быстрым выпуклым оптимизатором в Matlab или есть ли способ ускорить текущие решатели? Я использую CVX, но настойчиво решает проблему оптимизации, которую у меня есть. оптимизация у меня есть, чтобы решить

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

где размер А и В очень большой.

Есть ли способ решить это методом наименьшего квадрата, а затем передать его в версию ограничения, чтобы сделать ее быстрее?

Answer 1

Я не уверен, что означает x.d <= delta, но я просто предполагаю, что это должно быть x <= delta.

Вы можете решить эту проблему, используя проецируемый метод градиента или метод ускоренного проецируемого градиента (который представляет собой лишь небольшую модификацию проецируемого градиентного метода, который «магически» сходится намного быстрее). Вот некоторый код python, который показывает, как свести к минимуму .5 || Ax - b ||^2 при условии, что 0 < = x < = delta с использованием FISTA, который является ускоренным проецируемым градиентом. Более подробную информацию о проецированном градиентном методе и FISTA можно найти, например, у Boyd's 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)