Seveda bi v gravitacijsko simulacijo radi postavili raketo, ki si lahko spreminja hitrost. Prvo ponovimo raketno enačbo, ki jo lahko zelo hitro izpeljemo iz drugega Newtonovega zakona. Tu moramo upoštevati spremembo gibalne količine zaradi spremembe mase rakete in zaradi odmetavanja materiala s konstantno hitrostjo (to bomo označili z ) - ta sprememba je negativna saj odmetavamo v nasprotno smer gibanja. Silo izračunamo s primerjavo tlakov in to pomnožimo s površino motorja:
Tu smo uvedli efektivno hitrost motorja zaradi lažje oznake in v zadnji enakosti upoštevali, da se masa zaradi odmetavanja zmanjšuje. Dobili smo zelo preprosto diferencialno enačbo, katere rešitev dobimo s preprosto integracijo:
Spremembo v hitrosti torej lahko podamo v obliki spremembe hitrosti ali mase. Iz enačbe uvidimo, da večje razmerje začetne in končne mase pomeni večjo spremembo hitrost - dobro se je torej sproti znebiti odvečne mase (rakete imajo zato več stopenj). V simulatorju smo spisali razred, ki upošteva vse od naštetega, in nam dovoli premikanje po prostoru. Ker simulacija ni čisto natančna in za čudne orbite ne poznamo točne periode, raketo usmerjamo s tako imenovanimi manevrskimi vozlišči. Te postavimo v prostor in jim določimo radij in spremembo, ter jih povežemo s telesom. Ko se ta približa znotraj danega radija naredi spremembo in označi, da je manevrsko vozlišče bilo uporabljeno.
from simulatorv2 import *
import numpy as np
sun = Presets.sun(Vec2)
rocket = Rocket(
name="Apollo",
mass=1000 / Units.SUN_MASS,
dry_mass=500 / Units.SUN_MASS,
effective_velocity=4500 * Units.ASTRONOMICAL_UNIT / Units.YEAR,
radius=0.01,
color=(255, 0, 0),
position=Vec2(1, 0),
velocity=Vec2(0, 2*np.pi),
nodes=[
ManeuverNode(
position=Vec2(0, 1),
radius=0.00001,
eccentricity=0.8
)
]
)
sim = Simulator(
bodies=[sun, rocket],
steps=2000,
timestep=0.001,
diagnostics=[
RocketMass(),
TotalEnergy()
]
)
result = sim.run()
viewer = SimulationViewer(
result,
ViewerConfig(
camera_mode="follow",
follow_body="Apollo",
fixed_center=(0, 0),
trail_length=10,
view_size=0.5,
diagnostics=[
DiagnosticView(
"Celotna Energija",
["vsi"]
),
],
references=[
ReferenceConic(0.0, opacity=0.5),
ReferenceConic(0.8, opacity=0.5, rotation=np.pi*0.5)
]
),
)
viewer.show()