The Algorithms logo
The Algorithms

Basic Orbital Capture

from math import pow, sqrt

from scipy.constants import G, c, pi

These two functions will return the radii of impact for a target object
of mass M and radius R as well as it's effective cross sectional area sigma.
That is to say any projectile with velocity v passing within sigma, will impact the
target object with mass M. The derivation of which is given at the bottom
of this file.

The derivation shows that a projectile does not need to aim directly at the target
body in order to hit it, as  R_capture>R_target. Astronomers refer to the effective
cross section for capture as sigma=π*R_capture**2.

This algorithm does not account for an N-body problem.


def capture_radii(
    target_body_radius: float, target_body_mass: float, projectile_velocity: float
) -> float:
    Input Params:
    target_body_radius: Radius of the central body SI units: meters | m
    target_body_mass: Mass of the central body SI units: kilograms | kg
    projectile_velocity: Velocity of object moving toward central body
        SI units: meters/second | m/s
    >>> capture_radii(6.957e8, 1.99e30, 25000.0)
    >>> capture_radii(-6.957e8, 1.99e30, 25000.0)
    Traceback (most recent call last):
    ValueError: Radius cannot be less than 0
    >>> capture_radii(6.957e8, -1.99e30, 25000.0)
    Traceback (most recent call last):
    ValueError: Mass cannot be less than 0
    >>> capture_radii(6.957e8, 1.99e30, c+1)
    Traceback (most recent call last):
    ValueError: Cannot go beyond speed of light

    Returned SI units:
    meters | m

    if target_body_mass < 0:
        raise ValueError("Mass cannot be less than 0")
    if target_body_radius < 0:
        raise ValueError("Radius cannot be less than 0")
    if projectile_velocity > c:
        raise ValueError("Cannot go beyond speed of light")

    escape_velocity_squared = (2 * G * target_body_mass) / target_body_radius
    capture_radius = target_body_radius * sqrt(
        1 + escape_velocity_squared / pow(projectile_velocity, 2)
    return round(capture_radius, 0)

def capture_area(capture_radius: float) -> float:
    Input Param:
    capture_radius: The radius of orbital capture and impact for a central body of
    mass M and a projectile moving towards it with velocity v
        SI units: meters | m
    >>> capture_area(17209590691)
    >>> capture_area(-1)
    Traceback (most recent call last):
    ValueError: Cannot have a capture radius less than 0

    Returned SI units:
    meters*meters | m**2

    if capture_radius < 0:
        raise ValueError("Cannot have a capture radius less than 0")
    sigma = pi * pow(capture_radius, 2)
    return round(sigma, 0)

if __name__ == "__main__":
    from doctest import testmod



Let: Mt=target mass, Rt=target radius, v=projectile_velocity,
     r_0=radius of projectile at instant 0 to CM of target
     v_p=v at closest approach,
     r_p=radius from projectile to target CM at closest approach,
     R_capture= radius of impact for projectile with velocity v

(1)At time=0  the projectile's energy falling from infinity| E=K+U=0.5*m*(v**2)+0


(2)at time=0 the angular momentum of the projectile relative to CM target|


(3)The energy of the projectile at closest approach will be its kinetic energy
   at closest approach plus gravitational potential energy(-(GMm)/R)|


(4)The angular momentum of the projectile relative to the target at closest
   approach will be L_p=m*r_p*v_p*sin(Θ), however relative to the target Θ=90°

(5)Using conservation of angular momentum and energy, we can write a quadratic
   equation that solves for r_p|

    Ei=Ep-> 0.5*m*(v**2)=0.5*m*(v_p**2)-(G*Mt*m)/r_p-> v**2=v_p**2-(2*G*Mt)/r_p

    Li=Lp-> m*v*R_capture=m*r_p*v_p-> v*R_capture=r_p*v_p-> v_p=(v*R_capture)/r_p

   (c) b plugs int a|



   (d) Using the quadratic formula, we'll solve for r_p then rearrange to solve to

    r_p=(-2*G*Mt ± sqrt(4*G^2*Mt^2+ 4(v^4*R_c^2)))/(2*v^2)->

    r_p=(-G*Mt ± sqrt(G^2*Mt+v^4*R_c^2))/v^2->

    r_p<0 is something we can ignore, as it has no physical meaning for our purposes.->

    r_p=(-G*Mt)/v^2 + sqrt(G^2*Mt^2/v^4 + R_c^2)

   (e)We are trying to solve for R_c. We are looking for impact, so we want r_p=Rt

    Rt + G*Mt/v^2 = sqrt(G^2*Mt^2/v^4 + R_c^2)->

    (Rt + G*Mt/v^2)^2 = G^2*Mt^2/v^4 + R_c^2->

    Rt^2 + 2*G*Mt*Rt/v^2 + G^2*Mt^2/v^4 = G^2*Mt^2/v^4 + R_c^2->

    Rt**2 + 2*G*Mt*Rt/v**2 = R_c**2->

    Rt**2 * (1 + 2*G*Mt/Rt *1/v**2) = R_c**2->

    escape velocity = sqrt(2GM/R)= v_escape**2=2GM/R->

    Rt**2 * (1 + v_esc**2/v**2) = R_c**2->

    R_capture = Rt * sqrt(1 + v_esc**2/v**2)

Source: Problem Set 3 #8 c.Fall_2017|Honors Astronomy|Professor Rachel Bezanson

Source #2:
           8.8 Planetary Rendezvous: Pg.368