Is there a range of launch velocities that would get my satellite as close to Mars as possible?

Question

Introduction

Hi, I am trying to send a satellite which has a mass of 1062kg from Earth to Mars. My aim is to get it as close as possible to Mars but not orbit it. I am looking for a range of possible launch velocities for an Hohmann transfer orbit to Mars. For the snippets of relevant parts of my code I provide, you can assume that that my simulation of the solar system is working perfectly. The numerical algorithm that I am using is the Beeman with a timestep of 31577s.

def main():
    

    # get the velocity range
    velocity_range = np.linspace(0, 5000, 20)  # Example range of additional velocities
    min_distances = []
    velocities = []
    # run the simulation for additional velocities
    for additional_velocity in velocity_range:

        velocity_at_launch = sim_beeman.run(additional_velocity)
      
        
        velocities.append(velocity_at_launch)
        min_distance = sim_beeman.calculate_min_distance()
        min_distances.append(min_distance)


    print("velocities: ", velocities)

     
    # plot the velocity range vs minimum distance
    plt.plot(velocities, min_distances)
    plt.xlabel("satellite velocity at launch(m/s)")
    plt.ylabel("Minimum Distance (km)")
    plt.title("Minimum Distance vs satellite velocity at launch")
    plt.show()

def hohmann_transfer_parameters(G, mass_sun, radius_earth, radius_mars):
   

    # Semi-major axis of the transfer orbit
    a_transfer = (radius_earth + radius_mars) / 2  
    V_earth = np.sqrt(G * mass_sun / radius_earth)  # Earth's orbital velocity
    V_mars = np.sqrt(G * mass_sun / radius_mars)  # Mars' orbital velocity

    V_transfer_perihelion = np.sqrt(G * mass_sun * (2 / radius_earth - 1 / a_transfer))
    V_transfer_aphelion = np.sqrt(G * mass_sun * (2 / radius_mars - 1 / a_transfer))

    Delta_V1 = V_transfer_perihelion - V_earth  # ΔV to enter the transfer orbit


    return Delta_V1

The above equations are what I found to calculate the first burn velocity. Ignore the calculations for aphelion. The Hohmann function is not part of a class and all of the functions below are class functions.

 # function to run the simulation
def run(self, additional_velocity):
    # initialize simulation data for each run
    self.initialize_simulation_data()
    velocity_at_launch = 0
    self.satellite_launched = False
   
    # iterate over the number of iterations
    for i in range(self.iteration_counts):
        
        # update the position and velocity of each body in each iteration
        self.update_position()
        self.update_velocity()
    
        # Calculate Earth-Sun distance in the current iteration
        earth_distance_to_sun = np.linalg.norm(self.positions['earth'][-1])
        # get mars distance to sun
        mars_distance_to_sun = np.linalg.norm(self.positions['mars'][-1])

        # get the distance between earth and mars
        earth_to_mars_distance = np.linalg.norm(self.positions['earth'][-1] - self.positions['mars'][-1])
        #print("Current distance to sun: ", current_distance_to_sun)


        # check if the angle between the sun-earth and sun-mars vectors is about 44.3 degrees for launch
        if self.is_launch_time(earth_distance_to_sun, mars_distance_to_sun, earth_to_mars_distance ) and not self.satellite_launched: 
            # launch the satellite
            velocity_at_launch = self.launch_satellite(1.989e30, 147.1e9, 227.9e9, additional_velocity)
            self.satellite_launched = True
        


       
    return velocity_at_launch

def launch_satellite(self, mass_sun, radius_earth, radius_mars, additional_velocity):            
        print("Attempting to launch satellite near perihelion...")
        # Calculate Hohmann transfer parameters
        Delta_V1 = hohmann_transfer_parameters(self.G, mass_sun, radius_earth, radius_mars)
        print("Delta V1: ", Delta_V1)
        
        # get the earth's velocity vector at launch timestep
        earth_velocity = self.velocities["earth"][-1]  
        # get the earth's position vector at launch timestep
        earth_position = self.positions["earth"][-1]  # Earth's position vector
       

        # Calculate the unit vector in the direction of Earth's velocity
        earth_velocity_unit_vector = earth_velocity / np.linalg.norm(earth_velocity)
        earth_position_unit_vector = earth_position / np.linalg.norm(earth_position)

        # Calculate the satellite's velocity by applying ΔV in the direction of Earth's velocity
        satellite_velocity = earth_velocity + (Delta_V1 + additional_velocity) * earth_velocity_unit_vector + 1000

        # Calculate the satellite's position; placing it 200 km above Earth's surface
    
        satellite_position = earth_position + np.array([0, 6371e3 + 200e3])*earth_position_unit_vector

        # print parameters for launch
        print("Satellite velocity vector at launch:", satellite_velocity)
        print("Satellite position vector at launch:", satellite_position)
        
        print("satellite velocity at launch: ", np.linalg.norm(satellite_velocity))
        print("satellite position at launch : ", np.linalg.norm(satellite_position))

        print("Earth velocity vector at launch: ", earth_velocity)
        print("Earth position vector at launch: ", earth_position)

        print("Earth velocity magnitude at launch: ", np.linalg.norm(earth_velocity))
        print("Earth position magnitude at launch: ", np.linalg.norm(earth_position))

        print("\n\n")

        # Update satellite data
        self.velocities["satellite"].append(satellite_velocity)
        self.positions["satellite"].append(satellite_position)
        self.accelerations["satellite"].append(np.zeros(2))


       
        return np.linalg.norm(satellite_velocity)

        
          
def is_launch_time(self, sun_to_earth_distance, sun_to_mars_distance, Earth_to_mars_distance):
   # use the cosine rule to calculate the angle between the sun-earth and sun-mars vectors 
    angle_rads = math.acos((sun_to_earth_distance**2 + sun_to_mars_distance**2 - Earth_to_mars_distance**2) / (2 * sun_to_earth_distance * sun_to_mars_distance))
    angle_deg = np.degrees(angle_rads)
   
    # if the angle is between 44 and 46 degrees return true
    if 44 <= angle_deg <= 46:
        return True
    return False

Here are what my data structures for storing positions look like:

 {'sun': [array([0., 0.])], 'mercury': [array([5.78952e+10, 0.00000e+00])], 'venus': [array([1.081608e+11, 0.000000e+00])], 'earth': [array([1.496e+11, 0.000e+00])], 'mars': [array([2.279904e+11, 0.000000e+00])], 'jupiter': [array([7.7792e+11, 0.0000e+00])], 'saturn': [array([1.4212e+12, 0.0000e+00])], 'uranus': [array([2.87232e+12, 0.00000e+00])], 'neptune': [array([4.50296e+12, 0.00000e+00])], 'satellite': [array([0., 0.])]}

Above is the initial state of self.positions, where everything is in SI units. Velocities are also stored exactly the same way.

Problem

I read on the internet that the best time to launch a satellite for a hohmann transfer orbit is when the angle between earth and mars is about 44.3 degrees with sun as the origin. This is what I used to decide the timing of my launch. Since my aim was to discover the range of velocities, I ran the simulation multiple times and printed out the launch parameters in SI units, here are some of them:

Attempting to launch satellite near perihelion...
Delta V1:  3078.5906999237886
Satellite velocity vector at launch: [-31805.41409535  -1371.22945039]
Satellite position vector at launch: [-1.07767306e+10  1.49062981e+11]
satellite velocity at launch:  31834.95933690441
satellite position at launch :  149452033556.3857
Earth velocity vector at launch:  [-29734.83427786  -2149.2828756 ]
Earth position vector at launch:  [-1.07767306e+10  1.49056427e+11]
Earth velocity magnitude at launch:  29812.409939671335
Earth position magnitude at launch:  149445496725.23538



Attempting to launch satellite near perihelion...
Delta V1:  3078.5906999237886
Satellite velocity vector at launch: [-32067.88721996  -1390.20144085]
Satellite position vector at launch: [-1.07767306e+10  1.49062981e+11]
satellite velocity at launch:  32098.006959903447
satellite position at launch :  149452033556.3857
Earth velocity vector at launch:  [-29734.83427786  -2149.2828756 ]
Earth position vector at launch:  [-1.07767306e+10  1.49056427e+11]
Earth velocity magnitude at launch:  29812.409939671335
Earth position magnitude at launch:  149445496725.23538



Attempting to launch satellite near perihelion...
Delta V1:  3078.5906999237886
Satellite velocity vector at launch: [-32330.36034456  -1409.1734313 ]
Satellite position vector at launch: [-1.07767306e+10  1.49062981e+11]
satellite velocity at launch:  32361.05637597122
satellite position at launch :  149452033556.3857
Earth velocity vector at launch:  [-29734.83427786  -2149.2828756 ]
Earth position vector at launch:  [-1.07767306e+10  1.49056427e+11]
Earth velocity magnitude at launch:  29812.409939671335
Earth position magnitude at launch:  149445496725.23538



Attempting to launch satellite near perihelion...
Delta V1:  3078.5906999237886
Satellite velocity vector at launch: [-32592.83346917  -1428.14542175]
Satellite position vector at launch: [-1.07767306e+10  1.49062981e+11]
satellite velocity at launch:  32624.10754173499
satellite position at launch :  149452033556.3857
Earth velocity vector at launch:  [-29734.83427786  -2149.2828756 ]
Earth position vector at launch:  [-1.07767306e+10  1.49056427e+11]
Earth velocity magnitude at launch:  29812.409939671335
Earth position magnitude at launch:  149445496725.23538



Attempting to launch satellite near perihelion...
Delta V1:  3078.5906999237886
Satellite velocity vector at launch: [-32855.30659378  -1447.1174122 ]
Satellite position vector at launch: [-1.07767306e+10  1.49062981e+11]
satellite velocity at launch:  32887.160415209415
satellite position at launch :  149452033556.3857
Earth velocity vector at launch:  [-29734.83427786  -2149.2828756 ]
Earth position vector at launch:  [-1.07767306e+10  1.49056427e+11]
Earth velocity magnitude at launch:  29812.409939671335
Earth position magnitude at launch:  149445496725.23538



Attempting to launch satellite near perihelion...
Delta V1:  3078.5906999237886
Satellite velocity vector at launch: [-33117.77971838  -1466.08940266]
Satellite position vector at launch: [-1.07767306e+10  1.49062981e+11]
satellite velocity at launch:  33150.21495574151
satellite position at launch :  149452033556.3857
Earth velocity vector at launch:  [-29734.83427786  -2149.2828756 ]
Earth position vector at launch:  [-1.07767306e+10  1.49056427e+11]
Earth velocity magnitude at launch:  29812.409939671335
Earth position magnitude at launch:  149445496725.23538

Finally, I made a plot of Min distance from mars vs launch velocities that looks like:

The lowest value y value of the graph is about 10^6 km at about 37580m/s. However, the radius of Mars is nearly 3785 km this suggests that I could get even closer to Mars, about 600 km above the surface of Mars. I don't know what I am missing about calculating the correct parameters or time of launch. Thank you for your time and help.