Ричард Фейнман - Фейнмановские лекции по физике 1. Современная наука о природе, законы механики
Заметьте, что новое положение тела равно старому плюс интервал времени ?, умноженный на скорость. Но что это за скорость? В какой момент? В начале интервала одна скорость, а в конце она совсем другая. Прием состоит в том, чтобы брать скорость в середине интервала. Если известна скорость в настоящий момент и известно, что она меняется, как же можно надеяться получить удовлетворительный результат, считая, что тело все время движется с той же скоростью, что и в настоящий момент? Более разумно использовать какую–то среднюю скорость между началом и концом интервала. Те же рассуждения применимы к изменению самой скорости: для подсчета ее изменений нужно использовать ускорение в средней точке между двумя моментами времени, в которых необходимо найти скорость. Таким образом, реально мы будем пользоваться следующими уравнениями: положение в конце интервала равно положению в начале плюс интервал ?, умноженный на скорость в середине интервала. Эта скорость в свою очередь равна скорости в середине предыдущего интервала (т. е. на отрезок ? меньше) плюс ускорение в начале интервала, умноженное на ?.
Таким образом, мы будем пользоваться уравнениями
Остается еще один небольшой вопрос: что такое v (?/2)? Вначале у нас было v (0), а не v (-?/2). Но теперь, чтобы начать наши вычисления, необходимо использовать дополнительное уравнение v(?/2)=v (0)+( ?/2)а(0).
Таблица 9.1 • решение уравнения (dvx/dt)=-x Интервал ?=0,10 сек
Ну, а теперь все готово для расчетов. Для удобства можно их выполнить в виде таблицы, в столбцах которой стоят время, положение, скорость и ускорение, причем скорость пишется в промежутках между строками (табл. 9.1). Такая таблица есть, конечно, просто удобный способ записи результатов, полученных из уравнений (9.16), и фактически полностью заменяет их. Мы просто заполняем одно за другим свободные места в ней и получаем очень интересную картину движения: сначала грузик находится в покое, затем понемногу приобретает отрицательную скорость (вверх), а это приводит к уменьшению его расстояния от точки равновесия. При этом хотя ускорение и становится меньше, оно все еще «подгоняет» скорость. Однако по мере приближения к положению равновесия (х=0) ускорение становится все меньше и меньше, скорость нарастает все медленней и медленней, но все же еще нарастает вплоть до точки x=0, которая достигается примерно через 1,5 сек. Скажем по секрету, что произойдет дальше. Грузик, конечно, не остановится в точке х=0, а пойдет дальше, но теперь все пойдет наоборот: его положение х станет отрицательным, а ускорение – положительным. Скорость начнет уменьшаться. Интересно сравнить полученные нами числа с функцией cost. Результат этого сравнения представлен на фиг. 9.4.
Фиг. 9.4. График движения грузика на пружинке.
Оказывается, что в пределах точности наших расчетов (три знака после запятой) совпадение полное! Позднее вы узнаете, что функция cos t – точное решение нашего уравнения, так что у вас теперь есть наглядное представление о мощи численного анализа: столь простой расчет дает столь точный результат.
§ 6. Движение планет
Приведенный анализ очень подходит к движению осциллирующей пружинки с грузиком, но можно ли таким же путем вычислять движение планеты вокруг Солнца? Давайте посмотрим, можно ли при некоторых приближениях получить эллиптическую орбиту. Предположим, что Солнце бесконечно тяжелое в том смысле, что его движение не будет приниматься в расчет.
Допустим, что в известной точке планета начала свое движение и имеет определенную скорость. Она движется вокруг Солнца по какой–то кривой, и мы попытаемся определить с помощью уравнений движения Ньютона и его же закона всемирного тяготения, что это за кривая. Как это сделать? В некоторый момент времени планета находится в каком–то определенном месте, на расстоянии r от Солнца; в этом случае известно, что на нее действует сила, направленная по прямой к Солнцу, которая, согласно закону тяготения, равна определенной постоянной, умноженной на произведение масс планеты и Солнца и деленной на квадрат расстояния между ними. Чтобы рассуждать дальше, нужно выяснить, какое ускорение вызывает эта сила.
Однако в отличие от предыдущей задачи нам потребуются теперь компоненты ускорения в двух направлениях, которые мы назовем х и у. Положение планеты в данный момент будет определяться координатами х и у, поскольку третья координата z всегда равна нулю.
Действительно, координатная плоскость ху выбрана нами таким образом, что z–компоненты как силы, так и начальной скорости равны нулю, а поэтому нет никаких причин, которые бы заставили планету выйти из этой плоскости. Сила при этом будет направлена по линии, соединяющей планету с Солнцем, как это показано на фиг. 9.5.
Фиг. 9.5. Сила притяжения, действующая на планету.
Из этого рисунка видно, что горизонтальная компонента силы так относится к полной ее величине, как координата х относится к расстоянию r. Это сразу следует из подобия треугольников. Кроме того, если х положительна, то Fx отрицательна, и наоборот.
Таким образом, Fx?F?=-x/r, или Fя=-?F?xlr=-GM mx/r3 и соответственно Fy=-GMmy/r3. Теперь можно воспользоваться динамическими законами (9.7) и написать, что х–или y–компонента ускорения, умноженная на массу планеты, равна соответственно х–или y–компоненте силы:
Это именно та система уравнений, которую мы должны решить. Для того чтобы упростить вычисления, предположим, что либо единицы измерения времени или массы выбраны соответствующим образом, либо нам просто повезло, словом, получилось так, что GM=1. Для нашего случая предположим, что в начальный момент t=0 планета находилась в точке с координатами х=0,500 и у=0,000, а скорость ее в этот момент направлена параллельно оси у и равна 1,6300. Как же в этом случае делаются расчеты? Снова составляется таблица со столбцами для времени t, координаты х, x–компонент скорости vx и ускорения ах. Затем идут отделенные чертой три колонки: для координаты y, у–компонент скорости и ускорения. Однако, для того чтобы подсчитать ускорения, мы должны воспользоваться уравнением (9.17), согласно которому его компоненты равны –х/r3 и –у/r3, а r=?(x2+y2). Так что, получив х и у, мы должны где–то в сторонке провести небольшие вычисления – извлечь квадратный корень из суммы квадратов и получить расстояние. Удобно также отдельно вычислить и 1/r3.
После этого все готово, чтобы определить компоненты ускорения. Всю эту работу можно сильно облегчить, если пользоваться таблицами квадратов, кубов и обратных величин. На нашу долю останется тогда только умножение х на 1/r3, которое легко выполняется на логарифмической линейке.
Перейдем к дальнейшему. Возьмем интервал времени ?=0,100. В начальный момент t=0
x(0)=0,500,. у(0)=0,000,
vx(0) = 0,000, vy(0)=+1,630.
Отсюда находим
r(0)=0,500, 1/r3=8,000,
ax=-4,000, ау=0,000.
После этого можно вычислять компоненты vx (0,05) и vy (0,05):
vя (0,05)=0,000–4,000•0,050 = -0,200,
vy(0,05)=1,630+0,000–0,100=1,630.
А теперь начнем наш основной расчет:
и т. д.
В результате мы получим числа, приведенные в табл. 9.2, где приблизительно за 20 шагов прослежена половина пути нашей планеты вокруг Солнца. На фиг. 9.6 отложены координаты планеты х и y, приведенные в табл. 9.2.
Фиг. 9.6. График движения планеты вокруг Солнца.
Точки представляют собой последовательные положения планеты через каждую десятую долю выбранной нами единицы времени. Видно, что сначала она двигалась быстро, а затем – все медленней и медленней. Видна также и форма кривой движения планеты. Итак, вы теперь знаете, как реально можно вычислять движение планет!
Давайте посмотрим теперь, как вычислить движение Нептуна, Юпитера, Урана и остальных планет. Можно ли сделать подробные расчеты со множеством планет, учитывая к тому же и движение Солнца? Разумеется, можно. Найдем сначала силу, действующую на каждую данную планету, например на ту, которую мы обозначим номером i и координаты которой хi, yi и zi (i=1 может означать Солнце, i=2 – Меркурий, i=3 – Венеру и т. д.). Наша задача – найти координаты всех планет. По закону тяготения x–компонента силы, действующая на i–ю планету со стороны планеты номер j' с координатами хj уj, zj, будет равна –Gmimj(xi–xj)/r3jj . Если же учесть силы со стороны всех планет, то получим следующую систему уравнений:
где rij – расстояние между i–й и j–й планетами:
? означает суммирование по всем остальным планетам, r. е. по всем значениям j, за исключением, конечно, j = i. Таким образом, чтобы решить это уравнение, нужно лишь значительно увеличить количество столбцов в нашей таблице. Для движения Юпитера понадобится девять столбцов, для Сатурна – тоже девять и т. д. Если нам заданы все начальные положения и скорости, то из уравнения (9.18) можно подсчитать все ускорения, вычислив, конечно, предварительно по формуле (9.19) все расстояния rij,. А сколько же времени потребуется на все эти вычисления? Если вы будете делать их сами дома, то очень много! Однако сейчас уже имеются машины, неимоверно быстро выполняющие все арифметические расчеты. Сложение, например, такая машина выполняет за 1 мксек, т. е. за одну миллионную долю секунды, а умножение – за 10 мксек. Так что если один цикл расчетов состоит из 30 операций умножения, то это займет всего лишь 300 мксек, или за 1 сек можно сделать 3000 циклов. Если мы хотим считать с точностью до одной миллиардной, то для того, чтобы покрыть все время обращения планеты вокруг Солнца, требуется 4•105 циклов. (Оказывается, что ошибка в расчетах приблизительно пропорциональна квадрату ?. Если брать интервал в тысячу раз меньший, то ошибка уменьшится в миллион раз. Так что для обеспечения нашей точности нужно взять интервал в 10 000 раз меньше.) На машине это займет 130 сек, или около 2 мин. Всего лишь 2 мин, для того чтобы «прогнать» Юпитер вокруг Солнца и при этом еще с точностью до одной миллиардной учесть все возмущения от других планет!
