chapter07.qq

\chapter Консервативные системы с одной степенью свободы \label{chap:7:cons}

\section Отступление: уравнения в Гамильтоновой форме

\definition \label def:7:Ham
    Пусть $H(x,y)$ — какая-то гладкая функция, $H(x,y):\mathbb{ R}^2 \to
    \mathbb{ R}$. Система
    \equation \label{eq:7:Ham}
        \begin{cases}
            \dot x=\frac{\partial H}{\partial y}\\\\
            \dot y=-\frac{\partial H}{\partial x}
        \end{cases}
    называется системой \em{в гамильтоновой форме}. Функция $H$ называется
    \em{функцией Гамильтона}

\proposition
    Функция $H$ является \snref[первым интегралом] системы \ref{eq:7:Ham}

\proof 
    Посчитаем \snref[производную функции $H$ вдоль векторного
    поля][snip:Lie-deriv] 
    \eq
        v=\left(\frac{\partial H}{\partial y}, -\frac{\partial H}{\partial x}\right). 
    Имеем:

    \eq
        L_v H=\frac{\partial H}{\partial x} \dot x + \frac{\partial H}{\partial
            y}\dot y=\frac{\partial H}{\partial x}\frac{\partial H}{\partial
            y}-\frac{\partial H}{\partial y}\frac{\partial H}{\partial x}=0
    
    Иными словами, скорость изменения функции $H$ при движении вдоль решений
    дифференциального уравнения \ref{eq:7:Ham} равна нулю, то есть $H$ — первый
    интеграл. (См. \ref[утверждение][prop:6:crit] из
    \ref[главы][chap:6:firstint]).

Функция Гамильтона $H$ является аналогом полной энергии — которая, как известно, сохраняется. (Слово «консервативные» в заголовке лекции означает, что что-то сохраняется.) Конечно, не всякая система, обладающая первым интегралом, является гамильтоновой. Однако, если вам дали функцию и попросили придумать систему, для которой эта функция является первым интегралом, то вы можете смело использовать гамильтонову форму.

\section Уравнение Ньютона
Рассмотрим уравнение

\equation \label eq:7:newton
    \ddot x=F(x).
Оно задает поведение частицы с единичной массой, движущейся в одномерном фазовом
пространстве под действием силы $F$, которая зависит только от положения
частицы.

\figure \showcode \collapsed
    \pythonfigure \style max-width: 300px;
        plt.figure(figsize=(3,1))
        plt.axis([-6,6,-1,2])
        plt.axis('off')
        plt.arrow( -5, 0, 10, 0, fc="k", ec="k", head_width=.5, head_length=1  )
        plt.plot([0],[0],'o',linewidth=1)
        plt.arrow( 0, 1, -1, 0, fc="k", ec="k", head_width=.5, head_length=.5)
        plt.text(-1,1.6,"$F(x)$",fontsize=16)
        plt.text(4.5, 1, "$x$", fontsize=16)
    \caption На точечное тело с координатой $x$ действует сила $F(x)$

Перепишем уравнение \ref{eq:7:newton} в виде системы:

\equation \label eq:7:newton-system
    \begin{cases}
        \dot x=y\\\\
        \dot y=F(x)
    \end{cases}
Можно найти фазовые кривые этого уравнения
\align
    \item
        \frac{dy}{dx}& = \frac{F(x)}{y}
    \item 
        y\,dy &= F(x)\,dx
    \item
        \int y\, dy &= \int F(x)\, dx
    \item \label eq:7:phasecurves
        \frac{y^2}{2} & - \int F(x)\, dx = C
        
Введём следующие функции:
\itemize
    \item $U(x)=-\int F(x)\, dx$ — потенциальная энергия.
    \item $V(y)=\frac{y^2}{2}$ — кинетическая энергия.
    \item $H(x,y)=U(x)+V(y)$ — полная энергия.

\proposition 
    Функция $H$ является первым интегралом системы \ref{eq:7:newton-system}.

\proof 
    Во-первых, это следует из \ref{eq:7:phasecurves}. Во-вторых, можно заметить,
    что система \ref{eq:7:newton-system} — гамильтонова. Наконец, можно явно
    посчитать \snref[производную вдоль векторного поля][snip:Lie-deriv]: 
    \eq
        L_{(y,F(x))} H = \frac{\partial H}{\partial x}y+\frac{\partial
            H}{\partial y} F(x)=-F(x)y+yF(x)=0
    Утверждение доказано.

Итак, полная энергия системы является первым интегралом. Это позволяет нам рисовать фазовые кривые любых уравнений вида \ref{eq:7:newton}

\example \label ex:07:oscillator
    Уже знакомое нам уравнение осциллятора.
    \eq F(x)=-x.
    Его потенциальная энергия:
    $U(x)=\frac{x^2}{2}$, полная энергия $H(x,y)=\frac{x^2}{2}+\frac{y^2}{2}$.
    Построим график полной энергии:

    \figure \showcode 
        \plotly
            X, Y = np.mgrid[-2:2:100j, -2:2:100j]
            Z = X**2/2 + Y**2/2

            data = go.Surface(x = X, y = Y, z = Z)
            plot([data])
        \caption
            График полной энергии гармонического осциллятора


И его линии уровня (они же — фазовые кривые):

\figure \showcode \collapsed
    \pythonfigure \style max-width: 400px;
        plt.figure(figsize=(6,6))
        ob.axes4x4(labels=('x','y'))
        F = lambda x, y: y**2 / 2 + x**2 / 2
        ob.mcontour(np.linspace(-5,5,300),np.linspace(-5,5,300),F, 
                 levels=[0.01,-0.5,0.5,0,-3,3,6,-6,12],linewidths=2,
                 vmin=0,vmax=25,cmap=plt.cm.coolwarm)
        ob.vectorfield(np.arange(-4,4,0.5),np.arange(-4,4,0.5), lambda x, y: (y, -x))
    \caption
        Линии уровня функции полной энергии гармонического осциллятора

Заметим, что тут есть единственная особая точка (положение равновесия), где правая часть системы \ref{eq:7:newton-system} равна нулевому вектору — это точка $(0,0)$.

\example \label ex:07:inverse-pendulum
    «Перевернутый маятник» вблизи положения равновесия.
    \eq
        F(x)=x.
    Здесь сила пропорциональна отклонению и действует в сторону увеличения
    отклонения.
        
    Потенциальная энергия $U(x)=-\frac{x^2}{2}$, полная энергия
    $H(x,y)=\frac{y^2}{2}-\frac{x^2}{2}$. Как и в предыдущем примере, построим
    график и линии уровня полной энергии.

    \figure \showcode \collapsed
        \plotly
            X, Y = np.mgrid[-2:2:100j, -2:2:100j]
            Z = -X**2/2 + Y**2/2

            data = go.Surface(x = X, y = Y, z = Z)
            plot([data])

        \pythonfigure \style max-width: 400px;
            ob.axes4x4(labels=('x','y'))
            F=lambda x, y: y**2 / 2 - x**2 / 2
            ob.mcontour(np.linspace(-5, 5, 300), np.linspace(-5, 5, 300),
                     F, levels=[-0.5, 0.5, 0, -3, 3, 6, -6],
                     linewidths=2, cmap=plt.cm.coolwarm, vmin=-15, vmax=15)
            ob.vectorfield(np.arange(-4,4,0.5),np.arange(-4,4,0.5), lambda x, y: (y, x))

        \caption
            График и линии уровня полной энергии перевернутого маятника вблизи положения
            равновесия 

Фазовыми кривыми этого уравнения является семейства гипербол, а также пять специальных кривых: четыре прямолинейных луча, выходящих из начала координат (не содержащих при этом начало координат) и само начало координат.

\section Математический маятник
    \label sec:7:pendulum

\definition
    \em{Математическим маятником} называется механическая система, состоящая из
    стержня, один конец которого закреплён, а к другому прикреплен
    маленький грузик (материальная точка). Стержень может свободно вращаться
    вокруг закрепленного конца. Трение отсутствует. Масса стержня считается
    нулевой.

\figure \label fig:07:pendulum-forces
    \rawhtml
        <img src="https://upload.wikimedia.org/wikipedia/commons/6/66/Pendulum_gravity.svg" style='max-width:300px'>
    \caption 
        Силы, действующие на маятник. Автор изображения: \href[Krishnavedala /
        Wikimedia
        Commons][https://commons.wikimedia.org/wiki/User:Krishnavedala], лицензия
        \href[CC BY-SA][https://creativecommons.org/licenses/by-sa/3.0/deed.ru].
        \href[Оригинальное
        изображение][https://commons.wikimedia.org/wiki/File:Pendulum_gravity.svg].

Пусть маятник отклонён от положения равновесия на угол $\theta$. Рассмотрим
силы, действующие на грузик (см. рис.  \ref{fig:07:pendulum-forces}).  На него
действует сила тяжести $mg$ (здесь $g$ — ускорение свободного падения вблизи
поверхности Земли), направленная вертикально вниз, и сила реакции опоры
со стороны стержня, направленная вдоль стержня. Разложим силу тяжести в сумму
двух сил: действующей в направлении стержня ($mg \cos \theta$) и перпендикулярно
ему ($mg \sin \theta$). Сила реакции опоры заставляет грузик двигаться по
окружности, не позволяя ему отдаляться от центра вращения или приближаться к
нему: это означает, что она в точности компенсирует проекцию силы тяжести
на ось стержня, равную $mg\cos \theta$. Интересующая нас проекция силы тяжести
направлена перпендикулярно к направлению стержня и равна $mg \sin \theta$: её
никто не компенсирует и именно она приводит маятник в движение.

Теперь нам нужно записать соответствующее дифференциальное уравнение с помощью
второго закона Ньютона. В качестве фазовой координаты мы рассматриваем угол
$\theta$, но положение (а значит и скорость, и ускорение) грузика зависит не
только от $\theta$, но и от длины маятника $l$. Нетрудно видеть, что скорость и
ускорение грузика пропорциональны соответственно угловой скорости и ускорению
(то есть $d\theta/dt$ и $d^2\theta/dt^2$), а также длине маятника.
Действительно, если увеличить длину вдвое, при такой же угловой скорости
реальная скорость грузика увеличится также вдвое. Это означает, что реальное
ускорение равно $l\cdot d^2 \theta/d\theta^2$ и второй закон Ньютона в проекции
на направление, перпендикулярное к стержню, запишется в следующем виде:

\eq
    ml\frac{d^2 \theta}{dt^2}=-mg\sin \theta
Знак минус в правой части связан с тем фактом, что сила направлена в направлении
уменьшения $\theta$. Сокращая на $m$ и перенося $l$ в правую часть, имеем:

\eq
    \frac{d^2 \theta}{dt^2}=\frac{g}{l}\sin \theta

Будем считать для простоты, что у нас $l=g$.

\eq
    \frac{d^2 \theta}{dt^2}=-\sin \theta

Потенциальная энергия имеет вид:
\eq
    U(\theta)=-\int \sin \theta d\theta = - \cos \theta

Полная энергия:

\eq
    E(\theta, v) = \frac{v^2}{2} - \cos \theta.
Здесь мы считаем, что $m=1$.

Построим график полной энергии и её линии уровня.

\figure \showcode \collapsed
    \plotly
        x, y = np.mgrid[-2*np.pi:2*np.pi:300j, -2:2:300j]
        surface = go.Surface(
            x=x, y=y, z=-np.cos(x)+y**2/2, 
            contours=dict(
                y=dict(show=False), 
                x=dict(show=True), 
                z=dict(show=True, project=dict(z=True)))
        )
        plot([surface])
    

Построим также график функции потенциальной энергии и линии уровня полной
энергии один под другим (они связаны).

\figure \showcode \collapsed
    \pythonfigure \style max-width: 400px;
        fig = plt.figure(figsize=(6,12))
        ax=fig.add_subplot(211) #FIGURE 1
        plt.plot([3.14,3.14],[4,-4],'--k',[-3.14,-3.14],[4,-4],'--k',[0,0],[4,-4],'--k')
        ob.axes4x4(labels=('\\\\theta', 'U'))
        ob.mplot(np.linspace(-4, 4), lambda x: -np.cos(x))
        ax=fig.add_subplot(212) #FIGURE 2
        ob.axes4x4(labels=('\\\\theta','v'))
        ob.mcontour(np.linspace(-5, 5, 300), np.linspace(-5, 5, 300),
            lambda x, y: y**2/2-np.cos(x), 
            levels=list(np.arange(-1, 5, 0.5)) + [0.1, -1 + 0.1, -1 + 0.01])
        plt.plot([3.14,3.14],[4,-4],'--k',[-3.14,-3.14],[4,-4],'--k',[0,0],[4,-4],'--k')
    \caption Потенциальная энергия и линии уровня полной энергии

\question В какую сторону должны быть направлены стрелки на картинке в фазовом
    пространстве?

\subsection Интерпретация

Положениями равновесия системы Ньютона являются точки в фазовом пространстве, в которых сила равна нулю и скорость равна нулю (то есть это точки горизонтальной оси, в которых $F(x)=0$). Они же соответствуют точкам максимума и минимума для потенциальной энергии (что неудивительно, потому что последняя с точностью до знака равна интегралу от силы). Вблизи точек минимума потенциальной энергии система ведёт себя подобно
осциллятору (см. \ref[пример][ex:07:oscillator]), а вблизи точек максимума
потенциальной энергии — подобно системе из \ref[примера][ex:07:inverse-pendulum].

Положения равновесия, соответствующие минимумам потенциальной энергии, являются
устойчивыми — точнее, \em{устойчивыми по Ляпунову}. Точное определение этого
термина мы дадим позже, но его суть состоит в том, что малое отклонение от
положения равновесия не приводит к кардинальному изменению траектории: вместо
состояния покоя мы будем наблюдать колебания вблизи положения равновесия с
небольшой амплитудой, которые остаются в некоторой окрестности этого положения
равновесия.

Положения равновесия, соответствующие максимумам потенциальной энергии, являются
неустойчивыми: сколь угодно малое отклонение от такого равновесия чаще всего приводит к
уходу траектории далеко от равновесия.

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

Допустим теперь, что мы начинаем с положения, находящегося левее положения
равновесия, соответствующего вертикально стоящему перевернутому маятнику $\theta=\pi$. Если начальная скорость не слишком велика, маятник приблизится
к положению равновесия, но не достигнет его, затормозит (в момент пересечения
оси $\theta$) и станет снова падать влево. Если же начальная скорость большая,
маятник перелетит через положение равновесия, и упадёт направо. Наконец, есть
промежуточная траектория (она называется \em{сепаратрисой}) — вдоль неё маятник
будет обречен вечно приближаться к положению равновесия, но никогда его не
достигнет (если бы он его достиг, нарушилась бы \snref[теорема существования и
единственности]).

В окрестности точки $(0,0)$ (точнее, любой точки $(2\pi k,0)$ имеем минимум потенциальной энергии и точку устойчивого равновесия (фазовые кривые близки к окружностям), в окрестности $(\pi+2\pi k,0)$ — максимум потенциальной энергии и точку неустойчивого равновесия (фазовые кривые близки к гиперболам). Фазовые кривые, похожие на окружности, соответствуют колебательным движениям маятника. Фазовые кривые, похожие на синусоиды, соответствуют вращению маятника («солнышко» на качелях).

Эта картинка имеет ещё одну интерпретацию. Представим себе, что график
потенциальной энергии сделан из тонкой проволоки, по которой может свободно
скользить (без трения) маленькая бусинка. Движение бусинки в этом случае будет
описываться фазовым портретом, похожим на портрет математического маятника.
Например, если отпустить бусинку с нулевой начальной скоростью вблизи точек
минимума потенциальной энергии бусинка будет колебаться в «ямке», а вблизи точек
максимума — падать с «горки» влево или вправо. Аналогичная интерпретация
работает и в случае других функций полной энергии.

Вы можете поэкспериментировать с математическим маятником с помощью \href[интерактивной демонстрации][http://math-info.hse.ru/f/2015-16/nes-ode/pendulum.html].