Convex Inner Problem For SVI (As described by Zeliade)

convex_inner_svi.py

def optimize(y,max_vi,market_vols,tte,s,weights):
    y = y[market_vols>0]
    vols = market_vols[market_vols>0]
    w = weights[market_vols>0]

    n = len(y)
        
    R = np.concatenate([np.ones(n), np.sqrt(y**2+1), y]).reshape(3,n).T
    W = np.diag(w)
    
    var_t = vols**2*tte
    max_vi = np.nanmax(var_t)
    
    Q = 2*matrix(R.T.dot(W).dot(R))


    Y1  = np.sum(y)
    Y2  = np.sum(y**2)
    Y3  = np.sum(np.sqrt(y**2+1))
    Y4  = np.sum(y*np.sqrt(y**2+1))
    Y5  = np.sum(y**2+1)

    v   = np.sum(var_t)
    vsq = np.sum(var_t*var_t)
    vY  = np.sum(var_t*y)
    vY2 = np.sum(var_t*np.sqrt(y**2+1))

    q = -2*matrix(var_t.T.dot(W).dot(R))

    h = matrix([0.0, 4.0*s, 0.0, 0.0, 4.0*s, 4.0*s, 0.0, max_vi])

 
    # G is in column-major order
    G = matrix([[0.0,0.0,0.0,0.0,0.0,0.0,-1.0,1.0],[-1.0,1.0,-1.0,-1.0,1.0,1.0,0.0,0.0],[0.0,0.0,-1.0,1.0,-1.0,1.0,0.0,0.0]])

    cvxopt.solvers.options['show_progress'] = False
    sol = solvers.qp(Q,q,G,h)
 
    # I convert back to the notation in the paper here
    # to verify that the constraints are satisfied:
    a = sol['x'][0]
    c = sol['x'][1]
    d = sol['x'][2]
    
    if(c < 0 or c > 4*s):
        return [a,c,d,-1] 
    if(math.fabs(d) > c or math.fabs(d) > ((4*s) -c)):
        return [a,c,d,-1] 
    if(a < 0.0 or a > max_vi):
        return [a,c,d,-1] 
    #1 indicates sucess, -1 failure (in case you want to do something about that...)
    return [a,c,d,1]