diff --git a/.gitignore b/.gitignore index f8bee76..5e5df5b 100644 --- a/.gitignore +++ b/.gitignore @@ -1,3 +1,4 @@ pdf/ **/*.pdf julia/**/*.html +julia/test/ipopt.out diff --git a/julia/src/Thesis.jl b/julia/src/Thesis.jl index fd82220..7ca2c0d 100644 --- a/julia/src/Thesis.jl +++ b/julia/src/Thesis.jl @@ -16,6 +16,7 @@ module Thesis end include("./errors.jl") + include("./bodies.jl") include("./constants.jl") include("./spacecraft.jl") include("./mission.jl") diff --git a/julia/src/bodies.jl b/julia/src/bodies.jl new file mode 100644 index 0000000..6cf6f48 --- /dev/null +++ b/julia/src/bodies.jl @@ -0,0 +1,17 @@ +export state + +struct Body + μ::Float64 + r::Float64 # radius + color::String + id::Int # SPICE id +end + +function state(p::Body, t::DateTime, add_v∞::Vector{Float64}=[0., 0., 0.], add_mass::Float64=1e10) + time = utc2et(Dates.format(t,"yyyy-mm-ddTHH:MM:SS")) + [ spkssb(p.id, time, "ECLIPJ2000"); 0.0 ] + [ zeros(3); add_v∞; add_mass ] +end + +function state(p::Body, t::Float64, add_v∞::Vector{Float64}=[0., 0., 0.], add_mass::Float64=1e10) + [ spkssb(p.id, t, "ECLIPJ2000"); 0.0 ] + [ zeros(3); add_v∞; add_mass ] +end diff --git a/julia/src/constants.jl b/julia/src/constants.jl index 3156e6d..7f07e85 100644 --- a/julia/src/constants.jl +++ b/julia/src/constants.jl @@ -7,13 +7,6 @@ export Jupiter, Saturn, Uranus, Neptune, Pluto export G, AU, init_STM, hour, day, year, second export Pathlist -struct Body - μ::Float64 - r::Float64 # radius - color::String - id::Int # SPICE id -end - const Sun = Body(1.32712440018e11, 696000., "Electric", 10) const Mercury = Body(2.2032e4, 2439., "heat", 1) const Venus = Body(3.257e5, 6052., "turbid", 2) diff --git a/julia/src/errors.jl b/julia/src/errors.jl index 5949997..4b0a6b6 100644 --- a/julia/src/errors.jl +++ b/julia/src/errors.jl @@ -2,12 +2,21 @@ struct LaGuerreConway_Error <: Exception end struct ΔVsize_Error <: Exception end +struct HitPlanet_Error <: Exception end +Base.showerror(io::IO, e::HitPlanet_Error) = print(io, "spacecraft hit the planet...") + struct Convergence_Error <: Exception end Base.showerror(io::IO, e::Convergence_Error) = print(io, "NLP solver didn't converge...") struct GenOrbit_Error <: Exception end Base.showerror(io::IO, e::GenOrbit_Error) = print(io, "infinite Loop trying to generate the init orbit") +struct V∞_Error <: Exception + v∞_in::AbstractVector + v∞_out::AbstractVector +end +Base.showerror(io::IO, e::V∞_Error) = print(io, "v∞s weren't the same: ", e.v∞_in, " and ", e.v∞_out) + struct PropOne_Error <: Exception ΔV_unit::AbstractVector end diff --git a/julia/src/inner_loop/monotonic_basin_hopping.jl b/julia/src/inner_loop/monotonic_basin_hopping.jl index 2995c07..0eb4478 100644 --- a/julia/src/inner_loop/monotonic_basin_hopping.jl +++ b/julia/src/inner_loop/monotonic_basin_hopping.jl @@ -41,30 +41,35 @@ function gen_date(date_range::Tuple{DateTime, DateTime}) l0 + Millisecond(floor(rand()*(lf-l0).value)) end -""" -Returns a random amount of time in a range -""" -function gen_period(date_range::Tuple{DateTime, DateTime}) - l0, lf = date_range - Millisecond(floor(rand()*(lf-l0).value)) -end - """ Perturbs a valid mission with pareto-distributed variables, generating a mission guess """ function perturb(mission::Mission) - new_launch_date = mission.launch_date + Dates.Second(floor(7day * pareto_add())) - new_launch_v∞ = mission.launch_v∞ .* pareto(3) - new_phases = Vector{Phase}() - for phase in mission.phases - new_v∞_in = phase.v∞_in .* pareto(3) - new_δ = phase.δ * pareto() - new_tof = phase.tof * pareto() - new_thrust_profile = phase.thrust_profile .* pareto(size(phase.thrust_profile)) - push!(new_phases, Phase(phase.planet, new_v∞_in, new_δ, new_tof, new_thrust_profile)) + mission_guess = Bad_Mission("Starting point") + while typeof(mission_guess) == Bad_Mission + new_launch_date = mission.launch_date + Dates.Second(floor(7day * pareto_add())) + new_launch_v∞ = mission.launch_v∞ .* pareto(3) + new_phases = Vector{Phase}() + for phase in mission.phases + new_v∞_in = phase.v∞_in .* pareto(3) + temp_v∞_out = phase.v∞_out .* pareto(3) + new_v∞_out = norm(new_v∞_in) * temp_v∞_out/norm(temp_v∞_out) + new_tof = phase.tof * pareto() + new_thrust_profile = phase.thrust_profile .* pareto(size(phase.thrust_profile)) + push!(new_phases, Phase(phase.planet, new_v∞_in, new_v∞_out, new_tof, new_thrust_profile)) + end + # TODO: Mission_Guess.validate() + try + mission_guess = Mission_Guess(mission.sc, mission.start_mass, new_launch_date, new_launch_v∞, new_phases) + catch e + if isa(e, Mass_Error) ||isa(e,V∞_Error) || isa(e,HitPlanet_Error) + continue + else + rethrow() + end + end end - # TODO: Mission_Guess.validate() - Mission_Guess(mission.sc, mission.start_mass, new_launch_date, new_launch_v∞, new_phases) + return mission_guess end """ @@ -78,66 +83,63 @@ function mission_guess( flybys::Vector{Body}, max_v∞_in_mag::Float64, latest_arrival::DateTime, primary::Body=Sun ) - # TODO: Eventually I can calculate n more intelligently - n = 20 + mission_guess = Bad_Mission("Keep trying to generate a guess") + while mission_guess == Bad_Mission("Keep trying to generate a guess") + # TODO: Eventually I can calculate n more intelligently + n = 20 - # Determine the launch conditions - launch_date = gen_date(launch_window) - launch_v∞_normalized = rand(-1:0.0001:1, 3) - launch_v∞ = rand(0:0.0001:√max_C3_out) * launch_v∞_normalized/norm(launch_v∞_normalized) + # Determine the launch conditions + launch_date = gen_date(launch_window) + launch_v∞_normalized = rand(-1:0.0001:1, 3) + launch_v∞ = rand(0:0.0001:√max_C3_out) * launch_v∞_normalized/norm(launch_v∞_normalized) - # Determine the leg lengths - num_phases = length(flybys) - 1 - total_tof = 100year - max_tof = (latest_arrival - launch_date).value / 1000 - tofs = rand(30day : hour : 0.7max_tof, num_phases) - total_tof = sum(tofs) - while total_tof > max_tof + # Determine the leg lengths + num_phases = length(flybys) - 1 + total_tof = 100year + max_tof = (latest_arrival - launch_date).value / 1000 tofs = rand(30day : hour : 0.7max_tof, num_phases) total_tof = sum(tofs) - end - - # For each phase, determine the v∞_in and δ - phases = Vector{Phase}() - for i in 1:num_phases - flyby = flybys[i+1] - v∞_normalized = rand(-1:0.0001:1, 3) - v∞ = rand(0:0.0001:10) * v∞_normalized/norm(v∞_normalized) - δ = rand(0:0.0001:2π) - periapsis = (flyby.μ/(v∞ ⋅ v∞)) * ( 1/sin(δ/2) - 1 ) - while periapsis < flyby.r + 100. - δ = rand(0:0.0001:2π) - periapsis = (flyby.μ/(v∞ ⋅ v∞)) * ( 1/sin(δ/2) - 1 ) + while total_tof > max_tof + tofs = rand(30day : hour : 0.7max_tof, num_phases) + total_tof = sum(tofs) + end + + # For each phase, determine the v∞_in and δ + phases = Vector{Phase}() + for i in 1:num_phases + flyby = flybys[i+1] + v∞_in_normalized = rand(-1:0.0001:1, 3) + v∞_in = rand(0:0.0001:10) * v∞_in_normalized/norm(v∞_in_normalized) + v∞_out_normalized = rand(-1:0.0001:1, 3) + v∞_out = norm(v∞_in) * v∞_out_normalized/norm(v∞_out_normalized) + δ = acos( ( v∞_in ⋅ v∞_out ) / ( norm(v∞_in) * norm(v∞_out) ) ) + periapsis = (flyby.μ/(v∞_in ⋅ v∞_in)) * ( 1/sin(δ/2) - 1 ) + while periapsis < flyby.r + 100. + v∞_out_normalized = rand(-1:0.0001:1, 3) + v∞_out = norm(v∞_in) * v∞_out_normalized/norm(v∞_out_normalized) + δ = acos( ( v∞_in ⋅ v∞_out ) / ( norm(v∞_in) * norm(v∞_out) ) ) + periapsis = (flyby.μ/(v∞_in ⋅ v∞_in)) * ( 1/sin(δ/2) - 1 ) + end + thrusts = rand(-1:0.0001:1,(n,3)) + push!(phases, Phase(flyby, v∞_in, v∞_out, tofs[i], thrusts)) + end + + # Finally, determine the arrival v∞ + arrival_v∞_normalized = rand(-1:0.0001:1, 3) + arrival_v∞ = rand(0:0.0001:max_v∞_in_mag) * arrival_v∞_normalized/norm(arrival_v∞_normalized) + + # And we can construct a mission guess object with these values + try + mission_guess = Mission_Guess( sc, start_mass, launch_date, launch_v∞, phases ) + catch e + if isa(e, Mass_Error) ||isa(e,V∞_Error) || isa(e,HitPlanet_Error) + continue + else + rethrow() + end end - thrusts = rand(-1:0.0001:1,(n,3)) - push!(phases, Phase(flyby, v∞, δ, tofs[i], thrusts)) end - - # Finally, determine the arrival v∞ - arrival_v∞_normalized = rand(-1:0.0001:1, 3) - arrival_v∞ = rand(0:0.0001:max_v∞_in_mag) * arrival_v∞_normalized/norm(arrival_v∞_normalized) - - # And we can construct a mission guess object with these values - Mission_Guess( sc, start_mass, launch_date, launch_v∞, phases ) -end - -""" -A convenience function for calculating mass usage given a certain thrust profile -""" -function mass_consumption(sc::Sc, phase::Phase) - weighted_thrusting_time = 0.0 - n = size(phase.thrust_profile)[1] - for i in 1:size(phase.thrust_profile,1) - weighted_thrusting_time += norm(phase.thrust_profile[i,:]) * phase.tof/n - end - return weighted_thrusting_time*sc.mass_flow_rate -end - -""" -This attempts to determine v∞_out from v∞_in and the turning angle, assuming we are staying in the -plane of the three planets -""" -function calc_turn(p1::Body, p2::Body, p3::Body, v∞_in::Vector{Float64}, δ::Float64) + return mission_guess end """ @@ -145,25 +147,26 @@ Sequentially calls the NLP solver to attempt to solve based on the initial guess """ function inner_loop_solve(guess::Mission_Guess) v∞_out = guess.launch_v∞ - time = utc2et(Dates.format(guess.launch_date,"yyyy-mm-ddTHH:MM:SS")) + mass = guess.start_mass current_planet = Earth - start = [spkssb(Earth.id, time, "ECLIPJ2000"); 0.0] + [ zeros(3); guess.launch_v∞; guess.start_mass ] + time = utc2et(Dates.format(guess.launch_date,"yyyy-mm-ddTHH:MM:SS")) + start = state(current_planet, time, v∞_out, mass) corrected_phases = Vector{Phase}() for i in 1:length(guess.phases) phase = guess.phases[i] + current_planet = phase.planet time += phase.tof - goal = spkssb(phase.planet.id, time, "ECLIPJ2000") + [zeros(3); phase.v∞_in] + goal = state(current_planet, time, phase.v∞_in) result = solve_phase( start, goal, guess.sc, phase.tof, phase.thrust_profile) result.converged || return Bad_Mission(result.info) # Drop if it's not working - corrected_phase = Phase(phase.planet, phase.v∞_in, phase.δ, phase.tof, result.sol) + corrected_phase = Phase(phase.planet, phase.v∞_in, phase.v∞_out, phase.tof, result.sol) push!(corrected_phases, corrected_phase) mass_used = mass_consumption(guess.sc, corrected_phase) + mass -= mass_used if i != length(guess.phases) - v∞_out = calc_turn(current_planet, phase.planet, guess.phases[i+1].planet, phase.v∞_in, phase.δ) - current_planet = phase.planet - planet_state = [spkssb(current_planet.id, time, "ECLIPJ2000"); 0.0] - start = planet_state + [ zeros(3); v∞_out; start_mass - mass_used ] + v∞_out = phase.v∞_out + start = state(current_planet, time, v∞_out, mass) end end @@ -183,8 +186,13 @@ function cost(mission::Mission, max_C3::Float64, max_v∞::Float64) return 3mass_percent + C3_percent + v∞_percent end +cost(_::Bad_Mission) = Inf - +""" +This is the main monotonic basin hopping function. There's a lot going on here, but the general idea +is that hopefully you can provide mission parameters and a list of flybys and get the optimal +mission that uses those flybys. +""" function mbh( flybys::Vector{Body}, sc::Sc, start_mass::Float64, @@ -193,49 +201,64 @@ function mbh( flybys::Vector{Body}, max_v∞::Float64, latest_arrival::DateTime, primary::Body=Sun; - search_patience::Int=1_000, - drill_patience::Int=50) - # Let's pseudo-code this bitch - # - # First, we need a function (mission_guess below) that randomly generates the: - # - Launch date - # - Launch v∞_out - # - for each phase: - # - tof - # - v∞_in - # - turning angle - # - # Also need a function (inner_loop_solve below) that takes the generated decision vector guess and - # attempts to correct it with the NLP solver. It will either converge or not. - # - # Also need a costing function (may be provided potentially) - # - # Also need a perturb function - # - guess = mission_guess(flybys, sc, start_mass, launch_window, max_C3, max_v∞, latest_arrival, primary) - inner_loop_solve(guess) - # search_count = 0 - # x_current = "Bad" - # while search_count < search_patience - # search_count += 1 - # drill_count = 0 - # x_star = inner_loop_solve(mission_guess()) - # if x_star.converged - # if cost(x_star) < cost(x_current) - # x_current = x_star - # while drill_count < drill_patience - # x_star = inner_loop_solve(perturb(x_current)) - # if x_star.converged and cost(x_star) < cost(x_current) - # x_current = x_star - # drill_count = 0 - # else - # drill_count += 1 - # end - # end - # push!(archive, x_current) - # end - # end - # - # + search_patience::Int=10_000, + drill_patience::Int=50, + verbose::Bool=false) + + # A convenience function + random_guess() = mission_guess(flybys,sc,start_mass,launch_window,max_C3,max_v∞,latest_arrival,primary) + cost(m::Mission) = cost(m, max_C3, max_v∞) + status(s,d,l) = print("\r\t", "search: ", s, " drill: ", d, " archive length: ", l, " ") + + # Initialize stuff + search_count = 0 + x_current = Bad_Mission("Starting point") + archive = Vector{Mission}() + + # The main loop + if verbose println("Starting Monotonic Basin Hopper") end + while search_count < search_patience + + # Intialize an x_star, if it doesn't converge, hop on to the next basin + search_count += 1 + drill_count = 0 + if verbose status(search_count, drill_count, length(archive)) end + x_star = inner_loop_solve(random_guess()) + x_star.converged || continue + x_basin = x_star + + + # If it does, though, we check to see if it's better than the current + if cost(x_star) < cost(x_current) x_current = x_star end + + # Either way, we need to drill into this particular basin, since it's valid + while drill_count < drill_patience + + if verbose status(search_count, drill_count, length(archive)) end + + #Perturb to generate a slightly different x_star + x_star = inner_loop_solve(perturb(x_basin)) + + # If better than the best, then keep it as current + if x_star.converged && cost(x_star) < cost(x_current) + x_current = x_star + x_basin = x_star + drill_count = 0 + # If better than the best in this particular basin, keep it as basin + elseif x_star.converged && cost(x_star) < cost(x_basin) + x_basin = x_star + drill_count = 0 + # Otherwise, keep drilling + else + drill_count += 1 + end + end + + x_current in archive || push!(archive, x_current) + + end + + return x_current, archive + end diff --git a/julia/src/inner_loop/phase.jl b/julia/src/inner_loop/phase.jl index b4701c7..09f9aad 100644 --- a/julia/src/inner_loop/phase.jl +++ b/julia/src/inner_loop/phase.jl @@ -39,15 +39,17 @@ function solve_phase( start::Vector{Float64}, lower_x = -1 * ones(3n) upper_x = ones(3n) num_constraints = 6 + filename = "ipopt_"*string(rand(UInt)) ipopt_options = Dict("constr_viol_tol" => tol, "acceptable_constr_viol_tol" => 1e-4, "max_iter" => 10_000, "max_cpu_time" => 30., - "print_level" => 0) - # if verbose ipopt_options["print_level"] = 5 end + "print_level" => 0, + "output_file" => filename) options = Options(solver=IPOPT(ipopt_options), derivatives=ForwardAD()) x, _, info = minimize(f!, vec(x0), num_constraints, lower_x, upper_x, final[1:6], final[1:6], options) + rm(filename) if info in [:Solve_Succeeded, :Solved_To_Acceptable_Level] return Result(true, info, reshape(x,(n,3))) else diff --git a/julia/src/mission.jl b/julia/src/mission.jl index 0284a6a..250154b 100644 --- a/julia/src/mission.jl +++ b/julia/src/mission.jl @@ -1,17 +1,30 @@ export Phase, Mission_Guess, Mission, Bad_Mission export test_phase1, test_phase2 -export test_mission_guess, test_mission_guess_simple -export test_mission, test_mission_simple +export test_mg, test_mission struct Phase planet::Body v∞_in::Vector{Float64} - δ::Float64 + v∞_out::Vector{Float64} tof::Float64 thrust_profile::Matrix{Float64} end -const test_phase1 = Phase(Venus, [10.4321, -6.3015, -0.01978], 0.2, 1.30464e7, zeros(20,3)) -const test_phase2 = Phase(Jupiter, [0.3, 7.1, 0.2], 2π, 3.9year, zeros(20,3)) + +const test_phase1 = Phase(Venus, [9., 10., 0.], [10., 9, 0.], 1.30464e7, zeros(20,3)) +const test_phase2 = Phase(Jupiter, [0.3, 7.1, 0.2], [0.3, 7.1, 0.2], 3.9year, zeros(20,3)) +const test_phases = [test_phase1, test_phase2] + +""" +A convenience function for calculating mass usage given a certain thrust profile +""" +function mass_consumption(sc::Sc, phase::Phase) + weighted_thrusting_time = 0.0 + n = size(phase.thrust_profile)[1] + for i in 1:size(phase.thrust_profile,1) + weighted_thrusting_time += norm(phase.thrust_profile[i,:]) * phase.tof/n + end + return weighted_thrusting_time*sc.mass_flow_rate +end struct Mission_Guess sc::Sc @@ -21,17 +34,28 @@ struct Mission_Guess phases::Vector{Phase} converged::Bool end -Mission_Guess(args...) = Mission_Guess(args..., false) -const test_mission_guess = Mission_Guess( bepi, - 12_000., - DateTime(1992, 11, 19), - [-3.4, 1.2, 0.1], - [test_phase1, test_phase2] ) -const test_mission_guess_simple = Mission_Guess(bepi, - 12_000., - DateTime(1992, 11, 19), - [-3.4, 1.2, 0.1], - [test_phase1]) + +""" +Constructor for a mission guess. Generally mission guesses are not converged +""" +function Mission_Guess(sc::Sc, mass::Float64, date::DateTime, v∞::Vector{Float64}, phases::Vector{Phase}) + # First do some checks to make sure that it's valid + mass_used = 0 + for phase in phases + mass_used += mass_consumption(sc, phase) + mass - mass_used > sc.dry_mass || throw(Mass_Error(mass - mass_used)) + v∞_in, v∞_out = phase.v∞_in, phase.v∞_out + if phase != phases[end] + norm(v∞_in) ≈ norm(v∞_out) || throw(V∞_Error(v∞_in, v∞_out)) + δ = acos( ( v∞_in ⋅ v∞_out ) / ( norm(v∞_in) * norm(v∞_out) ) ) + periapsis = (phase.planet.μ/(v∞_in ⋅ v∞_in)) * ( 1/sin(δ/2) - 1 ) + periapsis > 1.1phase.planet.r || throw(HitPlanet_Error()) + end + end + Mission_Guess(sc, mass, date, v∞, phases, false) +end + +const test_mg = Mission_Guess(bepi, 12_000., DateTime(1992,11,19), [-3.4,1.2,0.1], test_phases) struct Mission sc::Sc @@ -41,20 +65,33 @@ struct Mission phases::Vector{Phase} converged::Bool end -Mission(args...) = Mission(args..., true) -const test_mission = Mission(bepi, - 12_000., - DateTime(1992, 11, 19), - [-3.4, 1.2, 0.1], - [test_phase1, test_phase2]) -const test_mission_simple = Mission(bepi, - 12_000., - DateTime(1992, 11, 19), - [4.2984, -4.3272668, 1.43752], - [test_phase1]) + +""" +Constructor for a mission. Generally mission guesses are converged +""" +function Mission(sc::Sc, mass::Float64, date::DateTime, v∞::Vector{Float64}, phases::Vector{Phase}) + # First do some checks to make sure that it's valid + mass_used = 0 + for phase in phases + mass_used += mass_consumption(sc, phase) + mass - mass_used > sc.dry_mass || throw(Mass_Error(mass - mass_used)) + v∞_in, v∞_out = phase.v∞_in, phase.v∞_out + if phase != phases[end] + norm(v∞_in) ≈ norm(v∞_out) || throw(V∞_Error(v∞_in, v∞_out)) + δ = acos( ( v∞_in ⋅ v∞_out ) / ( norm(v∞_in) * norm(v∞_out) ) ) + periapsis = (phase.planet.μ/(v∞_in ⋅ v∞_in)) * ( 1/sin(δ/2) - 1 ) + periapsis > 1.1phase.planet.r || throw(HitPlanet_Error()) + end + end + Mission(sc, mass, date, v∞, phases, true) +end + +const test_mission = Mission(bepi, 12_000., DateTime(1992,11,19), [-3.4,1.2,0.1], test_phases) struct Bad_Mission - message::Symbol + message::String converged::Bool end -Bad_Mission(s) = Bad_Mission(s,false) +Bad_Mission(s::String) = Bad_Mission(s,false) +Bad_Mission(s::Symbol) = Bad_Mission(String(s),false) + diff --git a/julia/test/inner_loop/monotonic_basin_hopping.jl b/julia/test/inner_loop/monotonic_basin_hopping.jl index 4d0641f..c1adbd5 100644 --- a/julia/test/inner_loop/monotonic_basin_hopping.jl +++ b/julia/test/inner_loop/monotonic_basin_hopping.jl @@ -5,42 +5,54 @@ println("Testing Monotonic Basin Hopper") # First we test the random mission guess generator + println("Testing guess generator function") flybys = [Earth, Venus, Jupiter] launch_window = ( DateTime(2021,12,25), DateTime(2025,12,25) ) max_C3 = 10. max_v∞ = 8. - latest_arrival = DateTime(2035,12,25) + latest_arrival = DateTime(2030,12,25) random_guess = Thesis.mission_guess(flybys, bepi, 12_000., launch_window, max_C3, max_v∞, latest_arrival) @test typeof(random_guess) == Mission_Guess # Then the perturb function + println("Testing perturb function") mission_guess = Thesis.perturb(test_mission) @test mission_guess.launch_date != test_mission.launch_date @test mission_guess.launch_v∞ != test_mission.launch_v∞ for i in 1:2 @test mission_guess.phases[i].v∞_in != test_mission.phases[i].v∞_in - @test mission_guess.phases[i].δ != test_mission.phases[i].δ + @test mission_guess.phases[i].v∞_out != test_mission.phases[i].v∞_out @test mission_guess.phases[i].tof != test_mission.phases[i].tof end @test !mission_guess.converged # # Then the inner loop builder function - mission = Thesis.inner_loop_solve(test_mission_guess) + println("Testing inner loop solver function") + mission = Thesis.inner_loop_solve(test_mg) @test !mission.converged # For the valid case we need to use a lambert's solver - # TODO: This is probably not acceptable for how close I have to be - leave = DateTime(1992,11,19) - arrive = DateTime(1993,4,1) - time_leave = utc2et(Dates.format(leave,"yyyy-mm-ddTHH:MM:SS")) - time_arrive = utc2et(Dates.format(arrive,"yyyy-mm-ddTHH:MM:SS")) - earth_state = [spkssb(Earth.id, time_leave, "ECLIPJ2000"); 0.0] - venus_state = [spkssb(Venus.id, time_arrive, "ECLIPJ2000"); 0.0] + leave, arrive = DateTime(1992,11,19), DateTime(1993,4,1) + earth_state = state(Earth, leave) + venus_state = state(Venus, arrive) v∞_out, v∞_in, tof = Thesis.lamberts(Earth, Venus, leave, arrive) - phase = Phase(Venus, 1.01v∞_in, 0.2, tof, 0.01*ones(20,3)) + phase = Phase(Venus, 1.01v∞_in, -1.01v∞_in, tof, 0.01*ones(20,3)) mission_guess = Mission_Guess(bepi, 12_000., leave, v∞_out, [phase]) mission = Thesis.inner_loop_solve(mission_guess) @test mission.converged if !mission.converged println(mission.message) end + # Now we're ready to test the whole thing! + println("Testing whole thing") + flybys = [ Earth, Venus ] + start_mass = 12_000. + launch_window = DateTime(1997, 10, 1), DateTime(1997, 10, 31) + max_C3 = 15. + max_v∞ = 20. + latest_arrival = DateTime(2001, 1, 1) + best, archive = mbh(flybys, bepi, start_mass, launch_window, max_C3, max_v∞, latest_arrival, + verbose=true) + @test typeof(best) == Mission + + end diff --git a/notes.md b/notes.md index 7740e7f..8774df3 100644 --- a/notes.md +++ b/notes.md @@ -33,4 +33,13 @@ Todo: 1. ~~Add interpolate function~~ 2. ~~Add NLP solver sun test~~ -3. Start on MBH +3. ~~Start on MBH~~ + +## Refactor Notes (2021-09-26) + +Todo: + +1. Add mission validate function +2. Add turning angle function +3. Finish mbh +4. Test it