Merge branch 'master' into use_momentum_changes
This commit is contained in:
commit
a9450179be
|
@ -10,7 +10,7 @@ c = np.array(
|
||||||
, 1 # *c*exp(
|
, 1 # *c*exp(
|
||||||
, -0.7 # c
|
, -0.7 # c
|
||||||
, 0 # (r - c))
|
, 0 # (r - c))
|
||||||
, 0 # + c*exp(
|
, +5 # + c*exp(
|
||||||
, -.1 # c
|
, -.1 # c
|
||||||
, 2 # (r - c))
|
, 2 # (r - c))
|
||||||
, -0 # + c*exp(
|
, -0 # + c*exp(
|
||||||
|
|
3
force.py
3
force.py
|
@ -4,9 +4,12 @@ import matplotlib.pyplot as plt
|
||||||
|
|
||||||
from coefficients import c
|
from coefficients import c
|
||||||
|
|
||||||
|
# This is the quite unreadable way to create
|
||||||
|
# UFuncs with given parameters. FIXME: add this to another module.
|
||||||
force_function = UFuncWrapper(0, c)
|
force_function = UFuncWrapper(0, c)
|
||||||
potential_function = UFuncWrapper(2, c)
|
potential_function = UFuncWrapper(2, c)
|
||||||
|
|
||||||
|
# Plot the force and potential.
|
||||||
r = np.arange(0, 100, 0.02, dtype=np.float16)
|
r = np.arange(0, 100, 0.02, dtype=np.float16)
|
||||||
f, = plt.plot(r, force_function(r), label="force")
|
f, = plt.plot(r, force_function(r), label="force")
|
||||||
p, = plt.plot(r, potential_function(r), label="potential")
|
p, = plt.plot(r, potential_function(r), label="potential")
|
||||||
|
|
19
particles.py
19
particles.py
|
@ -11,23 +11,34 @@ from coefficients import c
|
||||||
#force_function = UFuncWrapper(0, c)
|
#force_function = UFuncWrapper(0, c)
|
||||||
#interaction2D = UFuncWrapper(1, c)
|
#interaction2D = UFuncWrapper(1, c)
|
||||||
|
|
||||||
|
# Borders for both the plot and the boundary condition
|
||||||
|
# (the boundary condition might be deactivated, when creating
|
||||||
|
# the BrownIterator).
|
||||||
borders_x = [-100, 100]
|
borders_x = [-100, 100]
|
||||||
borders_y = [-100, 100]
|
borders_y = [-100, 100]
|
||||||
n_particles = 600
|
n_particles = 600
|
||||||
|
# Idk, seems to not do anyting.
|
||||||
frames = 100
|
frames = 100
|
||||||
|
# Only spawn in 1/x of the borders.
|
||||||
spawn_restriction = 3
|
spawn_restriction = 3
|
||||||
|
# Time resolution. Note that setting this to a too
|
||||||
|
# high value (i.e. low resolution) will lead to
|
||||||
|
# erratic behaviour, because potentials can be skipped.
|
||||||
dt = 0.1
|
dt = 0.1
|
||||||
c[-1] = dt
|
c[-1] = dt
|
||||||
|
|
||||||
|
# Initial positions.
|
||||||
x_coords = np.random.uniform(borders_x[0] / spawn_restriction, borders_x[1] / spawn_restriction, n_particles).astype(np.float16)
|
x_coords = np.random.uniform(borders_x[0] / spawn_restriction, borders_x[1] / spawn_restriction, n_particles).astype(np.float16)
|
||||||
y_coords = np.random.uniform(borders_y[0] / spawn_restriction, borders_y[1] / spawn_restriction, n_particles).astype(np.float16)
|
y_coords = np.random.uniform(borders_y[0] / spawn_restriction, borders_y[1] / spawn_restriction, n_particles).astype(np.float16)
|
||||||
|
|
||||||
|
|
||||||
|
# Initial momenta are 0.
|
||||||
x_momenta = np.zeros(n_particles, dtype=np.float16)
|
x_momenta = np.zeros(n_particles, dtype=np.float16)
|
||||||
y_momenta = np.zeros(n_particles, dtype=np.float16)
|
y_momenta = np.zeros(n_particles, dtype=np.float16)
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
# Prepare the plot, remove axis & stuff.
|
||||||
fig = plt.figure(figsize=(7, 7))
|
fig = plt.figure(figsize=(7, 7))
|
||||||
ax = fig.add_axes([0, 0, 1, 1], frameon=False)
|
ax = fig.add_axes([0, 0, 1, 1], frameon=False)
|
||||||
ax.set_xlim(*borders_x)
|
ax.set_xlim(*borders_x)
|
||||||
|
@ -35,22 +46,28 @@ ax.set_xticks([])
|
||||||
ax.set_ylim(*borders_y)
|
ax.set_ylim(*borders_y)
|
||||||
ax.set_yticks([])
|
ax.set_yticks([])
|
||||||
|
|
||||||
|
# Plot the initial values.
|
||||||
plot, = ax.plot(x_coords, y_coords, "b.")
|
plot, = ax.plot(x_coords, y_coords, "b.")
|
||||||
center_of_mass, = ax.plot(x_coords.mean(), y_coords.mean(), "r-")
|
center_of_mass, = ax.plot(x_coords.mean(), y_coords.mean(), "r-")
|
||||||
|
# Keep track of the center of mass.
|
||||||
center_of_mass_history_x = deque([x_coords.mean()])
|
center_of_mass_history_x = deque([x_coords.mean()])
|
||||||
center_of_mass_history_y = deque([y_coords.mean()])
|
center_of_mass_history_y = deque([y_coords.mean()])
|
||||||
|
|
||||||
brown = BrownIterator(-1, c
|
brown = BrownIterator(-1, c # Max iterations, simulation parameters.
|
||||||
, x_coords, y_coords
|
, x_coords, y_coords
|
||||||
, y_momenta, y_momenta
|
, y_momenta, y_momenta
|
||||||
|
# The boundary condition: reflect at the borders,
|
||||||
, borders_x, borders_y
|
, borders_x, borders_y
|
||||||
|
# or just let propagate to infinity.
|
||||||
#, [], []
|
#, [], []
|
||||||
|
# Let the border dampen the system, border_dampening < 1 => energy is absorbed.
|
||||||
, border_dampening=1
|
, border_dampening=1
|
||||||
, dt=dt)
|
, dt=dt)
|
||||||
|
|
||||||
u = iter(brown)
|
u = iter(brown)
|
||||||
|
|
||||||
def update(i):
|
def update(i):
|
||||||
|
# Get the next set of positions.
|
||||||
data = next(u)
|
data = next(u)
|
||||||
center_of_mass_history_x.append(x_coords.mean())
|
center_of_mass_history_x.append(x_coords.mean())
|
||||||
center_of_mass_history_y.append(y_coords.mean())
|
center_of_mass_history_y.append(y_coords.mean())
|
||||||
|
|
|
@ -42,8 +42,10 @@ class BrownIterator(object):
|
||||||
self.py = delta_py + self.py
|
self.py = delta_py + self.py
|
||||||
# XXX: We need the (-1)**i to make the problem
|
# XXX: We need the (-1)**i to make the problem
|
||||||
# symmetric.
|
# symmetric.
|
||||||
|
# FIXME: is this necessary?
|
||||||
self.px[np.isnan(self.px)] = self.speed_of_light * (-1)**self._i
|
self.px[np.isnan(self.px)] = self.speed_of_light * (-1)**self._i
|
||||||
self.py[np.isnan(self.py)] = self.speed_of_light * (-1)**self._i
|
self.py[np.isnan(self.py)] = self.speed_of_light * (-1)**self._i
|
||||||
|
|
||||||
self._reflect_at_borders()
|
self._reflect_at_borders()
|
||||||
|
|
||||||
self.x += self.dt * self.px
|
self.x += self.dt * self.px
|
||||||
|
|
Loading…
Reference in New Issue
Block a user