import numpy as np

def dual_upper(beta, A, B, b, c, d, l, r):
    # if L(x, lambda1, beta * lambda1) =
    # A: (m,m) b: (m), c: ()
    # return L(beta)
    AA = A + beta * B
    AAinvb = np.linalg.solve(AA, b)
    AAinvc = np.linalg.solve(AA, c)
    sqrterm = np.sum(AAinvc * c) * (np.sum(AAinvb * b) * 0.25 - l + d - beta * r * r)
    return 0.5 * (np.sqrt(sqrterm) - np.sum(AAinvc * b))

def finite_random_feature_Q(x, s, a):
    return

def infinite_Q(x, s, a):
    return

def finite_random_feature(SASR, w, phi, Ktilde):
    return

def infinite_dim_RKHS(lambda1, lambda2, REW, K, Ktilde, r, l):
    # K : (n+1, n+1)
    # Ktilde : (n,n)
    # REW : (n)
    n = REW.shape[0]
    k0 = K[:, 0]    # (n+1)
    K1 = K[:, 1:]   # (n+1, n)
    K1Ktilde = np.matmul(K1, Ktilde)    # (n+1, n)
    A = -lambda1 / n / n * np.matmul(K1Ktilde,K1.T) - lambda2 * K   # (n+1, n+1)
    b = k0 + 2 * lambda1 / n / n * np.matmul(K1Ktilde, REW)         # (n+1)
    Ktilder = np.matmul(Ktilde, REW)                                # (n)
    c = lambda2 * r * r + lambda1 * l - lambda1 /n /n np.sum(REW * Ktilder) # ()
