Commit d4fa1921 authored by Olaf's avatar Olaf

Start saving every n timesteps instead of every timestep

parent 59de56d0
This diff is collapsed.
......@@ -84,7 +84,7 @@ def calc(bodies, settings):
pass
else:
pass
positions, velocities, times, E_kin = dynamics(bodies, settings.dimensions, settings.timestep, settings.time, settings.barneshut, settings.theta, settings.approx)
positions, velocities, times, E_kin = dynamics(bodies, settings.dimensions, settings.timestep, settings.time, settings.barneshut, settings.theta, settings.approx, settings.timesave)
for i, body in enumerate(bodies):
if isinstance(body, Body):
body.add_trajectory((positions.transpose(1,0,2)[i]).tolist())
......
......@@ -2,6 +2,7 @@ import numpy as np
import warnings
from predef import *
from barnes_hut import *
import math
def extract_data(bodies, D):
......@@ -224,7 +225,7 @@ def force_approx(mass, positions, D, n):
return Force
def dynamics(bodies, D, h, t_max, barnes_hut, theta, approximation):
def dynamics(bodies, D, h, t_max, barnes_hut, theta, approximation, j):
"""
Calculates the dynamics (trajectories) of all bodies
......@@ -242,6 +243,10 @@ def dynamics(bodies, D, h, t_max, barnes_hut, theta, approximation):
states if the Barnes-Hut algorithm should be used or not
theta: float
ratio used in the Barnes-Hut algorithm
approximation: bool
states if only interaction on asteroids and not of the asteroids should be taken
j: int
every j_th timestep the data is saved
Returns:
--------
......@@ -264,35 +269,45 @@ def dynamics(bodies, D, h, t_max, barnes_hut, theta, approximation):
else:
force0 = force_calc(mass, positions, D)
all_pos = np.zeros((T+1, positions.shape[0], positions.shape[1]))
all_pos = np.zeros((math.ceil((T+1)/j), positions.shape[0], positions.shape[1]))
times = [0]
E_kin = []
all_pos[0] = positions
index = 0
if barnes_hut and theta != 0 and not approximation:
for i in range(1,T+1):
all_pos[i] = position(all_pos[i-1], velo, force0, h)
force1 = force_barneshut(all_pos[i], mass, theta)
positions = position(positions, velo, force0, h)
force1 = force_barneshut(positions, mass, theta)
velo = velocity(velo, force0, force1, h)
force0 = force1
E_kin.append(0.5*mass*(np.linalg.norm(velo, axis=1)**2))
times.append(h*i)
if i%j == 0:
index += 1
all_pos[index] = positions
E_kin.append(0.5*mass*(np.linalg.norm(velo, axis=1)**2))
times.append(h*i)
elif approximation and not barnes_hut:
for i in range(1,T+1):
all_pos[i] = position(all_pos[i-1], velo, force0, h)
force1 = force_approx(mass, all_pos[i], D, n)
positions = position(positions, velo, force0, h)
force1 = force_approx(mass, positions, D, n)
velo = velocity(velo, force0, force1, h)
force0 = force1
E_kin.append(0.5*mass*(np.linalg.norm(velo, axis=1)**2))
times.append(i*h)
if i%j == 0:
index += 1
all_pos[index] = positions
E_kin.append(0.5*mass*(np.linalg.norm(velo, axis=1)**2))
times.append(h*i)
else:
for i in range(1,T+1):
all_pos[i] = position(all_pos[i-1], velo, force0, h)
force1 = force_calc(mass, all_pos[i], D)
positions = position(positions, velo, force0, h)
force1 = force_calc(mass, positions, D)
velo = velocity(velo, force0, force1, h)
force0 = force1
E_kin.append(0.5*mass*(np.linalg.norm(velo, axis=1)**2))
times.append(h*i)
if i%j == 0:
index += 1
all_pos[index] = positions
E_kin.append(0.5*mass*(np.linalg.norm(velo, axis=1)**2))
times.append(h*i)
# return
return all_pos, velo, times, E_kin
......@@ -517,6 +517,11 @@ class Settings:
self.approx = sets['approx']
except:
self.approx = False
try:
self.timesave = sets['save_timestep']
except:
self.timesave = 1
......
......@@ -166,11 +166,12 @@ Resonance_belt_parameters = {
'N' : 10, # number of asteroids
'mass' : 0.012, # total mass of all asteroids combined
'R' : 3.1, # orbital radius from the Sun
'theta' : 0,
'colour' : 'lightgray',
'id' : 10,
'name' : 'Resonance Asteroid Belt'
}
Resonance_belt = Resonance_astroids(Resonance_belt_parameters); solar_system.append(Resonance_belt)
# Resonance_belt = Resonance_astroids(Resonance_belt_parameters); solar_system.append(Resonance_belt)
""" Quick writing """
All = 'all'
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment