Numerical method of the solution of a boundary problem for the oscillating membrane

Maier Robert Valerievich Doctor of Pedagogy Professor, the department of Physics and Didactics of Physics, Glazov State Pedagogical Training Institute 427628 Russia, The Udmurt Republic, Glazov, Kalinina Street 8A, unit #79 robert_maier@mail.ru Other publications by this author

Изучение дисциплин “Численные методы” и “Компьютерное моделирование” предусматривает обсуждение способов численного решения дифференциальных уравнений в частных производных, таких как уравнение теплопроводности, волновое уравнение, уравнение Пуассона [1–6]. При этом важно не только проанализировать теоретические основы используемых методов, но и предоставить студентам возможность набрать компьютерную программу и выполнить серию вычислительных экспериментов, с целью изучения поведения системы при различных параметрах, начальных условиях и внешних воздействиях [3–5]. В настоящей работе рассмотрен метод компьютерного моделирования волновых явлений в двумерной среде, предложены простые программы и проанализированы получающиеся результаты.

Моделирование колебательных и волновых процессов в двумерных средах сводятся к решению краевой задачи для колеблющейся мембраны, которая формулируется следующим образом. Имеется мембрана, которая в положении покоя занимает область O, ограниченную контуром G. На различные точки мембраны действуют вынуждающие силы, описывающиеся функцией F=F(x,y,t). Необходимо рассчитать движение мембраны, то есть найти такую функцию координат и времени, которая удовлетворяет волновому уравнению и заданным краевым условиям ^{[1, 2, 6]}:

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

Как видно из формул, чтобы найти смещения элементов пластины на временном слое t+1, необходимо знать смещение элементов на слоях t и t–1. Для решения задачи создается цикл, в котором перебираются все узлы пространственной сетки и рассчитываются их смещения из положений равновесия при t=1, 2, 3… При этом следует учесть, что некоторые точки двумерной среды закреплены и поэтому не колеблются.

В случаях, когда форма пластины, расположение источников колебаний и закрепленных точек симметричны относительно некоторой оси, получающееся решение волнового уравнения также обладает симметрией относительно той же оси. Поэтому достаточно промоделировать колебания не всей пластины, а ее половины или четверти, правильно задав граничные условия. Если края пластины свободны, то смещениям граничных элементов присваиваются смещения соседних элементов, удаленных от края пластины на величину шага по координате. Для этого используются операторы xi(i,1):=xi(i,2); xi(i,M):=xi(i,M-1); xi(1,j):=xi(2,j); xi(N,j):=xi(N-1,j);. Если края пластины или какие-то другие ее точки закреплены, то смещения соответствующих элементов приравниваются к нулю. Для моделирования распространения волны в двумерной среде используется программа 1.

Рассмотренная модель позволяет изучить образование и распространение волн в случаях, когда источник непрерывно колеблется или совершает несколько колебаний и останавливается. Во втором случае образуется несколько волн (цуг), расходящихся во все стороны от источника. Перечислим явления, которые можно промоделировать с помощью программы: 1) распространение и отражение волн от края пластины, выпуклой или вогнутой поверхности (рис. 1.1 и рис. 4); 2) изменение длины волны при переходе волны из одной среды в другую (рис. 2); 3) интерференция от двух точечных источников (рис. 1.2); 4) стоячая волна, возникающая в результате интерференции падающей и отраженной волн; 5) дифракция волн, огибание препятствий (рис. 1.3 и рис. 3).

Допустим, точечный источник двумерной волны расположен вблизи границы раздела двух сред, в которых волна распространяется с различными скоростями v_1 и v_2. Программа 1 моделирует прохождение волны из одной среды в другую. Чтобы получить границу раздела двух сред, в узлах сетки, для которых i < N/2, задают скорость волны v_1, а во всех остальных узлах –– v_2. Из рис. 2 видно, что при переходе волны в среду с большей скоростью распространения длина волны увеличивается.

На рис. 2 представлены результаты моделирования дифракции линейной и круговой волн на препятствии. Источник И моделируется одним или несколькими элементами, колеблющимися по одному и тому же закону, а препятствие П –– совокупностью неподвижных элементов. Так как длина волны сравнима с размером препятствия П, то волны их огибают и заходят в область геометрической тени (рис. 3). Аналогичным образом моделируется дифракция на отверстии.

После незначительных изменений программа 1 позволяет промоделировать отражение двумерной волны от вогнутого цилиндрического зеркала. Вблизи центра пластины создается источник гармонических волн (четыре соседних элемента совершают несколько гармонических колебаний). Вогнутое зеркало моделируется как препятствие: ее точки не должны совершать колебания. Для этого в процедуру Pr, отвечающую за задание препятствий, необходимо добавить условный оператор: If (sqr(50-i)+sqr(50-j)>2000)and(i<50) then xi[i,j]:=0;. Результат моделирования представлен на рис. 4. Видно, как образующаяся круговая волна доходит до вогнутого зеркала и, отражаясь от нее, превращается в волну с линейным фронтом.

Рассмотренный выше метод численного решения волнового уравнения позволяет промоделировать целый ряд явлений, связанных с колебаниями упругой мембраны (пластины). Они обычно изучаются на занятиях курса “Математические методы физики”, но при этом студенты знакомятся с аналитическим способом решения соответствующих задач. Предложенная компьютерная программа после незначительных изменений (программа 2) позволяет рассчитать: 1) вынужденные колебания упругой пластины; 2) свободные колебания упругой пластины произвольной формы; 3) автоколебания упругой пластины, возникающие на одной из ее собственных частот. В результате многократных отражений и интерференции волн возникает устойчивое перераспределение энергии колебаний.

Наиболее просто промоделировать колебания упругой пластины (мембраны) квадратной формы (рис. 5). Допустим, что источник гармонических колебаний (то есть элемент, на который действует вынуждающая сила) находится в центре, края мембраны закреплены или свободны. Частота колебаний источника должна быть такой, чтобы на длине стороны пластины укладывалось от одной до четырех длин волн. Через некоторое время после начала колебаний источника, когда излучаемая волна многократно отразится от края пластины, на экране монитора наблюдаются картины, подобные приведенным на рис. 6 (первый ряд). На каждом из рисунков изображена "моментальная фотография" пластины; при этом ее элементы, смещенные на примерно равные величины, окрашены в одинаковые цвета. Видно, что пластина совершает вынужденные колебания, причем смещения ее точек примерно симметричны относительно вертикальной, горизонтальной и диагональных осей симметрии.

Чтобы повысить быстродействие программы, эту задачу следует решить для четверти пластины (рис. 5), задав граничные условия соответствующим образом и внеся изменения в процедуру Draw, чтобы программа рисовала на экране все четверти пластины в правильном расположении. Получающиеся результаты (рис. 6, второй и третий ряды) отличаются строгой симметрией. Рис. 7, 8 и 9 получены тем же способом.

На рис. 7 представлены результаты расчетов вынужденных колебаний квадратной пластины в случае, когда на нее воздействуют четыре источника, расположенные в середине ее четвертей и колеблющиеся одинаковым образом.

Чтобы промоделировать автоколебания пластины на одной из собственных частот, необходимо организовать положительную обратную связь, при которой со стороны источника на конкретный элемент пластины будет оказываться воздействие, “раскачивающее” пластину. Допустим, два источника помещены в середины верхней и нижней сторон пластины; когда соответствующий элемент пластины удаляется от нас, проходя положение равновесия, на него действует сила, увеличивающая скорость. Результаты моделирования представлены на рис. 8. На рис. 9 показаны автоколебания, возникающие в случае, когда на пластину воздействуют четыре источника, расположенные в середине ее четвертей и колеблющиеся одинаковым образом.

В настоящей работе рассмотрен простой метод численного решения волнового уравнения в двумерном случае, представлены две программы в среде Free Pascal, моделирующие волновые явления и колебания упругой мембраны. Показано, что предлагаемый метод позволяет промоделировать следующие явления: 1) распространение и отражение волн; 2) изменение длины волны при ее переходе из одной среды в другую; 3) интерференция волн от нескольких источников; 4) образование стоячей волны, возникающей в результате интерференции падающей и отраженной волн; 5) огибание волной препятствий, дифракция волн; 6) вынужденные колебания упругой пластины; 7) свободные колебания упругой пластины произвольной формы; 8) автоколебания упругой пластины, возникающие на одной из ее собственных частот. Использование этих программ на занятиях по физике и компьютерному моделированию позволит на более высоком уровне объяснить физическую сущность рассматриваемых явлений, продемонстрировать возможности компьютерных моделей, основанных на численном решении дифференциальных уравнений.