import numpy as np
import cvxpy as cp


def svm_primal_optimization(X, y, C):

    m=X.shape[0]
    n=X.shape[1]

    w = cp.Variable(n)
    b = cp.Variable(1)
    xi = cp.Variable(m)

    obj = 0.5*cp.norm(w)**2 + C*cp.sum(xi)
    constr=[]
    for i in range(m):
        constr+=[y[i]*(w@X[i,:] + b)>=1-xi[i]]
        constr+=[xi[i]>=0]

    prob = cp.Problem(cp.Minimize(obj), constr)
    prob.solve(solver="MOSEK", verbose=False)

    wstar = w.value
    bstar = b.value
    xistar = xi.value

    return wstar, bstar, xistar

def svm_dual_optimization(X, y, C, kernel, d, coef0):
    m=X.shape[0]
    n=X.shape[1]

    if kernel =='linear':
        K = np.zeros((m,m))
        for i in range(m):
            for j in range(m):
                K[i,j] = y[i]*y[j]*X[i,:]@X[j,:]
    elif kernel =='polynomial':
        K = poly_kernel(X, y, d, coef0)

    alpha = cp.Variable(m)

    obj = cp.sum(alpha) -0.5*cp.quad_form(alpha, cp.psd_wrap(K))

    constr=[alpha @ y==0]
    for i in range(m):
        constr+=[alpha[i]>=0]
        constr+=[alpha[i]<=C]

    prob = cp.Problem(cp.Maximize(obj), constr)
    prob.solve(solver="MOSEK", verbose=True)

    alphastar=alpha.value

    return alphastar


def poly_kernel(X, y, d, c):
    m=X.shape[0]
    n=X.shape[1]

    K = np.zeros((m,m))
    for i in range(m):
        for j in range(m):
            K[i,j] = y[i]*y[j]*((X[i,:]@X[j,:] + c)**d)

    return K