В этом примере показано, как решить систему дифференциальных уравнений с частными производными, которая использует ступенчатые функции в начальных условиях.
Рассмотрите УЧП
Уравнения включают постоянные параметры , , , , и , и заданы для и . Эти уравнения возникают в математической модели первых шагов связанного с опухолью ангиогенеза [1]. представляет плотность ячейки эндотелиальных клеток, и концентрация белка они выпускают в ответ на опухоль.
Эта проблема имеет постоянное, устойчивое состояние когда
.
Однако анализ устойчивости предсказывает эволюцию системы к неоднородному решению [1]. Таким образом, ступенчатые функции используются в качестве начальных условий, чтобы встревожить устойчивое состояние и стимулировать эволюцию системы.
Граничные условия требуют, чтобы оба компонента решения имели нулевой поток в и .
Чтобы решить эту систему уравнений в MATLAB, необходимо закодировать уравнения, начальные условия и граничные условия, затем выбрать подходящую mesh решения прежде, чем вызвать решатель pdepe
. Вы любой может включать необходимые функции как локальные функции в конце файла (как сделано здесь) или сохранить их как отдельные, именованные файлы в директории на пути MATLAB.
Прежде чем можно будет закодировать уравнение, необходимо убедиться что именно в форме pdepe
решатель ожидает:
С тех пор существует два уравнения в системе УЧП, система УЧП может быть переписана как
Значения коэффициентов в уравнении затем
(только диагональные значения)
Теперь можно создать функцию, чтобы закодировать уравнение. Функция должна иметь подпись [c,f,s] = angiopde(x,t,u,dudx)
:
x
независимая пространственная переменная.
t
независимая переменная времени.
u
зависимая переменная, дифференцируемая относительно x
и t
. Это - двухэлементный вектор где u(1)
и u(2)
.
dudx
частичная пространственная производная . Это - двухэлементный вектор где dudx(1)
и dudx(2)
.
Выходные параметры c
F
, и s
соответствуйте коэффициентам в стандартной форме уравнения УЧП, ожидаемой pdepe
.
В результате уравнения в этом примере могут быть представлены функцией:
function [c,f,s] = angiopde(x,t,u,dudx) d = 1e-3; a = 3.8; S = 3; r = 0.88; N = 1; c = [1; 1]; f = [d*dudx(1) - a*u(1)*dudx(2) dudx(2)]; s = [S*r*u(1)*(N - u(1)); S*(u(1)/(u(1) + 1) - u(2))]; end
(Примечание: Все функции включены как локальные функции в конце примера.)
Затем запишите функцию, которая возвращает начальное условие. Начальное условие применяется в первой временной стоимости и вводит значение и для любого значения x. Используйте функциональную подпись u0 = angioic(x)
записать функцию.
Эта проблема имеет постоянное, устойчивое состояние когда
.
Однако анализ устойчивости предсказывает эволюцию системы к негомогенному решению [1]. Так, ступенчатые функции используются в качестве начальных условий, чтобы встревожить устойчивое состояние и стимулировать эволюцию системы.
Соответствующая функция
function u0 = angioic(x) u0 = [1; 0.5]; if x >= 0.3 && x <= 0.6 u0(1) = 1.05 * u0(1); u0(2) = 1.0005 * u0(2); end end
Теперь запишите функцию, которая оценивает граничные условия
Для проблем, созданных на интервале , граничные условия применяются ко всем и также или . Стандартная форма для граничных условий, ожидаемых решателем,
Для , уравнение граничного условия
Таким образом, коэффициенты:
.
Для граничные условия являются тем же самым, таким образом, и .
Граничная функция должна использовать функциональную подпись [pl,ql,pr,qr] = angiobc(xl,ul,xr,ur,t)
, где:
Входные параметры xl
и ul
соответствовать и для левого контура.
Входные параметры xr
и ur
соответствовать и для правильного контура.
t
независимая переменная времени.
Выходные параметры pl
и ql
соответствовать и для левого контура ( для этой проблемы).
Выходные параметры pr
и qr
соответствовать и для правильного контура ( для этой проблемы).
Граничные условия в этом примере представлены функцией:
function [pl,ql,pr,qr] = angiobc(xl,ul,xr,ur,t) pl = [0; 0]; ql = [1; 1]; pr = pl; qr = ql; end
Долговременный интервал требуется, чтобы видеть ограничивающее поведение уравнений, так используйте 10 точек в интервале . Кроме того, предельное распределение варьируется только на приблизительно 0,1% на интервале , таким образом, относительно прекрасная пространственная mesh с 50 точками является соответствующей.
x = linspace(0,1,50); t = linspace(0,200,10);
Наконец, решите уравнение с помощью симметрии , уравнение УЧП, начальное условие, граничные условия и сетки для и .
m = 0; sol = pdepe(m,@angiopde,@angioic,@angiobc,x,t);
pdepe
возвращает решение в трехмерном массиве sol
, где sol(i,j,k)
аппроксимирует k
компонент th решения оцененный в t(i)
и x(j)
. Извлеките компоненты решения в отдельные переменные.
n = sol(:,:,1); c = sol(:,:,2);
Создайте объемную поверхностную диаграмму компонентов решения и построенный в выбранной mesh указывает для и .
surf(x,t,c) title('c(x,t): Concentration of Fibronectin') xlabel('Distance x') ylabel('Time t')
surf(x,t,n) title('n(x,t): Density of Endothelial Cells') xlabel('Distance x') ylabel('Time t')
Теперь постройте только итоговые распределения решений в . Эти графики соответствуют рисункам 3 и 4 в [1].
plot(x,n(end,:))
title('Final distribution of n(x,t_f)')
plot(x,c(end,:))
title('Final distribution of c(x,t_f)')
[1] Humphreys, M.E. и M.A.J. Капеллан. "Математическая модель первых шагов связанного с опухолью ангиогенеза: Капиллярное формирование ростка и вторичное ответвление". Журнал IMA Математики, Прикладной в Medicine & Biology, 13 (1996), стр 73-98.
Перечисленный здесь локальные функции помощника, которые решатель УЧП pdepe вызывает, чтобы вычислить решение. В качестве альтернативы можно сохранить эти функции как их собственные файлы в директории на пути MATLAB.
function [c,f,s] = angiopde(x,t,u,dudx) % Equation to solve d = 1e-3; a = 3.8; S = 3; r = 0.88; N = 1; c = [1; 1]; f = [d*dudx(1) - a*u(1)*dudx(2) dudx(2)]; s = [S*r*u(1)*(N - u(1)); S*(u(1)/(u(1) + 1) - u(2))]; end % --------------------------------------------- function u0 = angioic(x) % Initial Conditions u0 = [1; 0.5]; if x >= 0.3 && x <= 0.6 u0(1) = 1.05 * u0(1); u0(2) = 1.0005 * u0(2); end end % --------------------------------------------- function [pl,ql,pr,qr] = angiobc(xl,ul,xr,ur,t) % Boundary Conditions pl = [0; 0]; ql = [1; 1]; pr = pl; qr = ql; end % ---------------------------------------------