Решите систему УЧП

В этом примере показано, как сформулировать, вычислите и постройте решение системы двух дифференциальных уравнений с частными производными.

Рассмотрите систему УЧП

u1t=0.0242u1x2-F(u1-u2),

u2t=0.1702u2x2+F(u1-u2).

(Функция F(y)=e5.73y-e-11.46y используется в качестве сокращения.)

Уравнение держится интервал 0x1 в течение многих времен t0. Начальные условия

u1(x,0)=1,

u2(x,0)=0.

Граничные условия

xu1(0,t)=0,u2(0,t)=0,xu2(1,t)=0,u1(1,t)=1.

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

Уравнение кода

Прежде чем можно будет закодировать уравнение, необходимо убедиться что именно в форме pdepe решатель ожидает:

c(x,t,u,ux)ut=x-mx(xmf(x,t,u,ux))+s(x,t,u,ux).

В этой форме коэффициенты УЧП с матричным знаком, и уравнение становится

[1001]t[u1u2]=x[0.024u1x0.170u2x]+[-F(u1-u2)F(u1-u2)].

Таким образом, значения коэффициентов в уравнении

m=0

c(x,t,u,ux)=[11] (только диагональные значения)

f(x,t,u,ux)=[0.024u1x0.170u2x]

s(x,t,u,ux)=[-F(u1-u2)F(u1-u2)]

Теперь можно создать функцию, чтобы закодировать уравнение. Функция должна иметь подпись [c,f,s] = pdefun(x,t,u,dudx):

  • x независимая пространственная переменная.

  • t независимая переменная времени.

  • u зависимая переменная, дифференцируемая относительно x и t. Это - двухэлементный вектор где u(1) u1(x,t) и u(2) u2(x,t).

  • dudx частичная пространственная производная u/x. Это - двухэлементный вектор где dudx(1) u1/x и dudx(2) u2/x.

  • Выходные параметры cF, и s соответствуйте коэффициентам в стандартной форме уравнения УЧП, ожидаемой pdepe.

В результате уравнения в этом примере могут быть представлены функцией:

function [c,f,s] = pdefun(x,t,u,dudx)
c = [1; 1];
f = [0.024; 0.17] .* dudx;
y = u(1) - u(2);
F = exp(5.73*y)-exp(-11.47*y);
s = [-F; F];
end

(Примечание: Все функции включены как локальные функции в конце примера.)

Начальные условия кода

Затем запишите функцию, которая возвращает начальное условие. Начальное условие применяется в первой временной стоимости и вводит значение u(x,t0) для любого значения x. Количество начальных условий должно равняться количеству уравнений, таким образом, для этой проблемы существует два начальных условия. Используйте функциональную подпись u0 = pdeic(x) записать функцию.

Начальные условия

u1(x,0)=1,

u2(x,0)=0.

Соответствующая функция

function u0 = pdeic(x)
u0 = [1; 0];
end

Граничные условия кода

Теперь запишите функцию, которая оценивает граничные условия

xu1(0,t)=0,u2(0,t)=0,xu2(1,t)=0,u1(1,t)=1.

Для проблем, созданных на интервале axb, граничные условия применяются ко всем t и также x=a или x=b. Стандартная форма для граничных условий, ожидаемых решателем,

p(x,t,u)+q(x,t)f(x,t,u,ux)=0.

Написанный в этой форме, граничных условиях для частных производных u должен быть выражен в терминах потока f(x,t,u,ux). Таким образом, граничные условия для этой проблемы

Для x=0, уравнение

[0u2]+[10][0.024u1x0.170u2x]=0.

Коэффициенты:

pL(x,t,u)=[0u2],

qL(x,t)=[10].

Аналогично, для x=1 уравнение

[u1-10]+[01][0.024u1x0.170u2x]=0.

Коэффициенты:

pR(x,t,u)=[u1-10],

qR(x,t)=[01].

Граничная функция должна использовать функциональную подпись [pl,ql,pr,qr] = pdebc(xl,ul,xr,ur,t):

  • Входные параметры xl и ul соответствовать u и x для левого контура.

  • Входные параметры xr и ur соответствовать u и x для правильного контура.

  • t независимая переменная времени.

  • Выходные параметры pl и ql соответствовать pL(x,t,u) и qL(x,t) для левого контура (x=0 для этой проблемы).

  • Выходные параметры pr и qr соответствовать pR(x,t,u) и qR(x,t) для правильного контура (x=1 для этой проблемы).

Граничные условия в этом примере представлены функцией:

function [pl,ql,pr,qr] = pdebc(xl,ul,xr,ur,t)
pl = [0; ul(2)];
ql = [1; 0];
pr = [ur(1)-1; 0];
qr = [0; 1];
end

Выберите Solution Mesh

Решение этой проблемы изменяется быстро когда t мал. Несмотря на то, что pdepe выбирает временной шаг, который является соответствующим, чтобы разрешить резкие изменения, видеть, что поведение в выходе строит вас, должен выбрать соответствующие выходные времена. Для пространственной mesh существуют пограничные слои в решении в обоих концах 0x1, таким образом, необходимо задать точки mesh там, чтобы разрешить резкие изменения.

x = [0 0.005 0.01 0.05 0.1 0.2 0.5 0.7 0.9 0.95 0.99 0.995 1];
t = [0 0.005 0.01 0.05 0.1 0.5 1 1.5 2];

Решите уравнение

Наконец, решите уравнение с помощью симметрии m, уравнение УЧП, начальные условия, граничные условия и сетки для x и t.

m = 0;
sol = pdepe(m,@pdefun,@pdeic,@pdebc,x,t);

pdepe возвращает решение в трехмерном массиве sol, где sol(i,j,k) аппроксимирует kкомпонент th решения uk оцененный в t(i) и x(j). Извлеките каждый компонент решения в отдельную переменную.

u1 = sol(:,:,1);
u2 = sol(:,:,2);

Постройте решение

Создайте объемные поверхностные диаграммы решений для u1 и u2 построенный в выбранной mesh указывает для x и t.

surf(x,t,u1)
title('u_1(x,t)')
xlabel('Distance x')
ylabel('Time t')

surf(x,t,u2)
title('u_2(x,t)')
xlabel('Distance x')
ylabel('Time t')

Локальные функции

Перечисленный здесь локальные функции помощника что решатель УЧП pdepe вызовы, чтобы вычислить решение. В качестве альтернативы можно сохранить эти функции как их собственные файлы в директории на пути MATLAB.

function [c,f,s] = pdefun(x,t,u,dudx) % Equation to solve
c = [1; 1];
f = [0.024; 0.17] .* dudx;
y = u(1) - u(2);
F = exp(5.73*y)-exp(-11.47*y);
s = [-F; F];
end
% ---------------------------------------------
function u0 = pdeic(x) % Initial Conditions
u0 = [1; 0];
end
% ---------------------------------------------
function [pl,ql,pr,qr] = pdebc(xl,ul,xr,ur,t) % Boundary Conditions
pl = [0; ul(2)];
ql = [1; 0];
pr = [ur(1)-1; 0];
qr = [0; 1];
end
% ---------------------------------------------

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

Похожие темы

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