diff --git a/src/nature.rs b/src/nature.rs
index 78b9284..360059e 100644
--- a/src/nature.rs
+++ b/src/nature.rs
@@ -28,9 +28,11 @@ pub const G: f64 = 6.6743015e-11; // Gravitational constant in Nm²/kg²
 pub const SOL_RADIUS: f64 = 696_300_000.0;
 pub const JUPITER_RADIUS: f64 = 71_492_000.0;
 pub const JUPITER_RING_RADIUS: f64 = 229_000_000.0;
+pub const EARTH_RADIUS: f64 = 6_371_000.0;
 
 pub const SOL_MASS: f64 = 1.9885e30;
 pub const JUPITER_MASS: f64 = 1.8982e27;
+pub const EARTH_MASS: f64 = 5.972168e24;
 
 // Each star's values: (x, y, z, magnitude, color index, distance, name)
 pub const STARS: &[(f32, f32, f32, f32, f32, f32, &str)] = &include!("data/stars.in");
@@ -180,6 +182,22 @@ pub fn orbital_velocity(coords: Vec3, mass: f64) -> Vec3 {
     return perpendicular * speed;
 }
 
+/// Calculates the acceleration towards a mass in m/s
+pub fn gravitational_acceleration(coords: Vec3, mass: f64) -> Vec3 {
+    let r_squared = coords.length_squared() as f64;
+    let acceleration_magnitude = G * mass / r_squared;
+    return -acceleration_magnitude as f32 * (coords / r_squared.sqrt() as f32);
+}
+
+#[test]
+fn test_gravitational_acceleration() {
+    let coords = Vec3::new(EARTH_RADIUS as f32, 0.0, 0.0);
+    let mass = EARTH_MASS;
+    let g = gravitational_acceleration(coords, mass);
+    let g_rounded = (g * 10.0).round() / 10.0;
+    assert_eq!(g_rounded, Vec3::new(-9.8, 0.0, 0.0));
+}
+
 pub fn phase_dist_to_coords(phase_radians: f64, distance: f64) -> DVec3 {
     return DVec3::new(
         distance * phase_radians.cos(),