Волновое уравнение
Волновое уравнение – это пример непрерывных многочастичных систем, в которых частицы связаны пружинами. Пример конечной системы связанных пружинами масс для имитации колебаний приведен при рассмотрении колебаний винного бокала под действием звуковых волн (см. параграф 4.16). Волновое уравнение играет важную роль в электромагнетизме, где оно выводится из уравнения Максвелла для электромагнитных полей.
Волновое уравнение содержит константу – квадрат скорости распространения волны (с). В классической механике (в задаче о поперечных колебаниях непрерывной струны) квадрат скорости равен натяжению струны, деленному на линейную плотность массы.
Для решения волнового уравнения в одном измерении дернутой щипком струны применим ряд Фурье. Ранее это делалось в задаче о ряде Фурье, а в дальнейшем понадобится при анализе одномерного уравнения теплопроводности. Пусть u(x, t) обозначает вертикальное (поперечное) смещение струны как функцию координаты и времени. Предполагается, что продольного движения нет, что возможно в пределе малой амплитуды.
Из уравнения Ньютона в механическое уравнение переходят две производные по времени. Две пространственные производные возникают вследствие гармонической силы, действующей на точечную массу от двух ее соседей: разность относительных смещений можно записать как разность первых производных и, следовательно, как вторую производную по координате:
- restart;
- WE:=diff(u(x,t),t$2) = c^2*diff(u(x,t),x$2);
Берем «единичную» систему, где 0 < x < 1, c = 1.
- sol:=pdsolve(WE);
Команда pdsolve «понимает», что решения волнового уравнения можно выразить в виде суммы решений для бегущих влево и вправо волн. Хотя это интересно само по себе (на струне волны могут пересекаться, т. е. возможны решения в виде стоячих волн), это не обязательно использовать непосредственно, чтобы найти решение.
- c:=1;
- pdsolve(WE,HINT=f(x)*g(t));
Если попросить Maple найти решение в форме произведения, он напомнит о том, как разделение переменных связывает пространственную и временную части в решении. Структура ОДУ идентична, но обычно есть краевая задача в пространстве и начальная задача во времени. _c1 – константа разделения.
- ODE_x:=diff(f(x),x$2)=c1*f(x);
Это задача на собственные значения ДУ, которая определяет базис Фурье. Условие фиксирования струны при x = 0 используется только для выбора синусоидальных решений, т. е. для удаления одной константы интегрирования. Другое граничное условие выбирает собственное значение c1. Амплитуда синусоидальных решений произвольна, поскольку задача однородна (задача о собственных значениях есть всегда). Собственные значения отрицательны (для решений типа синус/косинус) и соответствуют отрицательному значению квадрата волнового числа k в решениях sin(kx). Найдем их:
- sol_x:=dsolve({ODE_x,f(0)=0});
Это напоминает о том, что тригонометрические решения могут быть выражены через комплексные экспоненты.
- assume(c1<0);
- sol_x:=dsolve({ODE_x,f(0)=0});
Потребуем, чтобы решение исчезало на правой границе при x = 1:
- solve(sin(k*1)=0,k);
Это решение неинтересно, так как дает тривиальное: f(x) = 0. Как сделать так, чтобы Maple дал интересные для нас решения?
- _EnvAllSolutions := true:
- solve(sin(k*1)=0,k);
Это известно из средней школы. Значения, кратные π, годятся для k:
- k_n:=n->n*Pi;
- fB:=n->sin(k_n(n)*x);
Постоянная разделения:
- c_n:=n->-k_n(n)^2;
Нормированы ли эти базисные функции?
- IP:=(f1,f2)->int(expand(f1*f2),x=0..1);
- IP(fB(1),fB(1));
- IP(fB(2),fB(2));
- IP(fB(5),fB(5));
Получилось доказательство того, что необходимо переопределить базисные функции:
- fBn:=n->sqrt(2)*sin(k_n(n)*x);
- IP(fBn(1),fBn(3));
- IP(fBn(3),fBn(3));
Ортонормированность базиса, похоже, удовлетворена; проверим это, взяв разные значения для n1, n2.
Следующий вопрос: как эти базисные состояния развиваются во времени? Константа разделения определяет их поведение во времени:
- ODE_t:=diff(g(t),t$2)=c1*g(t);
- dsolve(ODE_t);
Два начальных условия требуются для определения временного поведения каждой моды (обычно это начальное смещение и начальная скорость).
Для случая ущипнутой струны можно утверждать, что имеется известная начальная конфигурация струны и нулевая производная по времени для этой конфигурации в начальный момент. Это подразумевает, что нас интересуют только временные решения косинусного типа. Другим начальным условием может быть равновесная конфигурация в начальном профиле скорости, созданной, например, ударом молотка по струне в какой-то момент времени. Это вызвало бы только решения синусоидального типа. Общий случай должен содержать суперпозицию обоих типов, как записано выше.
Каждая мода отличается значением c1 и, следовательно, колеблется с различной частотой.
Фактическая суперпозиция требуемых мод получается путем проекции начальной конфигурации на базис Фурье. Предположим, что в нулевой момент времени струна смещается так, что в середине (для x = 1/2) ее высота равна а и идет как прямая линия в промежуточной области (рис. 4.17.1, 4.17.2).
- a:=1/10;
- f0:=piecewise(x<1/2,2*a*x,x>1/2,2*(a-a*x));
- plot(f0,x=0..1);

Рис. 4.17.1
Фурье-представление для конечного базиса:
- Nb:=9;
- b:=[seq(IP(f0,fBn(j)),j=1..Nb)];
- f0_ap:=add(b[j]*fBn(j),j=1..Nb);
- plot([f0,f0_ap],x=0..1,color=[red,blue]);

Рис. 4.17.2
Аппроксимация начальной конфигурации в порядке. Расширьте размер базиса и посмотрите, что улучшится!
Теперь соберем полное решение волнового уравнения (результат представлен на рис. 4.17.3, 4.17.4):
- sol_WE:=add(b[j]*fBn(j)*cos(sqrt(-c_n(j))*t),j=1..Nb);
- plot3d(sol_WE,x=0..1,t=0..5,axes=boxed);

Рис. 4.17.3
Сделаем анимацию:
- with(plots):
- animate(sol_WE,x=0..1,t=0..5,frames=50);

Рис. 4.17.4
Измените место, где ущипнули струну (берете не среднюю, а другую точку). Что произойдет с коэффициентами Фурье в сравнении с рассмотренным случаем?
Измените начальное условие так, чтобы представить случай, когда струна находится в равновесии, но ударена молотком.
Постройте график временной эволюции первых трех мод. Найдите соотношение между их частотами.
