From 16826cf07c692390f965cd2fe9de5491810caa1d Mon Sep 17 00:00:00 2001
From: Aaron Fenyes <aaron.fenyes@fareycircles.ooo>
Date: Tue, 20 Feb 2024 22:35:24 -0500
Subject: [PATCH 1/2] Try out the Gram matrix approach

---
 engine-proto/gram-test/gram-test.jl   | 30 +++++++++++++++++++++++++++
 engine-proto/gram-test/gram-test.sage | 27 ++++++++++++++++++++++++
 2 files changed, 57 insertions(+)
 create mode 100644 engine-proto/gram-test/gram-test.jl
 create mode 100644 engine-proto/gram-test/gram-test.sage

diff --git a/engine-proto/gram-test/gram-test.jl b/engine-proto/gram-test/gram-test.jl
new file mode 100644
index 0000000..b053097
--- /dev/null
+++ b/engine-proto/gram-test/gram-test.jl
@@ -0,0 +1,30 @@
+using LinearAlgebra
+using AbstractAlgebra
+
+F, (a, b, c) = rational_function_field(Generic.Rationals{BigInt}(), ["a", "b", "c"])
+M = matrix_space(F, 5, 5)
+
+# three mutually tangent spheres which are all perpendicular to the x, y plane
+gram = M(F.([
+  -1  0  0  0  0;
+   0 -1  a  b  c;
+   0  a -1  1  1;
+   0  b  1 -1  1;
+   0  c  1  1 -1;
+]))
+
+r, p, L, U = lu(gram)
+solution = transpose(p * L)
+mform = U * inv(transpose(L))
+
+concrete = [evaluate(entry, [0, 1, -3//4]) for entry in solution]
+
+std_basis = [
+  0 0 0 1  1;
+  0 0 0 1 -1;
+  1 0 0 0  0;
+  0 1 0 0  0;
+  0 0 1 0  0
+]
+std_solution = M(F.(std_basis)) * solution
+std_concrete = std_basis * concrete
\ No newline at end of file
diff --git a/engine-proto/gram-test/gram-test.sage b/engine-proto/gram-test/gram-test.sage
new file mode 100644
index 0000000..a95ce97
--- /dev/null
+++ b/engine-proto/gram-test/gram-test.sage
@@ -0,0 +1,27 @@
+F = QQ['a', 'b', 'c'].fraction_field()
+a, b, c = F.gens()
+
+# three mutually tangent spheres which are all perpendicular to the x, y plane
+gram = matrix([
+  [-1, 0,  0,  0,  0],
+  [0, -1,  a,  b,  c],
+  [0,  a, -1,  1,  1],
+  [0,  b,  1, -1,  1],
+  [0,  c,  1,  1, -1]
+])
+
+P, L, U = gram.LU()
+solution = (P * L).transpose()
+mform = U * L.transpose().inverse()
+
+concrete = solution.subs({a: 0, b: 1, c: -3/4})
+
+std_basis = matrix([
+  [0, 0, 0, 1,  1],
+  [0, 0, 0, 1, -1],
+  [1, 0, 0, 0,  0],
+  [0, 1, 0, 0,  0],
+  [0, 0, 1, 0,  0]
+])
+std_solution = std_basis * solution
+std_concrete = std_basis * concrete
\ No newline at end of file

From 717e5a6200a61db02b3cc03ae6117cf394a4ad6f Mon Sep 17 00:00:00 2001
From: Aaron Fenyes <aaron.fenyes@fareycircles.ooo>
Date: Wed, 21 Feb 2024 02:50:22 -0500
Subject: [PATCH 2/2] Extend Gram matrix automatically

The signature of the Minkowski form on the subspace spanned by the Gram
matrix should tell us what the big Gram matrix has to look like
---
 engine-proto/gram-test/gram-test.jl | 73 ++++++++++++++++++++++++-----
 1 file changed, 60 insertions(+), 13 deletions(-)

diff --git a/engine-proto/gram-test/gram-test.jl b/engine-proto/gram-test/gram-test.jl
index b053097..c962967 100644
--- a/engine-proto/gram-test/gram-test.jl
+++ b/engine-proto/gram-test/gram-test.jl
@@ -1,23 +1,56 @@
 using LinearAlgebra
 using AbstractAlgebra
 
-F, (a, b, c) = rational_function_field(Generic.Rationals{BigInt}(), ["a", "b", "c"])
-M = matrix_space(F, 5, 5)
+function printgood(msg)
+  printstyled("✓", color = :green)
+  println(" ", msg)
+end
+
+function printbad(msg)
+  printstyled("✗", color = :red)
+  println(" ", msg)
+end
+
+F, gens = rational_function_field(Generic.Rationals{BigInt}(), ["a₁", "a₂", "b₁", "b₂", "c₁", "c₂"])
+a = gens[1:2]
+b = gens[3:4]
+c = gens[5:6]
 
 # three mutually tangent spheres which are all perpendicular to the x, y plane
-gram = M(F.([
-  -1  0  0  0  0;
-   0 -1  a  b  c;
-   0  a -1  1  1;
-   0  b  1 -1  1;
-   0  c  1  1 -1;
+gram = [
+  -1  1  1;
+   1 -1  1;
+   1  1 -1
+]
+
+eig = eigen(gram)
+n_pos = count(eig.values .> 0.5)
+n_neg = count(eig.values .< -0.5)
+if n_pos + n_neg == size(gram, 1)
+  printgood("Non-degenerate subspace")
+else
+  printbad("Degenerate subspace")
+end
+sig_rem = Int64[ones(2-n_pos); -ones(3-n_neg)]
+unk = hcat(a, b, c)
+M = matrix_space(F, 5, 5)
+big_gram = M(F.([
+  diagm(sig_rem) unk;
+  transpose(unk) gram
 ]))
 
-r, p, L, U = lu(gram)
-solution = transpose(p * L)
-mform = U * inv(transpose(L))
+r, p, L, U = lu(big_gram)
+if isone(p)
+  printgood("Found a solution")
+else
+  printbad("Didn't find a solution")
+end
+solution = transpose(L)
+mform = U * inv(solution)
 
-concrete = [evaluate(entry, [0, 1, -3//4]) for entry in solution]
+vals = [0, 0, 0, 1, 0, -3//4]
+solution_ex = [evaluate(entry, vals) for entry in solution]
+mform_ex = [evaluate(entry, vals) for entry in mform]
 
 std_basis = [
   0 0 0 1  1;
@@ -27,4 +60,18 @@ std_basis = [
   0 0 1 0  0
 ]
 std_solution = M(F.(std_basis)) * solution
-std_concrete = std_basis * concrete
\ No newline at end of file
+std_solution_ex = std_basis * solution_ex
+
+println("Minkowski form:")
+display(mform_ex)
+
+big_gram_recovered = transpose(solution_ex) * mform_ex * solution_ex
+valid = all(iszero.(
+  [evaluate(entry, vals) for entry in big_gram] - big_gram_recovered
+))
+if valid
+  printgood("Recovered Gram matrix:")
+else
+  printbad("Didn't recover Gram matrix. Instead, got:")
+end
+display(big_gram_recovered)
\ No newline at end of file