-- Extended Euclidean Algorithm
local function extended_gcd(a, b)
local old_r, r = a, b
local old_s, s = 1, 0
local old_t, t = 0, 1
while math.abs(r) > 1e-10 do
local quotient = math.floor(old_r / r)
old_r, r = r, old_r - quotient * r
old_s, s = s, old_s - quotient * s
old_t, t = t, old_t - quotient * t
end
return old_r, old_s, old_t
end
-- Format time into years, days, hours, minutes, seconds
local function format_time(seconds)
local years = math.floor(seconds / (365.25 * 24 * 3600))
seconds = seconds % (365.25 * 24 * 3600)
local days = math.floor(seconds / (24 * 3600))
seconds = seconds % (24 * 3600)
local hours = math.floor(seconds / 3600)
seconds = seconds % 3600
local minutes = math.floor(seconds / 60)
seconds = math.floor(seconds % 60)
return string.format("%d years, %d days, %d hours, %d minutes, %d seconds",
years, days, hours, minutes, seconds)
end
-- Find meeting time between two orbiting objects
function find_meeting_time(t1, t2, T1, T2)
-- Normalize times to be within one period
t1 = t1 % T1
t2 = t2 % T2
-- Calculate time difference
local delta_t = t2 - t1
-- Special handling for very close periods
if math.abs(T1 - T2) < 1 then
-- Calculate how many periods needed to make up the position difference
local period_diff = math.abs(1/T1 - 1/T2) -- Difference in frequencies
if period_diff < 1e-10 then
return nil -- Periods too close, effectively equal
end
local meeting_time = math.abs(delta_t / period_diff)
-- Ensure we get the first future meeting time
if meeting_time < t1 then
-- Calculate approximate LCM for close periods
local lcm = (T1 * T2) / math.min(T1, T2)
meeting_time = meeting_time + lcm
end
return math.max(meeting_time, t1)
end
-- Regular case for different periods
local gcd, x, y = extended_gcd(T1, T2)
-- Check if meeting is possible
if math.abs(delta_t % gcd) > 1e-6 then
return nil
end
-- Calculate number of orbits needed
local n = (delta_t * x) / gcd
local k = math.ceil(-n * gcd / T2)
n = n + k * (T2 / gcd)
-- Calculate meeting time
local meeting_time = t1 + n * T1
-- Ensure positive future time
while meeting_time <= 0 do
local lcm = (T1 * T2) / gcd
meeting_time = meeting_time + lcm
end
return meeting_time
end
-- Timing function
local function time_execution(func, ...)
local start = os.clock()
local result = func(...)
local end_time = os.clock()
return result, end_time - start
end
-- Example usage with timing:
local t1 = 0 -- time to intersection for orbit 1
local t2 = 2500 -- time to intersection for orbit 2
local T1 = 5000 -- period of orbit 1
local T2 = 5000.1 -- period of orbit 2
-- Run the calculation multiple times to get an average
local num_runs = 1000
local total_time = 0
for i = 1, num_runs do
local result, execution_time = time_execution(find_meeting_time, t1, t2, T1, T2)
total_time = total_time + execution_time
-- Only print the result for the first run
if i == 1 then
if result then
print(string.format("Objects will meet after %.2f seconds", result))
print(format_time(result))
else
print("Objects will never meet at the intersection point")
end
end
end
print(string.format("\nAverage execution time over %d runs: %.6f seconds", num_runs, total_time/num_runs))
print(string.format("Total time for %d runs: %.6f seconds", num_runs, total_time))
