var('a,b,c,d')
vs = [(0, a), (1, b), (2, c), (3, d)]
def new():
    m = matrix(QQ, [[1/2, 1/2, 0, 0], [0, 0, 1/2, 1/2], [1/2, 1/2, 0, 0], [0, 0, 1/2, 1/2]])
    v = vector(QQ, [0, 0, 0, 0])
    return [m, v]
def bad(mv, i):
    mv[0][i // 2, i % 4] = 0
    return mv
def good(mv, i):
    bad(mv, i)
    mv[1][i // 2] += 1/2
    return mv
def make(mv):
    return [sum([mv[0][i, j] * v for (j, v) in vs]) + mv[1][i] == w for (i, w) in vs]
def resolve(i, j):
    return solve(make(bad(good(new(), j), i)), a, b, c, d)

all([len(resolve(i, j)) == 1 for i in range(8) for j in range(8)])

all([[resolve(i, j)[0][k].lhs() - v for (k, v) in vs] == [0, 0, 0, 0] for i in range(8) for j in range(8)])

[max([(i, j, sum(resolve(i, j)[0][k].rhs() for k in range(4)) / 4) for j in range(8) if i != j], key=lambda x: x[2]) for i in range(8)]