dyna3/engine-proto/gram-test/tetrahedron-radius-ratio.jl

include("Engine.jl")

using SparseArrays
using Random

# initialize the partial gram matrix for a sphere inscribed in a regular
# tetrahedron
J = Int64[]
K = Int64[]
values = BigFloat[]
for j in 1:10
  for k in 1:10
    filled = false
    if j == 10
      if k <= 4
        push!(values, 0)
        filled = true
      end
    elseif k == 10
      if j <= 4
        push!(values, 0)
        filled = true
      end
    elseif j == k
      push!(values, j <= 5 ? 1 : 0)
      filled = true
    elseif j <= 4
      if k <= 4
        push!(values, -1/BigFloat(3))
        filled = true
      elseif k == 5
        push!(values, -1)
        filled = true
      elseif k <= 9 && k - j != 5
        push!(values, 0)
        filled = true
      end
    elseif k <= 4
      if j == 5
        push!(values, -1)
        filled = true
      elseif j <= 9 && j - k != 5
        push!(values, 0)
        filled = true
      end
    end
    if filled
      push!(J, j)
      push!(K, k)
    end
  end
end
gram = sparse(J, K, values)

# set initial guess
Random.seed!(99230)
guess = hcat(
  sqrt(1/BigFloat(3)) * BigFloat[
    1  1 -1 -1  0
    1 -1  1 -1  0
    1 -1 -1  1  0
    0  0  0  0 -1.5
    1  1  1  1 -0.5
  ] + 0.2*Engine.rand_on_shell(fill(BigFloat(-1), 5)),
  Engine.point([-1, -1, -1]),
  Engine.point([-1,  1,  1]),
  Engine.point([ 1, -1,  1]),
  Engine.point([ 1,  1, -1]),
  BigFloat[0, 0, 0, 0, 1]
)
frozen = vcat(
  [CartesianIndex(4, k) for k in 6:9],
  [CartesianIndex(j, 10) for j in 1:5]
)

# complete the gram matrix
L, success, history = Engine.realize_gram(gram, guess, frozen)
completed_gram = L'*Engine.Q*L
println("Completed Gram matrix:\n")
display(completed_gram)
if success
  println("\nTarget accuracy achieved!")
else
  println("\nFailed to reach target accuracy")
end
println("Steps: ", size(history.scaled_loss, 1))
println("Loss: ", history.scaled_loss[end], "\n")