Симплектический интегратор - Symplectic integrator

В математика, а симплектический интегратор (СИ) это схема численного интегрирования за Гамильтоновы системы. Симплектические интеграторы образуют подкласс геометрические интеграторы которые по определению канонические преобразования. Они широко используются в нелинейная динамика, молекулярная динамика, методы дискретных элементов, физика ускорителя, физика плазмы, квантовая физика, и небесная механика.

Вступление

Симплектические интеграторы предназначены для численного решения Уравнения Гамильтона, который читается

${\displaystyle {\dot {p}}=-{\frac {\partial H}{\partial q}}\quad {\mbox{and}}\quad {\dot {q}}={\frac {\partial H}{\partial p}},}$

где ${\displaystyle q}$ обозначает координаты положения, ${\displaystyle p}$ координаты импульса, и ${\displaystyle H}$ гамильтониан. Набор координат положения и импульса ${\displaystyle (q,p)}$ называются канонические координаты.(Видеть Гамильтонова механика для получения дополнительной информации.)

Временная эволюция Уравнения Гамильтона это симплектоморфизм, что означает, что он сохраняет симплектическую 2-форма ${\displaystyle dp\wedge dq}$. Численная схема является симплектическим интегратором, если она также сохраняет эту 2-форму.

Симплектические интеграторы обладают, как сохраняющейся величиной, гамильтонианом, который возмущенный от оригинала. Благодаря этим преимуществам схема SI широко применяется для расчетов долговременной эволюции хаотических гамильтоновых систем в диапазоне от Проблема Кеплера к классическому и полуклассическому моделированию в молекулярная динамика.

Большинство обычных численных методов, таких как примитив Схема Эйлера и классический Схема Рунге – Кутты, не являются симплектическими интеграторами.

Методы построения симплектических алгоритмов

Методы расщепления сепарабельных гамильтонианов

Широко распространенный класс симплектических интеграторов составляют методы расщепления.

Предположим, что гамильтониан сепарабелен, то есть его можно записать в виде

${\displaystyle H(p,q)=T(p)+V(q).\qquad \qquad (1)}$

Это часто случается в гамильтоновой механике, где Т будучи кинетическая энергия и V то потенциальная энергия.

Для простоты обозначений введем символ ${\displaystyle z=(q,p)}$ для обозначения канонических координат, включая координаты положения и импульса. Тогда система уравнений Гамильтона, приведенная во введении, может быть выражена одним выражением как

${\displaystyle {\dot {z}}=\{z,H(z)\},\qquad \qquad (2)}$

где ${\displaystyle \{\cdot ,\cdot \}}$ это Скобка Пуассона. Кроме того, введя оператор ${\displaystyle D_{H}\cdot =\{\cdot ,H\}}$, который возвращает Скобка Пуассона операнда с Гамильтониан, выражение уравнения Гамильтона можно упростить до

${\displaystyle {\dot {z}}=D_{H}z.}$

Формальное решение этой системы уравнений дается как матричная экспонента:

${\displaystyle z(\tau )=\exp(\tau D_{H})z(0).\qquad \qquad (3)}$

Обратите внимание на положительность ${\displaystyle \tau D_{H}}$ в матрице экспоненты.

Когда гамильтониан имеет форму ур. (1) решение (3) эквивалентно

${\displaystyle z(\tau )=\exp[\tau (D_{T}+D_{V})]z(0).\qquad \qquad (4)}$

Схема SI аппроксимирует оператор временной эволюции ${\displaystyle \exp[\tau (D_{T}+D_{V})]}$ в формальном решении (4) произведением операторов в виде

${\displaystyle {\begin{aligned}\exp[\tau (D_{T}+D_{V})]&=\prod _{i=1}^{k}\exp(c_{i}\tau D_{T})\exp(d_{i}\tau D_{V})+O(\tau ^{k+1})\\&=\exp(c_{1}\tau D_{T})\exp(d_{1}\tau D_{V})\dots \exp(c_{k}\tau D_{T})\exp(d_{k}\tau D_{V})+O(\tau ^{k+1}),\end{aligned}}\qquad \qquad (5)}$

где ${\displaystyle c_{i}}$ и ${\displaystyle d_{i}}$ настоящие числа, ${\displaystyle k}$ целое число, которое называется порядком интегратора, и где ${\displaystyle \sum _{i=1}^{k}c_{i}=\sum _{i=1}^{k}d_{i}=1}$. Обратите внимание, что каждый из операторов ${\displaystyle \exp(c_{i}\tau D_{T})}$ и ${\displaystyle \exp(d_{i}\tau D_{V})}$ обеспечивает симплектическая карта, поэтому их произведение, фигурирующее в правой части (5), также составляет симплектическое отображение.

С ${\displaystyle D_{T}^{2}z=\{\{z,T\},T\}=\{({\dot {q}},0),T\}=(0,0)}$ для всех ${\displaystyle z}$, можно сделать вывод, что

${\displaystyle D_{T}^{2}=0.\qquad \qquad (6)}$

Используя ряд Тейлора, ${\displaystyle \exp(aD_{T})}$ можно выразить как

${\displaystyle \exp(aD_{T})=\sum _{n=0}^{\infty }{\frac {(aD_{T})^{n}}{n!}},\qquad \qquad (7)}$

где ${\displaystyle a}$ - произвольное действительное число. Комбинируя (6) и (7) и используя те же рассуждения для ${\displaystyle D_{V}}$ как мы использовали для ${\displaystyle D_{T}}$, мы получаем

${\displaystyle \left\{{\begin{aligned}\exp(aD_{T})&=1+aD_{T},\\\exp(aD_{V})&=1+aD_{V}.\end{aligned}}\right.\qquad \qquad (8)}$

Конкретно, ${\displaystyle \exp(c_{i}\tau D_{T})}$ дает отображение

${\displaystyle {\begin{pmatrix}q\\p\end{pmatrix}}\mapsto {\begin{pmatrix}q+\tau c_{i}{\frac {\partial T}{\partial p}}(p)\\p\end{pmatrix}},}$

и ${\displaystyle \exp(d_{i}\tau D_{V})}$ дает

${\displaystyle {\begin{pmatrix}q\\p\end{pmatrix}}\mapsto {\begin{pmatrix}q\\p-\tau d_{i}{\frac {\partial V}{\partial q}}(q)\\\end{pmatrix}}.}$

Обратите внимание, что обе эти карты практически вычислимы.

Примеры

Упрощенная форма уравнений (в выполненном порядке):

${\displaystyle q_{i+1}=q_{i}+c_{i}{\frac {p_{i+1}}{m}}t}$

${\displaystyle p_{i+1}=p_{i}+d_{i}F(q_{i})t}$

После преобразования в лагранжевые координаты:

${\displaystyle x_{i+1}=x_{i}+c_{i}v_{i+1}t}$

${\displaystyle v_{i+1}=v_{i}+d_{i}a(x_{i})t}$

Где ${\displaystyle F(x)}$ вектор силы при ${\displaystyle x}$, ${\displaystyle a(x)}$ - вектор ускорения при ${\displaystyle x}$, и ${\displaystyle m}$ - скалярная величина массы.

Ниже приведены некоторые симплектические интеграторы. Наглядный способ их использования - рассмотреть частицу с положением ${\displaystyle q}$ и скорость ${\displaystyle p}$.

Чтобы применить временной шаг со значениями ${\displaystyle c_{1,2,3},d_{1,2,3}}$ к частице выполните следующие действия:

Итеративно:

• Обновить позицию ${\displaystyle i}$ частицы, добавив к ней ее (ранее обновленную) скорость ${\displaystyle i}$ умножается на ${\displaystyle c_{i}}$
• Обновите скорость ${\displaystyle i}$ частицы, добавив к ней ее ускорение (в обновленном положении), умноженное на ${\displaystyle d_{i}}$

Пример первого порядка

В симплектический метод Эйлера - интегратор первого порядка с ${\displaystyle k=1}$ и коэффициенты

${\displaystyle c_{1}=d_{1}=1.\ \ }$

Обратите внимание, что приведенный выше алгоритм не работает, если требуется обратимость по времени. Алгоритм должен быть реализован в двух частях: одна для положительных временных шагов, другая для отрицательных временных шагов.

Пример второго порядка

В Верле метод интегратор второго порядка с ${\displaystyle k=2}$ и коэффициенты

${\displaystyle c_{1}=0,\qquad c_{2}=1,\qquad d_{1}=d_{2}={\tfrac {1}{2}}.}$

С ${\displaystyle c_{1}=0}$, приведенный выше алгоритм является симметричным по времени. Алгоритм состоит из 3 шагов, и шаги 1 и 3 точно такие же, поэтому версия с положительным временем может использоваться для отрицательного времени.

Пример третьего порядка

Симплектический интегратор третьего порядка (с ${\displaystyle k=3}$) был обнаружен Рональдом Рутом в 1983 году.^[1]Одно из многих решений дает

${\displaystyle {\begin{aligned}c_{1}&=1,&c_{2}&=-{\tfrac {2}{3}},&c_{3}&={\tfrac {2}{3}},\\d_{1}&=-{\tfrac {1}{24}},&d_{2}&={\tfrac {3}{4}},&d_{3}&={\tfrac {7}{24}}.\end{aligned}}}$

Пример четвертого порядка

Интегратор четвертого порядка (с ${\displaystyle k=4}$) также был обнаружен Рут в 1983 году и в то время распространен частным образом среди разработчиков ускорителей частиц. Об этом было написано в живой обзорной статье Forest.^[2]Этот интегратор четвертого порядка был опубликован в 1990 году Форестом и Рут, а также независимо обнаружен двумя другими группами примерно в то же время.^[3][4][5]

${\displaystyle {\begin{aligned}c_{1}&=c_{4}={\frac {1}{2(2-2^{1/3})}},&c_{2}&=c_{3}={\frac {1-2^{1/3}}{2(2-2^{1/3})}},\\d_{1}&=d_{3}={\frac {1}{2-2^{1/3}}},&d_{2}&=-{\frac {2^{1/3}}{2-2^{1/3}}},\quad d_{4}=0.\end{aligned}}}$

Для определения этих коэффициентов Формула Бейкера – Кэмпбелла – Хаусдорфа может быть использован. Йошида, в частности, дает элегантный вывод коэффициентов для интеграторов высшего порядка. Позже Бланес и Стон^[6] доработанный разделенный Методы Рунге – Кутты для интегрирования систем с сепарабельными гамильтонианами с очень малыми константами ошибок.

Методы расщепления общих несепарабельных гамильтонианов

Общие несепарабельные гамильтонианы также могут быть проинтегрированы явно и симплектически.

Для этого Тао ввел ограничение, которое связывает две копии фазового пространства вместе, чтобы сделать возможным явное разделение таких систем.^[7]Идея в том, что вместо ${\displaystyle H(Q,P)}$, один моделирует ${\displaystyle {\bar {H}}(q,p,x,y)=H(q,y)+H(x,p)+\omega \left(\|q-x\|_{2}^{2}/2+\|p-y\|_{2}^{2}/2\right)}$, решение которой совпадает с решением ${\displaystyle H(Q,P)}$ в том смысле, что ${\displaystyle q(t)=x(t)=Q(t),p(t)=y(t)=P(t)}$.

Новый гамильтониан полезен для явного симплектического интегрирования, поскольку его можно разбить на сумму трех субгамильтонианов: ${\displaystyle H_{A}=H(q,y)}$, ${\displaystyle H_{B}=H(x,p)}$, и ${\displaystyle H_{C}=\omega \left(\|q-x\|_{2}^{2}/2+\|p-y\|_{2}^{2}/2\right)}$. Можно явно получить точные решения всех трех субгамильтонианов: оба ${\displaystyle H_{A},H_{B}}$ решения соответствуют сдвигам несовпадения положения и импульса, а ${\displaystyle H_{C}}$ соответствует линейному преобразованию. Чтобы симметрично моделировать систему, нужно просто составить эти карты решений.

Приложения

В физике плазмы

В последние десятилетия симплектический интегратор в физике плазмы стал активной темой исследований.^[8] потому что прямое применение стандартных симплектических методов не удовлетворяет потребности в крупномасштабном моделировании плазмы, обеспечиваемом вычислительной аппаратурой от пета-до экза-масштабов. Необходимо разработать специальные симплектические алгоритмы, использующие особые структуры исследуемой физической проблемы. Одним из таких примеров является динамика заряженных частиц в электромагнитном поле. При канонической симплектической структуре гамильтониан динамики имеет вид

$${\displaystyle H({\boldsymbol {p}},{\boldsymbol {x}})={\frac {1}{2}}\left({\boldsymbol {p}}-{\boldsymbol {A}}\right)^{2}+\phi ,}$$

чей ${\textstyle {\boldsymbol {p}}}$-зависимость и ${\textstyle {\boldsymbol {x}}}$-зависимость неотделима, и стандартные явные симплектические методы не применяются. Однако для крупномасштабного моделирования на массивно-параллельных кластерах предпочтительны явные методы. Чтобы преодолеть эту трудность, мы можем изучить конкретный способ, которым ${\textstyle {\boldsymbol {p}}}$-зависимость и ${\textstyle {\boldsymbol {x}}}$-зависимости запутаны в этом гамильтониане, и пытаются разработать симплектический алгоритм только для того или иного типа проблемы. Прежде всего отметим, что ${\textstyle {\boldsymbol {p}}}$-зависимость квадратичная, поэтому симплектический метод Эйлера первого порядка неявный в ${\textstyle {\boldsymbol {p}}}$ на самом деле явный. Это то, что используется в канонической симплектической частица в клетке (PIC) алгоритм.^[9] Чтобы построить явные методы высокого порядка, отметим, что ${\textstyle {\boldsymbol {p}}}$-зависимость и ${\textstyle {\boldsymbol {x}}}$-зависимость в этом ${\textstyle H({\boldsymbol {p}},{\boldsymbol {x}})}$ разделимы по продукту, и явные симплектические алгоритмы 2-го и 3-го порядков могут быть построены с использованием производящих функций.^[10]

Более элегантная и универсальная альтернатива - взглянуть на следующую неканоническую симплектическую структуру задачи:

$${\displaystyle i_{({\dot {\boldsymbol {x}}},{\dot {\boldsymbol {v}}})}\Omega =-dH,\ \ \ \Omega =d({\boldsymbol {v}}+{\boldsymbol {A}})\wedge d{\boldsymbol {x}},\ \ \ H={\frac {1}{2}}{\boldsymbol {v}}^{2}+\phi .}$$

Здесь ${\textstyle \Omega }$ - непостоянная неканоническая симплектическая форма. Общий симплектический интегратор для непостоянной неканонической симплектической структуры, явной или неявной, не известен. Однако для этой конкретной задачи семейство явных неканонических симплектических интеграторов высокого порядка может быть построено с использованием метода He-расщепления.^[11] Расщепление ${\textstyle H}$ на 4 части,

$${\displaystyle {\begin{aligned}H&=H_{x}+H_{y}+H_{z}+H_{\phi },\\H_{x}&={\frac {1}{2}}v_{x}^{2},\ \ H_{y}={\frac {1}{2}}v_{y}^{2},\ \ H_{z}={\frac {1}{2}}v_{z}^{2},\ \ H_{\phi }=\phi ,\end{aligned}}}$$

мы случайно обнаруживаем, что для каждой подсистемы, например,

$${\displaystyle i_{({\dot {\boldsymbol {x}}},{\dot {\boldsymbol {v}}})}\Omega =-dH_{x}}$$

и

$${\displaystyle i_{({\dot {\boldsymbol {x}}},{\dot {\boldsymbol {v}}})}\Omega =-dH_{\phi },}$$

карту решения можно записать явно и точно рассчитать. Тогда явные неканонические симплектические алгоритмы высокого порядка могут быть построены с использованием различных композиций. Позволять ${\textstyle \Theta _{x},\Theta _{y},\Theta _{z}}$ и ${\textstyle \Theta _{\phi }}$ обозначают карты точных решений для 4 подсистем. Симплектическая схема 1-го порядка - это

$${\displaystyle {\begin{aligned}\Theta _{1}\left(\Delta \tau \right)=\Theta _{x}\left(\Delta \tau \right)\Theta _{y}\left(\Delta \tau \right)\Theta _{z}\left(\Delta \tau \right)\Theta _{\phi }\left(\Delta \tau \right)~.\end{aligned}}}$$

Симметричная симплектическая схема 2-го порядка:

$${\displaystyle {\begin{aligned}\Theta _{2}\left(\Delta \tau \right)&=\Theta _{x}\left(\Delta \tau /2\right)\Theta _{y}\left(\Delta \tau /2\right)\Theta _{z}\left(\Delta \tau /2\right)\Theta _{\phi }\left(\Delta \tau \right)\\&\Theta _{z}\left(\Delta t/2\right)\Theta _{y}\left(\Delta t/2\right)\Theta _{x}\left(\Delta t/2\right)\!,\end{aligned}}}$$

которое представляет собой обычно модифицированное расщепление Стренга. А ${\textstyle 2(l+1)}$Схема -го порядка может быть построена из ${\textstyle 2l}$схема -го порядка с использованием метода тройного прыжка,

$${\displaystyle {\begin{aligned}\Theta _{2(l+1)}(\Delta \tau )&=\Theta _{2l}(\alpha _{l}\Delta \tau )\Theta _{2l}(\beta _{l}\Delta \tau )\Theta _{2l}(\alpha _{l}\Delta \tau )~,\\\alpha _{l}&=1/(2-2^{1/(2l+1)})~,\\\beta _{l}&=1-2\alpha _{l}~.\end{aligned}}}$$

Метод He-расщепления является одним из ключевых методов, используемых в сохраняющих структуру геометрических частица в клетке (PIC) алгоритмы ^[12][13][14].

Смотрите также

• Энергетический дрейф
• Мультисимплектический интегратор
• Вариационный интегратор
• Интеграция Верле

Рекомендации

1. Рут, Рональд Д. (август 1983 г.). «Техника канонической интеграции». IEEE Transactions по ядерной науке. НС-30 (4): 2669–2671. Bibcode:1983ITNS ... 30.2669R. Дои:10.1109 / TNS.1983.4332919.
2. Лес, Этьен (2006). «Геометрическая интеграция для ускорителей частиц». J. Phys. A: Математика. Gen. 39 (19): 5321–5377. Bibcode:2006JPhA ... 39.5321F. Дои:10.1088 / 0305-4470 / 39/19 / S03.
3. Forest, E .; Рут, Рональд Д. (1990). «Симплектическое интегрирование четвертого порядка» (PDF). Physica D. 43: 105–117. Bibcode:1990 ФИД ... 43..105F. Дои:10.1016 / 0167-2789 (90) 90019-Л.
4. Йошида, Х. (1990). «Построение симплектических интеграторов высших порядков». Phys. Lett. А. 150 (5–7): 262–268. Bibcode:1990ФЛА..150..262Л. Дои:10.1016/0375-9601(90)90092-3.
5. Candy, J .; Розмус, В. (1991). "Алгоритм симплектического интегрирования для отделимых гамильтоновых функций". J. Comput. Phys. 92 (1): 230–256. Bibcode:1991JCoPh..92..230C. Дои:10.1016 / 0021-9991 (91) 90299-Z.
6. Blanes, S .; Стон, П. С. (май 2002 г.). «Практические симплектические разбитые методы Рунге – Кутты и Рунге – Кутты – Нистрома». Журнал вычислительной и прикладной математики. 142 (2): 313–330. Bibcode:2002JCoAM.142..313B. Дои:10.1016 / S0377-0427 (01) 00492-7.
7. Тао, Молей (2016). «Явная симплектическая аппроксимация несепарабельных гамильтонианов: алгоритм и длительное время работы». Phys. Ред. E. 94 (4): 043303. arXiv:1609.02212. Bibcode:2016PhRvE..94d3303T. Дои:10.1103 / PhysRevE.94.043303. PMID  27841574.
8. Цинь, H .; Гуань, X. (2008). «Вариационный симплектический интегратор для движения направляющего центра заряженных частиц для длительного моделирования в общих магнитных полях» (PDF). Письма с физическими проверками. 100: 035006. Дои:10.1103 / PhysRevLett.100.035006. PMID  18232993.
9. Цинь, H .; Liu, J .; Сяо, Дж. (2016). «Канонический симплектический метод частиц в ячейках для длительного крупномасштабного моделирования уравнений Власова – Максвелла». Термоядерная реакция. 56 (1): 014001. arXiv:1503.08334. Bibcode:2016NucFu..56a4001Q. Дои:10.1088/0029-5515/56/1/014001.
10. Zhang, R .; Цинь, H .; Тан, Ю. (2016). «Явные симплектические алгоритмы на основе производящих функций для динамики заряженных частиц». Физический обзор E. 94 (1): 013205. arXiv:1604.02787. Дои:10.1103 / PhysRevE.94.013205. PMID  27575228.
11. Привет.; Цинь, H .; Солнце, Ю. (2015). «Гамильтоновы методы интегрирования для уравнений Власова-Максвелла». Физика плазмы. 22: 124503. arXiv:1505.06076. Дои:10.1063/1.4938034.
12. Xiao, J .; Цинь, H .; Лю, Дж. (2015). "Явные неканонические симплектические алгоритмы частиц в ячейках высокого порядка для систем Власова-Максвелла". Физика плазмы. 22 (11): 112504. arXiv:1510.06972. Bibcode:2015ФПЛ ... 22к2504Х. Дои:10.1063/1.4935904.
13. Краус, М; Корманн, К; Morrison, P .; Зоннендрукер, Э (2017). «GEMPIC: геометрические методы электромагнитных частиц в ячейках». Журнал физики плазмы. 83 (4): 905830401. arXiv:1609.03053. Дои:10.1017 / S002237781700040X.
14. Xiao, J .; Цинь, H .; Лю, Дж. (2018). "Сохраняющие структуру геометрические методы частиц в ячейках для систем Власова-Максвелла". Наука и технологии плазмы. 20 (11): 110501. arXiv:1804.08823. Дои:10.1088 / 2058-6272 / aac3d1.

• Леймкухлер, Бен; Райх, Себастьян (2005). Моделирование гамильтоновой динамики. Издательство Кембриджского университета. ISBN  0-521-77290-7.
• Хайрер, Эрнст; Любич, Кристиан; Ваннер, Герхард (2006). Геометрическое численное интегрирование: сохраняющие структуру алгоритмы для обыкновенных дифференциальных уравнений (2-е изд.). Springer. ISBN  978-3-540-30663-4.