97 lines
2.9 KiB
Python
97 lines
2.9 KiB
Python
from brown.interaction import UFuncWrapper
|
|
from brown.brown import BrownIterator
|
|
import numpy as np
|
|
from collections import deque
|
|
from copy import copy
|
|
import matplotlib.pyplot as plt
|
|
import matplotlib.animation as ani
|
|
|
|
from coefficients import c
|
|
|
|
#force_function = UFuncWrapper(0, 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_y = [-100, 100]
|
|
n_particles = 60
|
|
# Idk, seems to not do anyting.
|
|
frames = 100
|
|
# Only spawn in 1/x of the borders.
|
|
spawn_restriction = 2
|
|
# 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.01
|
|
c[-1] = dt
|
|
|
|
# Draw arrows, makes the simulation fucking slow.
|
|
draw_arrows = True
|
|
|
|
# Initial positions.
|
|
x_coords = np.random.uniform(borders_x[0] / spawn_restriction, borders_x[1] / spawn_restriction, n_particles).astype(np.float32)
|
|
y_coords = np.random.uniform(borders_y[0] / spawn_restriction, borders_y[1] / spawn_restriction, n_particles).astype(np.float32)
|
|
|
|
|
|
# Initial momenta are 0.
|
|
x_momenta = np.zeros(n_particles, dtype=np.float32)
|
|
y_momenta = np.zeros(n_particles, dtype=np.float32)
|
|
|
|
|
|
|
|
# Prepare the plot, remove axis & stuff.
|
|
fig = plt.figure(figsize=(7, 7))
|
|
ax = fig.add_axes([0, 0, 1, 1], frameon=False)
|
|
ax.set_xlim(*borders_x)
|
|
ax.set_xticks([])
|
|
ax.set_ylim(*borders_y)
|
|
ax.set_yticks([])
|
|
|
|
# Plot the initial values.
|
|
plot, = ax.plot(x_coords, y_coords, "b.")
|
|
center_of_mass, = ax.plot(x_coords.mean(), y_coords.mean(), "r-")
|
|
if(draw_arrows):
|
|
arrows = [ plt.Arrow(x, y, dx, dy) for x, y, dx, dy in zip(x_coords, y_coords, x_momenta, y_momenta)]
|
|
for arrow in arrows:
|
|
ax.add_patch(arrow)
|
|
|
|
|
|
# Keep track of the center of mass.
|
|
center_of_mass_history_x = deque([x_coords.mean()])
|
|
center_of_mass_history_y = deque([y_coords.mean()])
|
|
|
|
brown = BrownIterator(-1, c # Max iterations, simulation parameters.
|
|
, x_coords, y_coords
|
|
, y_momenta, y_momenta
|
|
# The boundary condition: reflect at the borders,
|
|
, 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
|
|
, dt=dt)
|
|
|
|
u = iter(brown)
|
|
|
|
def update(i):
|
|
global arrows
|
|
# Get the next set of positions.
|
|
data = next(u)
|
|
center_of_mass_history_x.append(x_coords.mean())
|
|
center_of_mass_history_y.append(y_coords.mean())
|
|
|
|
plot.set_data(*data)
|
|
center_of_mass.set_data(center_of_mass_history_x, center_of_mass_history_y)
|
|
if(draw_arrows):
|
|
for arrow in arrows:
|
|
ax.patches.remove(arrow)
|
|
arrows = [plt.Arrow(x, y, dx, dy) for x, y, dx, dy in zip(data[0], data[1], u.px, u.py)]
|
|
for arrow in arrows:
|
|
ax.add_patch(arrow)
|
|
|
|
animation = ani.FuncAnimation(fig, update, range(frames), interval=1)
|
|
plt.show()
|
|
#animation.save("animation.mp4", fps=30)
|