Продолжение из
этой темы. Рассмотрим частный пример движения частицы с зарядом
и массой
пусть
частица движется в плоскости, перпендикулярной вектору постоянного однородного магнитного поля
Орт декартовой системы координат
выберем вдоль
Тогда
радиус-вектор частицы есть
Уравнение Ньютона (в гауссовой системе единиц)
сводится к системе двух уравнений для координат
и
частицы:
где
При
частица вращалась бы по окружности постоянного радиуса с "циклотронной частотой"
Четыре постоянных интегрирования (начальные данные): координаты центра орбиты
и
, начальная фаза вращения
, радиус орбиты
Полагая для удобства, что начало координат помещено в центр орбиты (так что
можно один раз проинтегрировать уравнения с равными нулю постоянными интегрирования; такие же уравнения получатся для переменных
и
с произвольными
Поступим так же при
Вот получившиеся дифференциальные уравнения (ДУ) в записи с безразмерной переменной времени
и с безразмерным неотрицательным параметром радиационного трения
1. Сначала найдём точное решение, чтобы было с чем сравнивать приближённые ответы. "Методом комплексных амплитуд" ищем решение в виде
где комплексные постоянные
определяются системой уравнений
Такая система имеет не равные нулю решения
только тогда, когда её определитель равен нулю:
Решение этого уравнения для
можно выбрать так, что:
Обозначим:
- мнимая часть,
- вещественная часть частоты
т.е.
Положим
- это означает, что начальную фазу
выбираем равной нулю, а расстояние от начала координат до частицы при
выбираем за единицу длины. Тогда
Таким образом, вот точные выражения для траектории частицы c
Для перехода к исходным единицам достаточно умножить
и
на произвольный вещественный коэффициет
с размерностью длины и подставить
где
К аргументам синуса и косинуса можно добавить произвольное вещественное число
Если к получившимся
и
добавить произвольные значения координат центра орбиты
и
то в ДУ следует под координатами
и
понимать разности
и
В решении
уравнения
величина
имеет смысл частоты вращения, а
- параметр "затухания спирали"; оба параметра безразмерные, измеряются в долях от
Из уравнения для
следует, что
и
Видно, что в отсутствие трения, т.е. при
получается, как и должно быть,
а при небольшом трении, т.е. при
получается немного меньшее единицы значение
и маленькое положительное
Значит, с трением спираль будет не раскручиваться, а сжиматься, "затухать".
2. Теперь для упражнения попытаемся решать ту же систему ДУ приближённо.
Оказывается, простой метод конечных разностей ("метод Эйлера") даже в отсутствие трения плохо работает: вместо стационарной окружности у меня получалась постепенно раскручивающаяся спираль. Поэтому применяю более совершенный способ приближённых вычислений, "метод Рунге-Кутта". А именно:
Переходим к системе ДУ первого порядка c вспомогательными функциями
и
:
Шаг по безразмерному времени
обозначаю как
Приращения функций записываю в виде оборванного ряда по степеням
с точностью до членов порядка
включительно; например, для функции
имеем приближённое равенство:
где производные относятся к точке
Затем эти производные последовательно выражаем через сами функции в точке
с помощью наших ДУ; например,
дифференцируя это ещё раз и снова пользуясь ДУ, получаем выражение для
и так далее.
Значения функций в точке
пишем с индексом, например,
В точке
они тогда будут с индексом
То есть:
и так далее. Так получается система рекуррентных соотношений, позволяющих в цикле по
от нуля до какого-нибудь большого
вычислять интересующее нас множество значений координат
вместе со вспомогательными
имеющими смысл значений скорости
a) Посмотрим для начала, что получается при положительном шаге
с параметром трения
Пусть начальные координаты частицы:
При этом из приведённых выше формул точного решения (c
видно, что начальную скорость можно выбрать так:
где
и
определяются из решения уравнения для
Вычислить
поручаю Маткаду (но не знаю, насколько хорошо он с таким заданием справляется :). Вот маткадные фрагменты расчёта, рекуррентные формулы:
Ниже на графиках красная линия - точное решение, пунктир - результат приближённого вычисления по указанным рекуррентным формулам. Пунктир это способ изображения графика, а не сами вычисленные точки. Приближённых точек так же много, как и точных (
но если их все изобразить, то они полностью закрашивают собой предыдущий график. Результат при
и
Видно, что начало приближённого решения соответствует точному решению (пунктир идёт по красной линии), но затем быстро развивается неустойчивость приближённого расчёта. По физическому смыслу величина
должна быть малой. Но при меньших значениях
картина становится только хуже: неустойчивость проявляется ещё раньше.
б) Выполняем расчёт "назад во времени" - при отрицательном
с прежним значением параметра трения
Начальные значения координат и скорости в этом примере выбираю те, которые получаются из точного решения для конечного момента времени
т.е. для точки с номером
После завершения цикла вычислений (с прежними рекуррентными формулами, по
от нуля до
но теперь с отрицательным шагом
все получившиеся координаты
перенумеровываем в обратном порядке. И... ура, неустойчивость исчезла, результат соответствует точному решению:
"Затухающая спираль" на плоскости
в этом примере выглядит так:
в) Расчёт "назад во времени", при
с прежним значением параметра трения
но с произвольными начальными условиями:
Начальные условия теперь задаю как попало; например, пусть
По отношению к нумерации "вперёд во времени" эти условия конечные - они для точки
с номером
После завершения такого же, как и раньше, цикла вычислений по
от нуля до
с отрицательным шагом
получившиеся координаты
перенумеровываем в обратном порядке. В результате, получившиеся
относятся уже к точке
По ним определяем начальный радиус
и начальную фазу
для сравнения с точным решением
Оказывается, в этом случае результат приближённого вычисления тоже соответствует точному решению:
Лишь на первых шагах расчёта "назад во времени" (т.е. вблизи
в этом примере) заметно обусловленное произвольными начальными условиями большое отличие приближённого результата от точного решения:
г) Расчёт "назад во времени" при
с малым значением параметра трения:
Пусть максимальное
будет равно
так что
Начальные условия (т.е. условия для точки
с номером
задаю произвольно:
Сравнение с точным решением делается тем же способом, что и в предыдущем примере; только теперь
в пятьдесят раз больше. Цель этого теста - убедиться, что даже на длинном куске траектории приближённый расчёт не сбивается, согласуется с точным решением. Оказывается, так и есть. Вот часть траектории в начале временного интервала:
Последняя часть траектории:
Вся "спираль" на плоскости
в этом примере выглядит так:
Маткад сумел даже кое-как анимировать примерно 80 оборотов по этой траектории (здесь
ссылка на avi-файл, некоторые браузеры могут его показывать, но не знаю, все ли ).
О перенумеровке точек траектории в обратном порядке:
Есть варианты. Вычисленные по указанным в начале рекуррентным формулам (c
массивы
можно в цикле по
от нуля до
переписать в новые массивы, например,
а затем, чтобы вернуться к исходным обозначениям, в отдельном цикле записываем их на место исходных величин:
и т.д. Более простой вариант: в рекуррентных формулах можно дать массивам новые обозначения, например,
вместо
и т.д.; тогда для перенумеровки в обратном порядке достаточно в цикле по
от нуля до
записать
Наконец, самый простой вариант - вести рекуррентные вычисления сразу в обратном порядке, т.е. в цикле по
от
до
с шагом
вычислять
и т.д.; тогда не надо будет ничего перенумеровывать.