Как создать граничное условие отражения «бильярдный шар» в питоне?

Question

Согласно Эрвину Шрёдингеру («Что такое жизнь?»), Диффузия может быть полностью объяснена случайным движением частиц. Я хочу проверить это сам, создав программу, которая создает временную визуализацию диффузии "молекул газа" в закрытом контейнере. Начальные условия будут иметь два раздела, один с низким и один с высокой концентрацией. После t0 перегородка удаляется, и газу позволяют диффундировать. Единственный механизм, который я хочу использовать, - это добавление случайных векторов смещения к каждой молекуле. Начальные условия будут выглядеть следующим образом.

Часть проблемы, которую я не понимаю, состоит в том, как создать простое отражение типа бильярда, когда молекула попадает на ограничивающие поверхности. Я предполагаю простые симметричные отражения (угол in = угол вне на границах). Я вообще не запускал код, потому что я не знаю, как справиться с этой частью, хотя я знаю, как сделать остальную часть. Я знаю, что это больше математический вопрос, но как я могу создать эти граничные условия в Python? В идеале мне бы хотелось запрограммировать эту функцию самостоятельно, чтобы я мог ее понять, а не использовать предварительно созданный пакет, который мог бы это сделать. Это то, что я ищу, для любой данной молекулы.

В конечном счете, все, что мне действительно нужно, это следующее: учитывая начальное местоположение (x1, y2), векторную величину v, угол тета, а также размер и местоположение прямоугольника, каково окончательное положение покоя молекулы (x2, у2).

Answer 1

Вам не нужно вычислять угол отражения, просто разбейте задачу на две части: одну для x и одну для y. В обоих случаях вам нужно, чтобы частица «возвращалась», когда она выходит за границы.

Я сделал это время назад для изучения плотности частиц в жидкостях. Проще всего рассмотреть (0, 1) границу в обоих направлениях. Следующий код должен сделать это (совет: правильное использование abs создаст эквивалент отражения):

x0 = [.1, .9]
delta = [-0.2, 0.3]
x1 = [(1-abs(abs(xi + di)-1)) for xi, di in zip(x0, delta)]
print(x1)
# 0.1, 0.8
#or using numpy:
x1 = 1-np.abs(np.abs(np.asarray(x0) + np.asarray(delta))-1)
print(x1)
>> [0.09999999999999998, 0.8]
   array([0.1, 0.8])

Я предполагаю из вашего вопроса, что вы пренебрегаете столкновением частиц и частиц и "не-суперпозицией" частиц-частиц *

Answer 2

Вот простая реализация. Я изменяю векторы движения только каждый десятый шаг, так что можно визуально проверять граничные отражения. Частицы мигают красным при обновлении векторов движения.

Хитрость, как описано *ere - «развернуть» ограничивающий прямоугольник. Вместо этого мы позволяем частицам двигаться без ограничений, а затем складываем пространство в ограничительную рамку.

import numpy as np
import pylab
from matplotlib.animation import FuncAnimation

xy = np.random.uniform(-1, 1, (2, 200))
xy[0, :160] = np.abs(xy[0, :160])
xy[0, 160:] = -np.abs(xy[0, 160:])
xy += 1

f, a = pylab.subplots()
pxy, = pylab.plot(*xy, 'o')

def init():
    a.set_xlim(0, 2)
    a.set_ylim(0, 2)
    return pxy,

def update(frame):
    global inc, xy
    if frame % 1 < 0.01:
        inc = np.random.normal(0, 0.01, xy.shape)
        pxy.set_markerfacecolor('red')
    elif frame % 1 < 0.11:
        pxy.set_markerfacecolor('blue')        
    xy += inc
    fxy = np.abs((xy+2)%4-2)
    pxy.set_data(*fxy)
    return pxy,

anim = FuncAnimation(f, update, frames=np.arange(1200) / 10,
                     init_func=init, blit=True)

pylab.show()

Answer 3

Итак, пару вещей, которые нужно иметь в виду:

1. Вам нужен компонент трения, иначе частица будет продолжать двигаться вечно (сохранение энергии). В этом случае трение происходит как функция скорости, и трение также будет происходить при скачках.

2. Если это всего лишь одна частица, вы можете рассчитать ее, определив ограничивающие прямоугольники, например, x между 0 и 5, y, например, между 0 и 3. Затем вы можете вычислить точку пересечения со стеной, вставив значение x = 5, а затем решив для y в уравнении для линии.

Для одной частицы вам не нужно делать это параметрически с шагом t_0, вы можете вычислить перехваты и в основном увеличить их там. Для кратных вам нужно будет рассчитать межмолекулярную диффузию и силы столкновения ... что является гораздо более сложной задачей, которая должна выполняться параметрически.

Вы должны будете рассчитать столкновения, когда центры двух молекул находятся на расстоянии 2 * радиуса друг от друга, а затем выполнить столкновение, которое сохраняет импульс .