Skip to content

Commit a0b642c

Browse files
BigAlsCodegithub-actionspre-commit-ci[bot]cclauss
authored
Physics/basic orbital capture (#8857)
* Added file basic_orbital_capture * updating DIRECTORY.md * added second source * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * fixed spelling errors * accepted changes * updating DIRECTORY.md * corrected spelling error * Added file basic_orbital_capture * added second source * fixed spelling errors * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * applied changes * reviewed and checked file * added doctest * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * removed redundant constnant * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * added scipy imports * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * added doctests to capture_radii and scipy const * fixed conflicts * finalizing file. Added tests * Update physics/basic_orbital_capture.py --------- Co-authored-by: github-actions <${GITHUB_ACTOR}@users.noreply.github.com> Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com> Co-authored-by: Christian Clauss <cclauss@me.com>
1 parent e406801 commit a0b642c

File tree

2 files changed

+179
-0
lines changed

2 files changed

+179
-0
lines changed

Diff for: ‎DIRECTORY.md

+1
Original file line numberDiff line numberDiff line change
@@ -741,6 +741,7 @@
741741

742742
## Physics
743743
* [Archimedes Principle](physics/archimedes_principle.py)
744+
* [Basic Orbital Capture](physics/basic_orbital_capture.py)
744745
* [Casimir Effect](physics/casimir_effect.py)
745746
* [Centripetal Force](physics/centripetal_force.py)
746747
* [Grahams Law](physics/grahams_law.py)

Diff for: ‎physics/basic_orbital_capture.py

+178
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,178 @@
1+
from math import pow, sqrt
2+
3+
from scipy.constants import G, c, pi
4+
5+
"""
6+
These two functions will return the radii of impact for a target object
7+
of mass M and radius R as well as it's effective cross sectional area σ(sigma).
8+
That is to say any projectile with velocity v passing within σ, will impact the
9+
target object with mass M. The derivation of which is given at the bottom
10+
of this file.
11+
12+
The derivation shows that a projectile does not need to aim directly at the target
13+
body in order to hit it, as R_capture>R_target. Astronomers refer to the effective
14+
cross section for capture as σ=π*R_capture**2.
15+
16+
This algorithm does not account for an N-body problem.
17+
18+
"""
19+
20+
21+
def capture_radii(
22+
target_body_radius: float, target_body_mass: float, projectile_velocity: float
23+
) -> float:
24+
"""
25+
Input Params:
26+
-------------
27+
target_body_radius: Radius of the central body SI units: meters | m
28+
target_body_mass: Mass of the central body SI units: kilograms | kg
29+
projectile_velocity: Velocity of object moving toward central body
30+
SI units: meters/second | m/s
31+
Returns:
32+
--------
33+
>>> capture_radii(6.957e8, 1.99e30, 25000.0)
34+
17209590691.0
35+
>>> capture_radii(-6.957e8, 1.99e30, 25000.0)
36+
Traceback (most recent call last):
37+
...
38+
ValueError: Radius cannot be less than 0
39+
>>> capture_radii(6.957e8, -1.99e30, 25000.0)
40+
Traceback (most recent call last):
41+
...
42+
ValueError: Mass cannot be less than 0
43+
>>> capture_radii(6.957e8, 1.99e30, c+1)
44+
Traceback (most recent call last):
45+
...
46+
ValueError: Cannot go beyond speed of light
47+
48+
Returned SI units:
49+
------------------
50+
meters | m
51+
"""
52+
53+
if target_body_mass < 0:
54+
raise ValueError("Mass cannot be less than 0")
55+
if target_body_radius < 0:
56+
raise ValueError("Radius cannot be less than 0")
57+
if projectile_velocity > c:
58+
raise ValueError("Cannot go beyond speed of light")
59+
60+
escape_velocity_squared = (2 * G * target_body_mass) / target_body_radius
61+
capture_radius = target_body_radius * sqrt(
62+
1 + escape_velocity_squared / pow(projectile_velocity, 2)
63+
)
64+
return round(capture_radius, 0)
65+
66+
67+
def capture_area(capture_radius: float) -> float:
68+
"""
69+
Input Param:
70+
------------
71+
capture_radius: The radius of orbital capture and impact for a central body of
72+
mass M and a projectile moving towards it with velocity v
73+
SI units: meters | m
74+
Returns:
75+
--------
76+
>>> capture_area(17209590691)
77+
9.304455331329126e+20
78+
>>> capture_area(-1)
79+
Traceback (most recent call last):
80+
...
81+
ValueError: Cannot have a capture radius less than 0
82+
83+
Returned SI units:
84+
------------------
85+
meters*meters | m**2
86+
"""
87+
88+
if capture_radius < 0:
89+
raise ValueError("Cannot have a capture radius less than 0")
90+
sigma = pi * pow(capture_radius, 2)
91+
return round(sigma, 0)
92+
93+
94+
if __name__ == "__main__":
95+
from doctest import testmod
96+
97+
testmod()
98+
99+
"""
100+
Derivation:
101+
102+
Let: Mt=target mass, Rt=target radius, v=projectile_velocity,
103+
r_0=radius of projectile at instant 0 to CM of target
104+
v_p=v at closest approach,
105+
r_p=radius from projectile to target CM at closest approach,
106+
R_capture= radius of impact for projectile with velocity v
107+
108+
(1)At time=0 the projectile's energy falling from infinity| E=K+U=0.5*m*(v**2)+0
109+
110+
E_initial=0.5*m*(v**2)
111+
112+
(2)at time=0 the angular momentum of the projectile relative to CM target|
113+
L_initial=m*r_0*v*sin(Θ)->m*r_0*v*(R_capture/r_0)->m*v*R_capture
114+
115+
L_i=m*v*R_capture
116+
117+
(3)The energy of the projectile at closest approach will be its kinetic energy
118+
at closest approach plus gravitational potential energy(-(GMm)/R)|
119+
E_p=K_p+U_p->E_p=0.5*m*(v_p**2)-(G*Mt*m)/r_p
120+
121+
E_p=0.0.5*m*(v_p**2)-(G*Mt*m)/r_p
122+
123+
(4)The angular momentum of the projectile relative to the target at closest
124+
approach will be L_p=m*r_p*v_p*sin(Θ), however relative to the target Θ=90°
125+
sin(90°)=1|
126+
127+
L_p=m*r_p*v_p
128+
(5)Using conservation of angular momentum and energy, we can write a quadratic
129+
equation that solves for r_p|
130+
131+
(a)
132+
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
133+
134+
(b)
135+
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
136+
137+
(c) b plugs int a|
138+
v**2=((v*R_capture)/r_p)**2-(2*G*Mt)/r_p->
139+
140+
v**2-(v**2)*(R_c**2)/(r_p**2)+(2*G*Mt)/r_p=0->
141+
142+
(v**2)*(r_p**2)+2*G*Mt*r_p-(v**2)*(R_c**2)=0
143+
144+
(d) Using the quadratic formula, we'll solve for r_p then rearrange to solve to
145+
R_capture
146+
147+
r_p=(-2*G*Mt ± sqrt(4*G^2*Mt^2+ 4(v^4*R_c^2)))/(2*v^2)->
148+
149+
r_p=(-G*Mt ± sqrt(G^2*Mt+v^4*R_c^2))/v^2->
150+
151+
r_p<0 is something we can ignore, as it has no physical meaning for our purposes.->
152+
153+
r_p=(-G*Mt)/v^2 + sqrt(G^2*Mt^2/v^4 + R_c^2)
154+
155+
(e)We are trying to solve for R_c. We are looking for impact, so we want r_p=Rt
156+
157+
Rt + G*Mt/v^2 = sqrt(G^2*Mt^2/v^4 + R_c^2)->
158+
159+
(Rt + G*Mt/v^2)^2 = G^2*Mt^2/v^4 + R_c^2->
160+
161+
Rt^2 + 2*G*Mt*Rt/v^2 + G^2*Mt^2/v^4 = G^2*Mt^2/v^4 + R_c^2->
162+
163+
Rt**2 + 2*G*Mt*Rt/v**2 = R_c**2->
164+
165+
Rt**2 * (1 + 2*G*Mt/Rt *1/v**2) = R_c**2->
166+
167+
escape velocity = sqrt(2GM/R)= v_escape**2=2GM/R->
168+
169+
Rt**2 * (1 + v_esc**2/v**2) = R_c**2->
170+
171+
(6)
172+
R_capture = Rt * sqrt(1 + v_esc**2/v**2)
173+
174+
Source: Problem Set 3 #8 c.Fall_2017|Honors Astronomy|Professor Rachel Bezanson
175+
176+
Source #2: http://www.nssc.ac.cn/wxzygx/weixin/201607/P020160718380095698873.pdf
177+
8.8 Planetary Rendezvous: Pg.368
178+
"""

0 commit comments

Comments
 (0)