From fa7fa0c09e29a4ade90fb59e18d5e58675365520 Mon Sep 17 00:00:00 2001
From: ysyapa <ysyapa@ncsu.edu>
Date: Wed, 30 Aug 2023 21:23:26 +0000
Subject: [PATCH] Almost correct

---
 common.jl | 28 +++++++++++++++-------------
 1 file changed, 15 insertions(+), 13 deletions(-)

diff --git a/common.jl b/common.jl
index 6a95d8a..8be2ddf 100644
--- a/common.jl
+++ b/common.jl
@@ -9,23 +9,25 @@ struct system{T}
     n::Int
     N::Int
     L::T
-    μs::Vector{Int}
-    invU::Matrix{Int}
+    μs::Vector{T}
+    invU::Matrix{T}
 
     function system{T}(d::Int, n::Int, N::Int, L::Real) where {T<:Float}
-        μs = [Int((coord + 1)^2 * reducedMass(coord, 1)) for coord in 1:(n - 1)]
+        μs = [1/((coord + 1)^2 * reducedMass(coord, 1)) for coord in 1:(n - 1)]
 
         # TODO: Optimize
-        invU = Matrix{Int}(undef, n, n - 1)
-        for i in CartesianIndices(invU)
-            if i[1] - 1 == i[2]
-                invU[i] = -i[2]
-            elseif i[1] > i[2]
-                invU[i] = 0
+        U = Matrix{T}(undef, n, n)
+        for i in CartesianIndices(U)
+            if i[1] + 1 == i[2]
+                U[i] = -i[1]
+            elseif i[1] >= i[2]
+                U[i] = 1
             else
-                invU[i] = 1
+                U[i] = 0
             end
         end
+        U[n, :] .= 1/n
+        invU = inv(U)[:, 1:(n - 1)]
         
         return new{T}(d, n, N, convert(T, L), μs, invU)
     end
@@ -44,7 +46,7 @@ end
 which_index(s::system, dim::Int, coord::Int)::Int = (dim - 1) * (s.n - 1) + coord
 
 "Get the distance to the nearest image of the particle"
-function nearest(s::system, Δk::Int)::Int
+function nearest(s::system, Δk)
     # TODO: Optimize
     while true
         if Δk >= s.N ÷ 2
@@ -67,9 +69,9 @@ function calculate_Vs(s::system{T}, V_twobody::Function, ϕ::T, n_image::Int)::A
         rs = s.invU * xs
         for p1 in 1:s.n
             for p2 in 1:(p1 - 1)
-                Δk = Array{Int}(undef, s.d)
+                Δk = Array{T}(undef, s.d)
                 for dim in 1:s.d
-                    Δk_temp = Int(rs[p1, dim] - rs[p2, dim])
+                    Δk_temp = rs[p1, dim] - rs[p2, dim]
                     Δk[dim] = nearest(s, Δk_temp)
                 end
                 for image in images