Skip to main content

Planet not following orbit in Solar System Simulation

I'm trying to simulate a solar system, but the planets don't seem to follow the orbit, instead, the distance between a planet and the sun increases. I can't figure out what's wrong with this code. Here is an MRE.

import numpy as np
import matplotlib.pyplot as plt


GRAVITATIONAL_CONSTANT = 6.674e-11
EARTH_MASS = 5.972 * (10**24)
SUN_MASS = 332954.355179 * EARTH_MASS
MERCURY_MASS = 0.06 * EARTH_MASS
AU2M = 1.495979e11
AUD2MS = AU2M / 86400


class SolarSystemBody:
    def __init__(self, name, mass, position, velocity):
        self.name = name
        self.mass = mass
        self.position = np.array(position, dtype=float) * AU2M
        self.velocity = np.array(velocity, dtype=float) * AUD2MS
        self.acceleration = np.zeros(3, dtype=float)

    def update_position(self, dt):
        self.position += self.velocity * dt + 0.5 * self.acceleration * dt**2

    def update_velocity(self, dt):
        self.velocity += self.acceleration * dt

    def gravitational_force(self, sun):
        r = sun.position - self.position
        distance = np.linalg.norm(r)
        direction = r / distance
        force_magnitude = GRAVITATIONAL_CONSTANT * self.mass * sun.mass / distance**2
        return force_magnitude * direction

    def calculate_acceleration(self, sun):
        force = self.gravitational_force(sun)
        self.acceleration = force / self.mass


mercury_position = [0.1269730114807624, 0.281031132701101, 0.01131924496141172]
mercury_velocity = [-0.03126433724097832, 0.01267637703164289, 0.00390363008183905]

sun = SolarSystemBody("Sun", SUN_MASS, [0, 0, 0], [0, 0, 0])
mercury = SolarSystemBody("Mercury", MERCURY_MASS, mercury_position, mercury_velocity)

dt = 3600 * 24
total_time = 365 * dt

pos = np.zeros((365, 3), dtype=float)
i = 0
for t in np.arange(0, total_time, dt):
    print(np.linalg.norm(sun.position - mercury.position))
    pos[i, :] = mercury.position
    mercury.calculate_acceleration(sun)
    mercury.update_velocity(dt)
    mercury.update_position(dt)
    i += 1


fig = plt.figure()
ax = fig.add_subplot(111, projection="3d")
ax.plot(pos[:, 0], pos[:, 1], pos[:, 2])
plt.show()

Getting the initial data using astroquery's HorizonsClass

from astropy.time import Time
from astroquery.jplhorizons import Horizons
import numpy as np

sim_start_date = "2023-01-01"

data = []
planet_id = 199  # Mercury

obj = Horizons(id=planet_id, location="@sun", epochs=Time(sim_start_date).jd, id_type='id').vectors()
name = obj["targetname"].data[0].split('(')[0].strip()
r = [np.double(obj[xi]) for xi in ['x', 'y', 'z']]
v = [np.double(obj[vxi]) for vxi in ['vx', 'vy', 'vz']]

Update:

Apparently mercury was slingshotting using the sun thats why the sudden increase in it's distance to sun.The issue seems to be that dt was too big so acceleration wasn't getting updated quickly. Decreasing dt fixes the issue.

dt = 360
total_time = 365 * 240 * dt


source https://stackoverflow.com/questions/76133729/planet-not-following-orbit-in-solar-system-simulation

Comments

Popular posts from this blog

ValueError: X has 10 features, but LinearRegression is expecting 1 features as input

So, I am trying to predict the model but its throwing error like it has 10 features but it expacts only 1. So I am confused can anyone help me with it? more importantly its not working for me when my friend runs it. It works perfectly fine dose anyone know the reason about it? cv = KFold(n_splits = 10) all_loss = [] for i in range(9): # 1st for loop over polynomial orders poly_order = i X_train = make_polynomial(x, poly_order) loss_at_order = [] # initiate a set to collect loss for CV for train_index, test_index in cv.split(X_train): print('TRAIN:', train_index, 'TEST:', test_index) X_train_cv, X_test_cv = X_train[train_index], X_test[test_index] t_train_cv, t_test_cv = t[train_index], t[test_index] reg.fit(X_train_cv, t_train_cv) loss_at_order.append(np.mean((t_test_cv - reg.predict(X_test_cv))**2)) # collect loss at fold all_loss.append(np.mean(loss_at_order)) # collect loss at order plt.plot(np.log(al...

Sorting large arrays of big numeric stings

I was solving bigSorting() problem from hackerrank: Consider an array of numeric strings where each string is a positive number with anywhere from to digits. Sort the array's elements in non-decreasing, or ascending order of their integer values and return the sorted array. I know it works as follows: def bigSorting(unsorted): return sorted(unsorted, key=int) But I didnt guess this approach earlier. Initially I tried below: def bigSorting(unsorted): int_unsorted = [int(i) for i in unsorted] int_sorted = sorted(int_unsorted) return [str(i) for i in int_sorted] However, for some of the test cases, it was showing time limit exceeded. Why is it so? PS: I dont know exactly what those test cases were as hacker rank does not reveal all test cases. source https://stackoverflow.com/questions/73007397/sorting-large-arrays-of-big-numeric-stings

How to load Javascript with imported modules?

I am trying to import modules from tensorflowjs, and below is my code. test.html <!DOCTYPE html> <html lang="en"> <head> <meta charset="UTF-8"> <title>Document</title </head> <body> <script src="https://cdn.jsdelivr.net/npm/@tensorflow/tfjs@2.0.0/dist/tf.min.js"></script> <script type="module" src="./test.js"></script> </body> </html> test.js import * as tf from "./node_modules/@tensorflow/tfjs"; import {loadGraphModel} from "./node_modules/@tensorflow/tfjs-converter"; const MODEL_URL = './model.json'; const model = await loadGraphModel(MODEL_URL); const cat = document.getElementById('cat'); model.execute(tf.browser.fromPixels(cat)); Besides, I run the server using python -m http.server in my command prompt(Windows 10), and this is the error prompt in the console log of my browser: Failed to loa...