Двухслойная акустическая схема

Лекция №13

Волновое уравнение

Схема “крест”

Гиперболическое волновое уравнение обрисовывает колебания струны, движение сжимаемого газа, распространение электрических волн и ряд других явлений.

Обычной одномерной задачей является задачка описания малых колебаний натянутой струны с распределенной по длине нагрузкой f(t,x):

; (1)

; (2)

. (3)

Уравнение колебания струны (1) дополняется исходными критериями (2), которые, в отличие от параболического уравнения требуют Двухслойная акустическая схема 2-ух критерий: исходного смещения относительно положения равновесия u0(x) и исходную скорость движения u1(x). Не считая того, задаются краевые условия (3), которые именуются краевыми критериями первого рода, они обрисовывают смещение концов струны относительно положений равновесия. Краевые условия могут быть и другого рода.

Составим легкую, но эффективную разностную схему для численного Двухслойная акустическая схема решения задачки (1) — (3), выбирая для простоты равномерные по t и x сетки. В качестве шаблона разностной схемы возьмем, представленный на рис.1 шаблон в форме “креста”. Аппроксимируя производные в (1) конечными разностями, получим трехслойную схему последующего вида

, (4)

с граничными критериями

. (4¢)

Рис.1. Разностная схема (4) в форме креста

По форме шаблона схему (4) именуют схемой “крест Двухслойная акустическая схема”. Исследуем схему крест.

Процедура вычисления решения смотрится последующим образом. На нулевом слое решение понятно из исходного условия

. (5)

На первом слое решение можно также вычислить по исходным данным. Простой метод получения решения на первом слое состоит в аппроксимации второго исходного условия в (2) согласно представлению:

. (6)

Аппроксимация (6) имеет 1-ый порядок точности, тогда как Двухслойная акустическая схема в схеме крест разностный оператор 2-ой производной по времени имеет 2-ой порядок точности. Потому для более четкой аппроксимации нужно учитывать очередной член разложения, т.е.

. (6¢)

Подставляя в (6¢) utt из (1), найдем

. (6¢¢)

Схема крест (4) является очевидной, так как позволяет выразить решение на последующем слое через значения yn и на Двухслойная акустическая схема 2-ух прошлых слоях. Потому, зная значения решения на нулевом (5) и первом слое (6¢¢), можно вычислить решения на всех следующих слоях. Итак, разностное решение существует и единственно.

Для исследования аппроксимации схемы крест, разложим четкое решение по формуле Тейлора с центром в узле (tm,xn), считая все появляющиеся в разложении Двухслойная акустическая схема производные непрерывными:

Используя данные разложения, просто можно отыскать невязку схемы крест:

(7)

также невязку исходного условия (6):

,

либо более четкого исходного условия (6¢¢):

(7¢)

Исходные (2) и краевые условия (3) аппроксимируются точно. В конечном итоге разностная схема (4) с исходными критериями (5), (6¢¢) имеет порядок аппроксимации , а та же разностная схема Двухслойная акустическая схема с исходным условием (6¢) имеет несколько худший порядок — .

Устойчивость схемы крест исследуем способом разделения переменных, полагая

. (8)

Подставляя представление (8) в (4), для множителя роста гармоник находим квадратное уравнение

. (9)

По аксиоме Виета произведение корней квадратного уравнения (9) равно единице, т.е. . Условие стойкости может таким макаром быть выполнено при . Это значит, что корешки уравнения (9) должны создавать комплексно Двухслойная акустическая схема сопряженную пару. Для этого нужно, чтоб дискриминант уравнения (9) был неположительным для всех гармоник, т.е.

. (10)

Для выполнения неравенства (10) нужно и довольно соблюдение условия Куранта:

ct < h. (11)

В (11) стоит серьезное неравенство, т.к., если правильно равенство Куранта ct = h для неких гармоник, то возникает слабенькая неустойчивость из-за того Двухслойная акустическая схема, что корень квадратного уравнения (9) становится кратным.

В конечном итоге разностная схема крест (4) с исходными критериями (5), (6¢) при выполнении условия Куранта (11) сходится с порядком аппроксимации . Из наших рассуждений следует, что сходимость имеет место в норме . Схема (4) обеспечивает неплохую точность расчета решения u(t,x), имеющих четвертую непрерывную производную. Схема позволяет Двухслойная акустическая схема также рассчитывать и наименее гладкие решения.

Изучим аппроксимацию и сходимость схемы крест на численном примере решения уравнения (1) с правой частью, исходными и граничными критериями вида:

(12)

Несложно проверить, что решение задачки (1), (12) имеет последующее аналитическое решение:

. (13)

На листинге_№1 приведен код программки численного решения задачки (1), (12) при помощи разностной схемы крест (4) и исходными Двухслойная акустическая схема критериями (5), (6¢¢), т.е. с порядком аппроксимации .

Листинг_№1

%Программка решения волнового уравнения (1), (12) с

%помощью схемы крест (4)

function cross

global c

%Определяем габариты области интегрирования по времени T

%и месту a, также параметр c

T=1; a=2*pi; c=1;

%Определяем наибольшее число удвоений числа узлов по

%времени и месту smax

smax=6; N=2;

%Организуем основной цикл расчетов Двухслойная акустическая схема с разными сетками

for s=1:smax

%Определяем число узлов по месту и двойное

%число узлов по времени

N=2*N; M=2*N;

%Определяем шаги по времени и месту

tau=T/(M-1); h=a/(N-1);

r=(c^2*tau^2)/h^2;

%Определяем сетки по времени и месту

t=0:tau:T; x=0:h:a;

%Используем исходные данные из (12) для

%определения Двухслойная акустическая схема численного решения на первом и

%втором слое по формуле (6'')

for n=1:N

y(1,n)=-x(n)*sin(x(n));

y(2,n)=y(1,n)+tau*c*x(n)*cos(x(n))+...

0.5*tau^2*(c^2*(-2*cos(x(n))+...

x(n)*sin(x(n)))+f(0,x(n)));

end

%Определяем левое и правое граничные условия,

%взятые из Двухслойная акустическая схема (12)

for m=1:M

y(m,1)=0; y(m,N)=a*sin(c*t(m)-a);

end

%Применяем схему крест (4) к нашей задачке

for m=2:(M-1)

for n=2:(N-1)

y(m+1,n)=2*y(m,n)-y(m-1,n)+r*(y(m,n+1)-...

2*y(m,n)+y(m,n-1))+tau^2*f(t(m),x(n));

end

end

%Оцениваем зависимость предстепенной Двухслойная акустическая схема константы

%const=||y-u||/h^2 ошибки численного решения на

%последнем шаге по времени в норме C от h

for n=1:N

z1(n)=абс(y(M,n)-u(t(M),x(n)));

end

const(s)=max(z1)/h^2;

step(s)=h;

end

mi=0;

for m=1:5:M

mi=mi+1; ti(mi)=t(m); ni Двухслойная акустическая схема=0;

for n=1:5:N

ni=ni+1; xi(ni)=x(n); z(mi,ni)=y(m,n);

end

end

%Рисуем трехмерную поверхность решения в координатах

%время-пространство

subplot(1,2,1); surf(xi,ti,z);

%Рисуем график зависимости предстепенной константы от

%шага сетки

subplot(1,2,2); loglog(step,const);

%Определяем функцию правой части

function y=f(t,x)

global c

y=2*c^2*cos(c*t-x);

%Определяем аналитическое Двухслойная акустическая схема решение

function y=u(t,x)

global c

y=x*sin(c*t-x);

Рис.2. Численное решение волнового уравнения (1), (12) при помощи
схемы крест (4), (5), (6¢¢)

На рис.2 приведен результат работы кода программки листинга_№1. На левом рисунке приведен внешний облик численного решения , m = 1,…,M. На правом рисунке приведена зависимость предстепенной константы const в Двухслойная акустическая схема представлении зависимости ошибки численного решения от шага сетки, при всем этом шаг по времени выбирался порядка шага по месту (при выполнении условия Куранта), т.е. t ~ h. Видно, что по мере уменьшения шага h величина const вправду выходит на некое неизменное значение, что подтверждает 2-ой порядок аппроксимации Двухслойная акустическая схема схемы крест (4) с исходными критериями (5), (6¢¢).

Неявная схема

В схеме крест условие Куранта (11), обеспечивающее устойчивость, может быть довольно обременительным. По этой причине построим неявную схему для задачки (1) — (3), которая будет непременно устойчивой.

На рис.3 приведен шаблон разыскиваемой схемы. Составим по шаблону рис.3 последующую разностную схему с весами:

(14)

Рис.3. Шаблон разностной схемы Двухслойная акустическая схема (14)

Веса схемы (14) веса неотрицательны при 0 £ s £ ½. При s = 0 неявная схема (14) перебегает в схему крест. Граничные и исходные условия можно определять по типу (4¢) и (5), (6¢¢) соответственно. На 3-ем и следующих слоях разностная схема представляет собой линейную систему уравнений с трехдиагональной матрицей, в какой имеет место диагональное Двухслойная акустическая схема доминирование, т.е. решение системы уравнений существует, единственно и рассчитывается способом прогонки.

Аналогично схеме крест (4) можно показать, что неявная схема (14) на решениях с 4-ыми непрерывными производными аппроксимирует задачку (1) — (3) с порядком при всех значениях параметра s.

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

(15)

Проводя рассуждения подобные тем, которые проводились при исследовании уравнения (9), можно прийти к выводу о том, что условие стойкости производится, когда корешки уравнения (15) комплексно сопряженные, что может быть только при . Из последнего неравенства следует условие стойкости разностной схемы (14):

. (16)

Из неравенства (16) следует, что при s ³ &frac Двухслойная акустическая схема14; схема (14) непременно устойчива, когда s < ¼ схема (14) условно устойчива при .

В конечном итоге, разностная схема (14) при выборе параметра веса ¼ £ s £ ½, непременно сходится с точностью .

Решим численно задачку (1), (12), (13) при помощи неявной разностной схемы (14). На листинге_№2 приведен код соответственной программки.

Листинг_№2

%Программка решения волнового уравнения (1), (12), (13)

%при помощи неявной схемы (14)

function implicit

global Двухслойная акустическая схема c

%Определяем габариты области интегрирования по времени T

%и месту a, также параметр c

T=10; a=2*pi; c=1;

%Определяем наибольшее число удвоений числа узлов по

%времени и месту smax

smax=6; N=4;

%Определяем весовую константу

sigma=0.25;

%Организуем основной цикл расчетов с разными сетками

for s=1:smax

%Определяем число узлов по месту и равное

%число узлов по времени

N Двухслойная акустическая схема=2*N; M=N;

%Определяем шаги по времени и месту

tau=T/(M-1); h=a/(N-1);

r=(sigma*c^2*tau^2)/h^2;

r2=((1-2*sigma)*c^2*tau^2)/h^2;

%Определяем сетки по времени и месту

t=0:tau:T; x=0:h:a;

%Используем исходные данные из (12) для

%определения численного решения на первом и

%втором слое по формуле (6'')

for Двухслойная акустическая схема n=1:N

y(1,n)=-x(n)*sin(x(n));

y(2,n)=y(1,n)+tau*c*x(n)*cos(x(n))+...

0.5*tau^2*(c^2*(-2*cos(x(n))+...

x(n)*sin(x(n)))+f(0,x(n)));

end

%Определяем левое и правое граничные условия,

%взятые из (12)

for m=1:M

y(m,1)=0; y(m,N)=a*sin(c Двухслойная акустическая схема*t(m)-a);

end

%Применяем неявную схему (14) к нашей задачке

for m=2:(M-1)

alpha(2)=0; beta(2)=y(m+1,1);

for n=2:(N-1)

alpha(n+1)=r/(1+r*(2-alpha(n)));

beta(n+1)=(r2*y(m,n-1)+(2-2*r2)*y(m,n)+...

r2*y(m,n+1)+r*y(m-1,n-1)-(1+2*r)*...

y(m-1,n)+r*y(m Двухслойная акустическая схема-1,n+1)+tau^2*...

f(t(m),x(n))+r*beta(n))/...

(1+r*(2-alpha(n)));

end

for n=N:-1:2

y(m+1,n-1)=alpha(n)*y(m+1,n)+beta(n);

end

end

%Оцениваем зависимость предстепенной константы

%const=||y-u||/h^2 ошибки численного решения на

%последнем шаге по времени в норме C от h

for n=1:N

z Двухслойная акустическая схема1(n)=абс(y(M,n)-u(t(M),x(n)));

end

const(s)=max(z1)/h^2;

step(s)=h;

end

mi=0;

for m=1:10:M

mi=mi+1; ti(mi)=t(m); ni=0;

for n=1:10:N

ni=ni+1; xi(ni)=x(n); z(mi,ni)=y(m,n);

end

end

%Рисуем трехмерную поверхность решения в координатах

%время-пространство

subplot(1,2,1); surf(xi Двухслойная акустическая схема,ti,z);

%Рисуем график зависимости предстепенной константы от

%шага сетки

subplot(1,2,2); loglog(step,const);

%Определяем функцию правой части

function y=f(t,x)

global c

y=2*c^2*cos(c*t-x);

%Определяем аналитическое решение

function y=u(t,x)

global c

y=x*sin(c*t-x);

На рис.4 приведен результат работы кода программки листинга Двухслойная акустическая схема_№2. На левом рисунке приведен внешний облик численного решения , m = 1,…,M. Обратим внимание на то, что внедрение неявной непременно устойчивой схемы (14) позволило при том же примерно количестве узлов прирастить отрезок интегрирования по времени в 10 раз по сопоставлению со схемой крест, которая является условно устойчивой. На правом рисунке приведена зависимость Двухслойная акустическая схема предстепенной константы const в представлении зависимости ошибки численного решения от шага сетки, при всем этом шаг по времени выбирался порядка шага по месту, т.е. t ~ h. Видно, что по мере уменьшения шага h величина const вправду выходит на некое неизменное значение, что подтверждает 2-ой порядок аппроксимации неявной схемы (14) с исходными Двухслойная акустическая схема критериями (5), (6¢¢), (12).

Рис.4. Решение волнового уравнения (1), (12), (13) при помощи неявной
разностной схемы (14)

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

, (17)

где k(t,x) º c2(t,x) > 0. Будем считать, что функции k(t,x), f(t,x) в Двухслойная акустическая схема (17) кусочно-непрерывны, при этом разрывы недвижны, т.е. лежат на линиях x = const. Подразумевается, что на разрывах производится условие сопряжения [u] = [kux] = 0, т.е. решение u и поток kux являются непрерывными функциями.

Выберем по времени равномерную сетку, а по месту специальную неравномерную сетку, у которой все Двухслойная акустическая схема точки разрыва коэффициентов являются узлами. Построим аналог лучшей ограниченной схемы (11.19) — (11.21¢), разобранной в лекции №11:

, (18)

где

(18¢)

Известно[1], что при изготовленных догадках и при довольно гладких исходных и граничных критериях схема (18), (18¢) умеренно сходится со скоростью при выполнении условия стойкости (16).

Проведем численный расчет по разностной схеме (18), (18¢) задачки (17), (2), (3) в прямоугольной области G Двухслойная акустическая схема(t,x) = [0,T]´[0,a] при условии, что

(19)

(20)

где k1, k2, k3, f1, f2, f3 — некие неизменные величины. С учетом положений разрывов в (19), (20) определим равномерную сетку по месту xn = nh, n = 1,…,N, где N = 4l + 1, l = 1,2,… Определим также равномерную сетку по времени, т.е. tm = t m, 1,…,M. С учетом изготовленных догадок Двухслойная акустическая схема разглядим последующую схему аппроксимации для интегралов в (18¢):

. (21)

В качестве исходных и граничных критерий положим

. (22)

На листинге_№3 приведен код программки численного решения задачки (17), (19), (20), (22) при помощи лучшей разностной схемы (18), (18¢), (21).

Листинг_№3

%Программка решения волнового уравнения (17), (19),

%(20) при помощи разностной схемы (18), (18'), (21)

function best_scheme

global a k1 k2 k3 f1 f2 f Двухслойная акустическая схема3

%Определяем габариты области интегрирования по

%времени T и месту a

T=5; a=1;

%Определяем константы k1,k2,k3,f1,f2,f3

k1=1; k2=10; k3=1; f1=1; f2=10; f3=1;

%Определяем весовую константу

sigma=0.25;

%Определяем число узлов по месту и времени

M=201; N=161;

%Определяем шаги по времени и месту

tau=T/(M-1); h=a/(N-1);

r=(sigma*tau Двухслойная акустическая схема^2)/h^2;

r2=((1-2*sigma)*tau^2)/h^2;

%Определяем сетки по времени и месту

t=0:tau:T; x=0:h:a;

%Используем исходные данные из (22) для

%определения численного решения на первом и

%втором слое по формуле (6'')

for n=1:N

y(1,n)=0;

y(2,n)=y(1,n)+0.5*tau^2*f(x(n));

end

%Определяем левое и правое граничные условия,

%взятые из (22)

for m Двухслойная акустическая схема=1:M

y(m,1)=0; y(m,N)=0;

end

%Применяем неявную лучшую схему (18), (18')

%к нашей задачке

for m=2:(M-1)

%Определяем коэффициенты прогонки

alpha(2)=0; beta(2)=y(m+1,1);

for n=1:(N-1)

kapa(n)=0.5*(k(x(n))+k(x(n+1)));

end

for n=2:(N-1)

alpha(n+1)=(r*kapa(n))/(1+r*(kapa(n)+...

kapa(n-1)*(1-alpha(n))));

beta(n Двухслойная акустическая схема+1)=(r2*kapa(n-1)*y(m,n-1)+...

(2-r2*(kapa(n-1)+kapa(n)))*y(m,n)+...

r2*kapa(n)*y(m,n+1)+r*kapa(n-1)*...

y(m-1,n-1)-(1+r*(kapa(n-1)+kapa(n)))*...

y(m-1,n)+r*kapa(n)*y(m-1,n+1)+0.25*...

tau^2*(f(x(n-1))+2*f(x(n Двухслойная акустическая схема))+f(x(n+1)))+...

r*kapa(n-1)*beta(n))/(1+r*(kapa(n)+...

kapa(n-1)*(1-alpha(n))));

end

for n=N:-1:2

y(m+1,n-1)=alpha(n)*y(m+1,n)+beta(n);

end

end

%Рисуем трехмерную поверхность решения в координатах

%время-пространство

surf(x,t,y);

%Определяем функцию квадрата скорости звука k(t,x)

function y=k(x)

global Двухслойная акустическая схема a k1 k2 k3

if (x>=0)&(x<=0.25*a)

y=k1;

end

if (x>0.25*a)&(x<0.75*a)

y=k2;

end

if (x>=0.75*a)&(x<=a)

y=k3;

end

%Определяем функцию правой части

function y=f(x)

global a f1 f2 f3

if (x>=0)&(x<=0.25*a)

y=f1;

end

if (x>0.25*a)&(x<0.75*a)

y=f2;

end

if (x>=0.75*a)&(x<=a)

y=f Двухслойная акустическая схема3;

end

Рис.5. Численное решение волнового уравнения (17), (19), (20), (22) с
помощью лучшей разностной схемы (18), (18¢), (21)

Результат работы кода программки листинга_№3 приведен на рис.5. На графике представлено решение в виде волны, что, конечно, типично для волнового уравнения (17). Плоскости разрывов x = a/4 и x = 3a/4 просматриваются на поверхности решения.

Двухслойная акустическая схема

Волновое уравнение второго порядка может Двухслойная акустическая схема быть представлено в виде пары уравнений первого порядка. Разглядим систему уравнений

(23)

Покажем, что система (23), состоящая из пары уравнений первого порядка, сводится к волновому уравнению (1) второго порядка. Продифференцируем по x 2-ое уравнение в (23) и подставим vx из первого уравнения, тогда получим

.

Полное совпадение с волновым уравнением (1) наступит при условии, что vx = ut Двухслойная акустическая схема и Fx = f. Если проинтегрировать последние два уравнения по x, то

. (24)

Величины (24) именуются потенциалами скорости и правой части.

Исходные данные (2) с учетом (24) можно переписать в виде:

, (25)

где V — некая константа.

Граничные условия (3) остаются без конфигурации, т.е.

. (26)

В ряде всевозможных случаев задачка акустики (23) — (26) оказывается более комфортной для Двухслойная акустическая схема численного решения, чем волновое уравнение (1). Построим и исследуем неявную разностную схему для решения уравнений акустики.

Рис.6. Шаблон разностной схемы (27), (27¢)

Возьмем в общем случае неравномерную сетку по месту 0 = x0 < x1 < … < xN = a, в узлах которой определим величины , n = 0,1,…,N, а в серединах интервалов — величины , n = 0,1,…,N-1. Выберем шаблон разностной Двухслойная акустическая схема схемы, представленный на рис.6 и составим по нему разностную схему с весами:

, (27)

, (27¢)

где , и считается, что s1, s2 Î [0,1]. Шаблон схемы (27) выделен на рис.6 красноватым цветом (сплошные полосы), шаблон схемы (27¢) — зеленоватым цветом (точечные полосы).

Если положить , то схема (27) будет симметричной по переменной x и иметь порядок аппроксимации . Если Двухслойная акустическая схема положить s1 = s2 = ½ и , то порядок аппроксимации по времени схемы (27), (27¢) станет вторым, т.е. .

Как обычно устойчивость схемы (27), (27¢) изучим при помощи способа разделения переменных, представляя возмущения функций y и z в виде гармоник последующего вида:

. (28)

Гармоники для y и z в (28) взяты с одной и той же частотой и Двухслойная акустическая схема множителем роста, но с разными амплитудами. Подставляя (28) в (27), (27¢) и полагая , для амплитуд a и b получим однородную систему уравнений

(29)

Система линейных уравнений (29) имеет нетривиальное решение, когда ее определитель равен нулю. Это условие дает квадратное уравнение для определения множителя роста r:

, (30)

где

(31)

Любой из 2-ух корней уравнения (30) меньше по модулю единицы Двухслойная акустическая схема, если и только если

. (32)

1-ое из неравенств в (32) следует из аксиомы Виета , 2-ое можно обосновать при помощи легких, но несколько массивных выкладок. Неравенства (32) производятся для всех гармоник, когда

, (33)

неравенства (33) являются критериями стойкости для разностной схемы (27), (27¢).

Из неравенства (33) следует, что если s1 ³ ½ и s2 ³ &frac Двухслойная акустическая схема12;, то схема (27), (27¢) является непременно устойчивой. Если s1 + s2 ³ 1, но один из весов меньше ½, то схема условно устойчива . При s1 + s2 < 1 схема непременно неустойчива.

Особо выделим симметричную схему, которая реализуется при s1 = s2 = ½. В данном случае разностная схема (27), (27¢) является непременно устойчивой и сходится со скоростью . Данная схема является Двухслойная акустическая схема одной из наилучших схем для задач акустики.

Разностная схема (27), (27¢) при случайных значениях весов s1 и s2 решается методом определения из уравнения (27), т.е.

(34)

и подстановки и в (27¢), что приводит к уравнению

(34¢)

Согласно (34¢) для получения значений можно применить способ прогонки. Соответственная трехдиагональная матрица обладает свойством диагонального Двухслойная акустическая схема доминирования, что обеспечивает существование и единственность решения задачки акустики при помощи разностных схем (27), (27¢) либо, что тоже, уравнений (34), (34¢).

Численно изучим разностную схему (34), (34¢) на примере решения задачки акустики (24) при

, (35)

также при исходных и граничных значениях вида:

(36)

Как несложно проверить уравнения акустики (23) при условии (35), (36) имеют последующее аналитическое решение:

. (37)

Аналитические решения Двухслойная акустическая схема (37) получены с помощь преобразования (24) при v(t,0) = c×cos(ct), исходя из решений задачки (1), (12), (13).

На листинге_№4 приведен код программки численного решения задачки акустики (23), (35), (36) при помощи разностной схемы (34), (34¢).

Листинг_№4

%Программка численного решения задачки акустики

%(23) - (26) при помощи разностной схемы (34), (34')

function acoustics

global c

%Определяем габариты области интегрирования

%G=[0,T]x[0,a] и скорость Двухслойная акустическая схема звука c

T=8; a=2*pi; c=1;

%Избираем для анализа симметричную схему, у которой

%весовые коэффициенты

sigma1=0.5; sigma2=0.5;

%Количество удвоений smax числа узлов сетки

smax=6; N=4;

%Организуем цикл расчетов по разностной схеме

%(34), (34') с разным числом узлов по

%месту и времени

for s=1:smax

%Удваиваем число узлов по месту

N=2*N;

%Считаем, что число узлов по времени равно

%числу узлов по Двухслойная акустическая схема месту

M=N;

%Определяем шаги по времени и месту

tau=T/(M-1); h=a/(N-1);

%Определяем сетки по времени и месту

t=0:tau:T; x=0:h:a;

r=(sigma1*sigma2*c^2*tau^2)/h^2;

r2=((1-sigma1)*sigma2*c^2*tau^2)/h^2;

%Определение исходного рассредотачивания для u-y

%согласно (36)

for n=1:N

y(1,n)=u(t(1),x Двухслойная акустическая схема(n));

end

%Определение исходного рассредотачивания для v-z

%согласно (36)

for n=1:(N-1)

z(1,n)=v(t(1),x(n)+0.5*h);

end

%Определение граничных критерий для u-y

%согласно (36)

for m=1:M

y(m,1)=u(t(m),x(1));

y(m,N)=u(t(m),x(N));

end

%Организуем цикл расчетов по времени

for m=1:(M-1)

%Находим коэффициенты прогонки Двухслойная акустическая схема для вычисления u-y

alpha(2)=0; beta(2)=y(m+1,1);

for n=2:(N-1)

alpha(n+1)=r/(1+r*(2-alpha(n)));

beta(n+1)=(r2*y(m,n-1)+(1-2*r2)*y(m,n)+...

r2*y(m,n+1)+(tau/h)*(z(m,n)-z(m,n-1))+...

((sigma2*tau^2)/h)*(F(t(m)+0.5*tau,...

x(n)+0.5*h)-F(t Двухслойная акустическая схема(m)+0.5*tau,x(n-1)+...

0.5*h))+r*beta(n))/(1+r*(2-alpha(n)));

end

for n=N:-1:2

y(m+1,n-1)=alpha(n)*y(m+1,n)+beta(n);

end

%Вычисляем функцию v-z

for n=1:(N-1)

z(m+1,n)=z(m,n)+((c^2*tau)/h)*(sigma1*...

(y(m+1,n+1)-y(m+1,n))+(1-sigma1)*...

(y Двухслойная акустическая схема(m,n+1)-y(m,n)))+tau*F(t(m)+...

0.5*tau,x(n)+0.5*h);

end

end

%Находим ошибку численного решения в норме C

for m=1:M

for n=1:N

ery(m,n)=абс(y(m,n)-u(t(m),x(n)));

end

end

for m=1:M

for n=1:(N-1)

erz(m,n)=абс(z(m,n)-v(t Двухслойная акустическая схема(m),x(n)+0.5*h));

end

end

%Делим ошибку численного решения на h^2 и находим

%предстепенную константу в зависимости ошибки

%численного решения от шага сетки

const(s)=max(max(max(ery)),max(max(erz)))/h^2;

step(s)=h;

end

mi=0;

for m=1:10:M

mi=mi+1; ni=0;

for n=1:10:N

ni=ni+1; yi(mi,ni)=y(m,n Двухслойная акустическая схема);

end

end

%Рисуем численное решение функции u-y

subplot(1,3,1); surf(x([1:10:N]),t([1:10:M]),yi);

mi=0;

for m=1:10:M

mi=mi+1; ni=0;

for n=1:10:(N-1)

ni=ni+1; zi(mi,ni)=z(m,n);

end

end

%Рисуем численное решение функции v-z

subplot(1,3,2); surf(x([1:10:(N-1)]),t([1:10:M]),zi);

%Рисуем зависимость предстепенной константы

%const = max/h^2 от шага сетки Двухслойная акустическая схема h

subplot(1,3,3); semilogx(step,const);

%Определяем функцию правой части

function y=F(t,x)

global c

y=-2*c^2*sin(c*t-x);

%Определяем аналитическое решение для u

function y=u(t,x)

global c

y=x*sin(c*t-x);

%Определяем аналитическое решение для v

function y=v(t,x)

global c

y=-c*x*sin(c*t Двухслойная акустическая схема-x)+c*cos(c*t-x);

На рис.7 приведен результат работы кода программки листинга_№4. На левом рисунке изображено численное решение u = u(tm,xn) » y(m,n), на среднем рисунке — v = v(tm,xn) » z(m,n). На правом рисунке изображена динамика предстепенной константы const в Двухслойная акустическая схема оценке ошибки численного решения симметричной схемы (s1 = s2 = ½) от шага по месту h, при всем этом полагалось, что t ~ h. Видно, что при уменьшении шага сетки h предстепенная константа вправду выходит на некое неизменное значение.

Рис.7. Численное решение задачки акустики (24), (35), (36) при помощи
разностной схемы (34), (34¢)

Инварианты. Перепишем уравнения акустики при Двухслойная акустическая схема помощи инвариантов:

. (38)

Умножим 1-ое уравнение в системе уравнений акустики (23) на c и сначала сложим оба уравнения вместе, а потом вычтем из второго уравнения 1-ое, тогда получим

(39)

С учетом (25) несложно отыскать исходные данные для инвариантов:

. (40)

С учетом (26) определим краевые условия для инвариантов:

. (41)

Инвариант r удовлетворяет уравнению переноса на право (c > 0), инвариант s Двухслойная акустическая схема — уравнению переноса на лево. Для однородной задачки, когда F = 0, m1 = m2 =0, величины r, s переносятся по своим чертам без конфигурации, с чем и связано их наименование.

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

. (42)

Схема (42) является схемой бегущего счета. Несложно показать, что при выполнении условия Куранта ct £ h схема устойчива и сходится с порядком точности на два раза безпрерывно дифференцируемых функциях.

Изучим разностную схему для инвариантов на примере численного решения задачки (39) при Двухслойная акустическая схема помощи разностной схемы (42), выбирая правую часть, исходные и граничные условия согласно (35), (36). В данном случае при учете (40), (41), находим последующие выражения для правой части, исходных и граничных критерий:

(43)

Несложно проверить, что аналитическим решением задачки (39), (43) являются выражения вида:

. (44)

На листинге_№5 приведен код программки расчета задачки акустики в инвариантах (39) со специальной правой частью, исходными Двухслойная акустическая схема и граничными критериями (43).

Листинг_№5

%Программка решения уравнений акустики (39) в

%инвариантах с правой частью, исходными и

%граничными критериями вида (43)

function invariants

global c

%Определяем габариты области интегрирования по

%времени T и месту a, также определяем

%величину скорости c

T=2*pi; a=2*pi; c=1;

%Определяем число удвоений числа узлов сетки по

%времени и месту dmax

dmax=5; N=4;

%Организуем цикл расчетов Двухслойная акустическая схема на разных сетках

for d=1:dmax

%Определяем число узлов сетки по месту

%и времени

N=2*N; M=N;

%определяем шаги по времени и месту

tau=T/(M-1); h=a/(N-1); kur=(c*tau)/h;

%Определяем сетки по времени и месту

t=0:tau:T; x=0:h:a;

%Определяем исходные условия (43) для

%инвариантов r и s

for n=1:N

r Двухслойная акустическая схема(1,n)=2*c*x(n)*sin(x(n))+c*cos(x(n));

s(1,n)=c*cos(x(n));

end

%Организуем цикл расчетов по времени

for m=1:(M-1)

%Осуществляем бегущий расчет инварианта r,

%осуществляемый слева вправо

for n=2:N

r(m+1,n)=(1-kur)*r(m,n)+kur*r(m,n-1)+...

tau*F(t(m),x Двухслойная акустическая схема(n));

end

%Учитываем правое граничное условие (43)

s(m+1,N)=r(m+1,N)+2*c*a*sin(c*t(m+1)-a);

%Осуществляем бегущий расчет инварианта s,

%осуществляемый справа влево

for n=(N-1):-1:1

s(m+1,n)=(1-kur)*s(m,n)+kur*s(m,n+1)+...

tau*F(t(m),x(n));

end

%Учитываем левое граничное условие (43)

r(m Двухслойная акустическая схема+1,1)=s(m+1,1);

end

%Ошибку численного расчета инвариантов r и s

%оцениваем как разность численного и

%аналитического решений в норме C. Делим

%полученную ошибку на шаг по месту h и

%находим предстепенную константу скорости

%сходимости схемы с инвариантами

for m=1:M

for n=1:N

er_r(m,n)=абс(r(m,n)-ra(t(m),x(n Двухслойная акустическая схема)));

er_s(m,n)=абс(s(m,n)-sa(t(m),x(n)));

end

end

const(d)=max(max(max(er_r)),max(max(er_s)))/h;

step(d)=h;

end

mi=0;

for m=1:5:M

mi=mi+1; ni=0;

for n=1:5:N

ni=ni+1;

ri(mi,ni)=r(m,n);

si(mi,ni)=s(m,n Двухслойная акустическая схема);

end

end

%Рисуем численное решение инварианта r

subplot(1,3,1); surf(x([1:5:N]),t([1:5:M]),ri);

%Рисуем численное решение инварианта s

subplot(1,3,2); surf(x([1:5:N]),t([1:5:M]),si);

%Рисуем график зависимости предстепенной

%константы от шага сетки h

subplot(1,3,3); semilogx(step,const);

%Определяем функцию правой части

function y=F(t,x)

global c

y=-2*c^2*sin(c*t Двухслойная акустическая схема-x);

%Определяем аналитический вид инварианта r

function y=ra(t,x)

global c

y=-2*c*x*sin(c*t-x)+c*cos(c*t-x);

%Определяем аналитический вид инварианта s

function y=sa(t,x)

global c

y=c*cos(c*t-x);

Рис.8. Решение задачки акустики (39), (43) в инвариантах

На рис.8 приведен результат работы Двухслойная акустическая схема кода программки листинга_№5. На левом рисунке приведено изображение численного решения инварианта r, на среднем рисунке — инварианта s. Правый график показывает зависимость const от h в представлении для ошибки численного решения последующего вида:

,

где r(t,x), s(t,x) — аналитические решения (44). Из графика видно, что по мере уменьшения шага сетки h величина Двухслойная акустическая схема const вправду выходит на некое неизменное значение, при всем этом шаг по времени выбирался порядка шага по месту, т.е. t ~ h. Это подтверждает также то, что скорость сходимости разностной схемы (42) — .

Многомерные схемы

Разглядим волновое уравнение в p-мерном пространстве последующего вида:

. (45)

Решение задачки (45) подразумевает определение исходных и краевых Двухслойная акустическая схема критерий:

(46)

Схема “крест” строится аналогично одномерной схеме (4) на шаблоне, вид которого для варианта 2-ух измерений показан на рис.9. При случайном числе измерений разностная схема крест может быть записана в виде:

, (47)

при всем этом в случае переменных коэффициентов ka операторы La определяются также как и для лучшей схемы Двухслойная акустическая схема (18), (18¢).

Рис.9. Шаблон схемы крест в двумерном случае

Трехслойная схема (47) является очевидной. Вычисления при помощи схемы (47) идиентично ординарны для хоть какого числа измерений. Просто проверить, что схема имеет порядок аппроксимации .

Исследуем устойчивость схемы (47) способом разделения переменных, считая коэффициенты ka неизменными. Для этого подставим в (47) многомерную гармонику, имеющую последующий вид на 3-х Двухслойная акустическая схема временных слоях:

.

В конечном итоге для множителя роста r получим квадратное уравнение

. (48)

Анализ корней уравнения (48) указывает, что разностная схема (47) устойчива при условии, что

. (49)

Условие стойкости (49) является аналогом условия Куранта (11).


dzhasharovich-pedagogicheskij-potencial-islama-v-svetskih-obrazovatelnih-praktikah-13-00-01-obshaya-pedagogika-istoriya-pedagogiki-i-obrazovaniya-13-00-08-teoriya-i-metodika-professionalnogo-obrazovaniya.html
dzhastin-sdelal-odin-vdoh-zatem-drugoj-tretij-chetvertij-i-nakonec-sedmoj-prezhde-chem-nachal.html
dzhataka-o-krokodile-muzej-spyashih-hippi.html