Prove Braikenridge–Maclaurin theorem by analytic geometry

Question

Braikenridge–Maclaurin theorem is the converse to Pascal's theorem.

I put point I onto the origin of Cartesian coordinates, and rotate GIH onto y-axis (set GI = g and IH = -h), then get

$\begin{cases}AB:y=jx+g\\DE:y=kx+g\\BC:y=mx+h\\EF:y=nx+h\\AF:y=px\\CD:y=qx\end{cases}$

My task is to prove 6 points

$\begin{cases}A \left( \frac{g}{- j + p}, \ \frac{g p}{- j + p}\right)\\ B \left( \frac{- g + h}{j - m}, \ \frac{- g m + h j}{j - m}\right)\\ C \left( \frac{h}{- m + q}, \ \frac{h q}{- m + q}\right)\\ D \left( \frac{g}{- k + q}, \ \frac{g q}{- k + q}\right)\\ E \left( \frac{- g + h}{k - n}, \ \frac{- g n + h k}{k - n}\right)\\ F \left( \frac{h}{- n + p}, \ \frac{h p}{- n + p}\right)\end{cases}$

lie on a conic.

Here I assume the conic doesn't go through origin I, so these points match the conic

$ax^2+bxy+cy^2+dx+ey+1=0$

According to this rule, the rest thing is to prove

$\det\left[\begin{matrix}\frac{g^{2}}{\left(- j + p\right)^{2}} & \frac{g^{2} p}{\left(- j + p\right)^{2}} & \frac{g^{2} p^{2}}{\left(- j + p\right)^{2}} & \frac{g}{- j + p} & \frac{g p}{- j + p} & 1\\\frac{\left(- g + h\right)^{2}}{\left(j - m\right)^{2}} & \frac{\left(- g + h\right) \left(- g m + h j\right)}{\left(j - m\right)^{2}} & \frac{\left(- g m + h j\right)^{2}}{\left(j - m\right)^{2}} & \frac{- g + h}{j - m} & \frac{- g m + h j}{j - m} & 1\\\frac{h^{2}}{\left(- m + q\right)^{2}} & \frac{h^{2} q}{\left(- m + q\right)^{2}} & \frac{h^{2} q^{2}}{\left(- m + q\right)^{2}} & \frac{h}{- m + q} & \frac{h q}{- m + q} & 1\\\frac{g^{2}}{\left(- k + q\right)^{2}} & \frac{g^{2} q}{\left(- k + q\right)^{2}} & \frac{g^{2} q^{2}}{\left(- k + q\right)^{2}} & \frac{g}{- k + q} & \frac{g q}{- k + q} & 1\\\frac{\left(- g + h\right)^{2}}{\left(k - n\right)^{2}} & \frac{\left(- g + h\right) \left(- g n + h k\right)}{\left(k - n\right)^{2}} & \frac{\left(- g n + h k\right)^{2}}{\left(k - n\right)^{2}} & \frac{- g + h}{k - n} & \frac{- g n + h k}{k - n} & 1\\\frac{h^{2}}{\left(- n + p\right)^{2}} & \frac{h^{2} p}{\left(- n + p\right)^{2}} & \frac{h^{2} p^{2}}{\left(- n + p\right)^{2}} & \frac{h}{- n + p} & \frac{h p}{- n + p} & 1\end{matrix}\right]=0$

It doesn't look too complicated. I use SymPy to do this task but it is running 6 hours and still doesn't return.

Here is my code:

from sympy import *

def L(a, b):
    return Eq(y, a * x + b)

def intersect(L1, L2):
    P = solve([L1, L2], (x, y))
    return simplify(P[x]), simplify(P[y])

g, h, j, k, m, n, p, q, x, y = symbols('g, h, j, k, m, n, p, q, x, y')
AB, DE, BC, EF, AF, CD = L(j, g), L(k, g), L(m, h), L(n, h), L(p, 0), L(q, 0)
A, E, C = intersect(AB, AF), intersect(DE, EF), intersect(BC, CD)
D, B, F = intersect(CD, DE), intersect(AB, BC), intersect(AF, EF)

mat = []
for P in [A, B, C, D, E, F]:
    x, y = P[0], P[1]
    mat.append([x * x, x * y, y * y, x, y, 1])
print(Matrix(mat).det())

Are there any ways to simplify this matrix to get the result quickly?

UPDATE

Applied with brainjam's solution and the trick, it can be run out quickly. Here is the improved code:

from sympy import *

def L(a, b):
    return Eq(y, a * x + b)

def intersect(L1, L2):
    P = solve([L1, L2], (x, y))
    return simplify(P[x]), simplify(P[y])

g, h, j, k, m, n, p, q, x, y = symbols('g, h, j, k, m, n, p, q, x, y')
AB, DE, BC, EF, AF, CD = L(j, g), L(k, g), L(m, h), L(n, h), L(p, 0), L(q, 0)
A, E, C = intersect(AB, AF), intersect(DE, EF), intersect(BC, CD)
D, B, F = intersect(CD, DE), intersect(AB, BC), intersect(AF, EF)
print('A:', A)
print('B:', B)
print('C:', C)
print('D:', D)
print('E:', E)
print('F:', F)

points = [A, B, C, D, E, F]
denominators = [j - p, j - m, q - m, q - k, k - n, p - n]
coefficients = [x * x, x * y, y * y, x, y, 1]
subs = {}
mat = []
for s in range(6):
    row = []
    denominator = denominators[s] ** 2
    subs_xy = {}
    subs_xy[x] = points[s][0]
    subs_xy[y] = points[s][1]
    for t in range(6):
        rst = symbols('r' + str(s) + str(t))
        subs[rst] = expand(simplify(N(coefficients[t] * denominator, subs = subs_xy)))
        row.append(rst)
    mat.append(row)
print('M =', N(Matrix(mat), subs = subs))
print('det M =', expand(N(Matrix(mat).det(), subs = subs)))

Answer 1

I found simply Matrix(mat).det(method='domain-ge') works. See here for details.

Full code:

from sympy import *

def L(a, b):
    return Eq(y, a * x + b)

def intersect(L1, L2):
    P = solve([L1, L2], (x, y))
    return simplify(P[x]), simplify(P[y])

g, h, j, k, m, n, p, q, x, y = symbols('g, h, j, k, m, n, p, q, x, y')
AB, DE, BC, EF, AF, CD = L(j, g), L(k, g), L(m, h), L(n, h), L(p, 0), L(q, 0)
A, E, C = intersect(AB, AF), intersect(DE, EF), intersect(BC, CD)
D, B, F = intersect(CD, DE), intersect(AB, BC), intersect(AF, EF)

mat = []
for P in [A, B, C, D, E, F]:
    x, y = P[0], P[1]
    mat.append([x*x, x*y, y*y, x, y, 1])
print('M =', Matrix(mat))
print('det M =', Matrix(mat).det(method='domain-ge'))
```

Answer 2

If you search for information like "sympy determinant slow" you get lots of hits. The gist of the problem is that the det() function seems to be trying to simplify as it goes, and can get bogged down once the matrices are $4\times 4$ or larger.

One problem with your code is that you have lots of rational expressions, and this is probably making things worse.

Since you just want to confirm that the determinant of your matrix is zero, you can multiply each row by its least common denominator. So, for example, the first row

$$ \frac{g^{2}}{\left(- j + p\right)^{2}} \quad \frac{g^{2} p}{\left(- j + p\right)^{2}} \quad \frac{g^{2} p^{2}}{\left(- j + p\right)^{2}} \quad \frac{g}{- j + p} \quad \frac{g p}{- j + p} \quad 1 $$

becomes

$$ g^{2} \quad g^{2} p \quad g^{2} p^{2} \quad g(- j + p) \quad g p(- j + p) \quad \left(- j + p\right)^{2} $$

With this modification I ran your code on replit.com and got the desired result (zero) after 15 minutes.

This is still a long time (but better than 6 hours or never). You could try something like this trick to eliminate the simplify-as-you-go problem to get a raw determinant, and then simplify the raw determinant.