Квант

2

Среди многоих применений компьютеров в науке и техническе одно из самых важных и распространенных — моделирование. Создатся компьютерные модели самолетов и морских течений, движения планет и экономических связей между странами…

В этой статье пойдет речь об одной из простейших модели газа. Это «идеальный бильярд»: абсолотно гладкие и упругие шары играют роль молекул газа, бортики бильярдного стола — роль стенок сосуда, поверхность стола абсолотно гладкая. Такого идеального бильярда не существует, однако можно «создать» его с помощью компьютера и наблюдать «жизнь» шаров на экране графического дисплея.

Возможно, «идеальный бильярд» показется вам очень прощеной модели газа. В принципие, на компьютере можно «сконструировать» и более реализную модель — перейти к движению шаров в пространстве, а не на плоскости, ввести взаимодействие шаров-молекул и т. п. Оdnako и «бильярд» — объект не только учебных, студенческих, но и серверных научных исследований.

Ясно, что движение идеальных шаров полностью определяется законами механики. Между тем движение части, составляющих газ, может рассматриваться какxaотическая, случайное*).

В этой статье мы опишем эксперименты с моделью «идеальный бильярд», которые позволяют наблюдать, как движение становитсяxaотическая, и попытаться выяснить при-

*) Картину возникновения «молекулярногоxaosa» в газе шаров исследовал советский физик Н. С. Крылов (1950 г.). Советский математик Я. Г. Синай провел подобное исследование на уровень строности, удовлетворяющим математиков (1963 г.). Научных работ, использованных компьютерную модель «газа шаров», очень много.

Рис. 1.

3

ра его о стенку) — это будут новые начальные данные — и повторим расчет.

Положения шаров будем задаать с помощью декартовых координат их центров, а скорости — с помощью их проекций на те же оси (рис. 1):

r⃗=(x,y),v⃗=(vx,vy).\vec{r} = (x, y), \quad \vec{v} = (v_x, v_y).

Формулы для вычисления оказываются совсем несложными. Пусть ось Y совпадает с левым бортиком бильярда. Тогда время до стокновения шара с этим бортиком —

t=−(x−a)/vx(a — радиус шара),t = -(x-a)/v_x \quad (a \text{ — радиус шара}),

причем стокновение состоит только если vx<0v_x < 0. Пусть это стокновение и есть само решение. При стокновении проекция vyv_y остается неизменной, а проекция vxv_x меняет знак:

vx→−vx,vy→vy.v_x \rightarrow -v_x, \quad v_y \rightarrow v_y.

Перемещение каждого шара за время tt —

r⃗it=r⃗i+v⃗it,\vec{r}_{it} = \vec{r}_i + \vec{v}_i t,

т. е. xit=xi+vixt,yit=yi+viytx_{it} = x_i + v_{ix} t, \quad y_{it} = y_i + v_{iy} t. Рассеты для случая стокновения шаров, хотя и несколько сложнее, тоже выполня «школьные». Предлагаем вам проверсти их самостояльно (см. задачи в конце статьи).

Описание полной программы, моделирующей движение NN шаров на графикском дисплее, приведено в статье «Газ на дисплее» (см. с. 61).

Какова скорость расчетов по подобной схеме? Компьютеры, использемые в студенческом практике физического факультета Новосибирского государственного университета, успешют считать с такой скоростью, что кадры на дисплее меняются так же быстро, как в мультфильме (если шаров не более десяти; с увеличением числа шаров время расчета каждого кадра «мультфильма» быстро расчет)*.

Все вычисления наш компьютер ведет с 7—8 значцами цифрами.

*) Счет идет на ЭВМ «Электронка-60», включенных в «Терминальный руководство системы» (ТЕВУС — разработка Новосибирского государственного университета).

Рис. 2.

4

ры — сначала один, а затем и остальные — «сбились» с пути и пришли совсем не туда, откуда вышили. Таким образом, движение шаров оказывается «неправильным». Почему же?

Весь расчет поизводится по точным формулам, поэтому единственная возможность причина «неправильности» — округление чисел в процессе счета. Но отклонение, связанное с округлением, имеет относительной величины 10−7−10−810^{-7} - 10^{-8} (как мы говорили, компьютер работает с 7—8 значениями цифрами), а на экране к моменту времени 2t02t_0 шары совсем «забыли» свои исходные положения. Можно ли эту «забывчность» объяснить округлением?

Оценки отклонений

Чтобы получить количественные оценки, представим себе, что из-за округления расчет направление движения одного из шаров отклонилось отправильно на угол φ0∼10−7\varphi_0 \sim 10^{-7}, а дальнейшее движение расчитывается точно.

Сначала рассмотрим случай, когда имеется всего один шар и отклонение произходит при ударе о стенку. Траектория центра шара изображена на рисунке 3, и там же приведено несколько построения, облегчающей расчет. Квадрат со стороны (b−2a)(b - 2a) — область, доступная центру шара (b−dлина стононы бильярного стола)(b - \text{dлина стононы бильярного стола}). Длина траектории 1→2→3→4→5→…1 \rightarrow 2 \rightarrow 3 \rightarrow 4 \rightarrow 5 \rightarrow \dots равна длине луча 1→2→3′→4′→5′→…1 \rightarrow 2 \rightarrow 3' \rightarrow 4' \rightarrow 5' \rightarrow \dots. Если в точке 1 в результате округления значения vxv_x направление скорости шара после удара отклонилось от точного направления на малый угол φ0\varphi_0, то на пути длиной LL «накопится» отклонение центра шара на расстояние Δ∼Lφ0\Delta \sim L\varphi_0. Ограничившись временем TT, за которое произойдет v1=1000v_1 = 1000 столковенный со стенкой, мы получим отклонение Δ1∼v1bφ0∼10−4b\Delta_1 \sim v_1 b\varphi_0 \sim 10^{-4}b. Для наших графических дисплеев толщина линии на экране составляет примерно 0,004b0,004b. Так что за время TT накопившееся отклонение не может быть замещенным — оно меньше толщины линии на экране. На заметное расстояние, равное, скажем, диаметр шара 2a2a, центр шара отклонится через v2∼2a/bφ0v_2 \sim 2a/b\varphi_0 столковенный. Время, за которое произойдет такое количество ударов о стенки, оценивается в 100 (!) часов.

Очевидно, округлением при расчету ударов о стенки нельзя обязывать «неправильное» движение шаров в нашем эксперименте.

Посмотрим, как сказывается округление чисел при расчете соударений шаров. Вновь применем, что направление движения одного из шаров отклонилось отправильно на угол φ0∼10−7\varphi_0 \sim 10^{-7}, и исследуем, как будет меньше угол отклонения при соударениях. При этом ограничится самыми группы оценками. Будем иметь в виду, что между соударениями шар движется по прямой и проходит в среднем пути ll, много больший радиуса шара (l≫a)(l \gg a). За время между двумя последовательными столковениями центр шара сместится на расстояние OO′∼lφ0OO' \sim l\varphi_0 (рис. 4). Тогда же порядка и смесение точки соприкосновения шаров при ударе. Участок поверхности второго шара в окрестности точки касания играет при ударе роль «зеркала». Это «зеркало» само движалось до удара, при ударе оно отскакивает, поэтому направление, в каком отскакивает наш шар, не определяется правильом «угол падения равен углу отражения». Тем не менее поворот «зеркала» на малый

5

угодл α ведет к изменению направления движения отскочившего шара на угол φ₁ ~ α (может быть, в полтора-два раза больше или меньше — такой уровень точности в данном случае насSTRAИВАЕТ). Так как α ~ (AA’/a) ~ (lφ₀/a), получаем

φ₁ ~ (l/a) φ₀ >> φ₀

— при ударе отклонение скорости резко возросло. После k стольновейный φₖ ~ (l/a) k φ₀.

Если (l/a) ~ 10, то достаточно 8—10 стольновейный, чтобы стало φₖ ~ 1 и направление движения шара перестало иметь какое бы то ни было относение к правильному.

КонтROLНЫЙ ЭКСПЕРИМЕНТ

Проделанные оценки выглядят убедительно, и тем не менее подвергнем их проверке с помощью компьютерного эксперимента. Для этого проведем два «запуска» шаров с почти одинаковыми начальными условиями, лишь слегка изменив во второй раз начальную скорость (или начальное положение) одного из шаров.*) «Закажем» компьютеру запоминать координаты какого-то шара через равные интервалы времени Δt; при втором запуске будем отмечать координаты того же шара в те же моменты tₖ = k · Δt, а затем определим, как изменяется со временем «расхождения» траекторий этого шара

Δrₖ = √((xₖ - x’ₖ)² + (yₖ - y’ₖ)²),

где xₖ, yₖ — значения координат в моменты tₖ при первом запуске, а x’ₖ, y’ₖ — при втором. Учитывая ожидаемую зависимость, удобно воспользоваться логарифмическим масштабом по оси Δrₖ, т. е. откладывать по этой оси значения lg (Δrₖ). Зависимость вида φₖ ~ (l/a) k φ₀ должна привести к прямолинейному графику

lg (Δrₖ) ~ lg φ₀ + (t/tₖ) lg (l/a).

Результаты экспериментов приведены на рисунке 5. На одном графике изображены точки, полученные при разных отклонениях в начальных условиях. Ожидаемый прямолинейный участок графика виден довольно отчетливо. Спустя время ℭ величина Δrₖ вырастает до предела, задаваемого размером стола, и далее изменяется хаотическая, но, в общем-то, не в миллионы раз (как на участке t < ℭ), так что в логарифмическом масштабе точки при t > ℭ не отклоняются даже от прямой lg (Δrₖ) ~ lg l.

Этот комплютерный эксперимент свидетельствует о разумности сделанных ранее оценок.

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