import numpy as np
import random, math
from scipy.optimize import minimize
import matplotlib.pyplot as plt

class Kernel:
    def __init__(self, ker, ps):
        self.ker = ker
        self.ps = ps
        return

    def __getitem__(self, input):
        if self.ker == 'linear': return np.dot(input[0],input[1])
        if self.ker == 'poly': return (np.dot(input[0],input[1])+1)**self.ps
        if self.ker == 'rbf': return math.exp(-(np.dot(np.subtract(input[0],input[1]) , np.subtract(input[0],input[1]))) / (2*self.ps)**2)


#---------------------------------------------------------------------------------

class Objective:
    def __init__(self, x, t, a, N, C=None, ker='linear', kerPs = 2):
        self.kernel = Kernel(ker, kerPs)
        p = np.zeros((N, N))
        for i in range(N):
            for j in range(N):
                p[i,j] = t[i]*t[j]*self.kernel[[x[i],x[j]]]

        self.p = p
        self.x = x
        self.t = t
        self.alpha = a
        self.N = N
        self.bounds = [(0, C) for i in range(self.N)]


    def obj1(self, a):
        val = np.dot(a, np.dot(a, self.p))/2 - sum(a)
        return val

    def zerofun(self, a):
        return np.dot(a,self.t)

    def indicator(self, s1, s2):
        s = [s1, s2]
        k = [self.kernel[[s, xi]] for xi in self.x]
        return np.dot(np.multiply(self.alpha, self.t), k) - self.b

#---------------------------------------------------------------------------------

classA = np.random.randn(20, 2)  + [0.0, 0.0]
classB = np.random.randn(20, 2)  + [0.0, 0.0]


inputs = np.concatenate((classA, classB))
targets = np.concatenate(
    (np.ones(classA.shape[0]),
    -np.ones(classB.shape[0])))

N = inputs.shape[0]

permute = list(range(N))
random.shuffle(permute)
inputs = inputs[permute, :]
targets = targets[permute]


O = Objective(inputs, targets, np.zeros(N), N, ker='rbf', kerPs = 2)

ret = minimize(O.obj1, O.alpha, bounds = O.bounds, constraints={'type':'eq', 'fun':O.zerofun})
alpha = ret['x']

res = np.array([[alpha[i], O.x[i], O.t[i]] for i in range(len(alpha)) if alpha[i] > 1e-5])
#print(res[:,0])