Problema de valor inicial

vanderpoldemo es una función que define la ecuación de van der Pol

type vanderpoldemo

function dydt = vanderpoldemo(t,y,Mu)
%VANDERPOLDEMO Defines the van der Pol equation for ODEDEMO.

% Copyright 1984-2014 The MathWorks, Inc.

dydt = [y(2); Mu*(1-y(1)^2)*y(2)-y(1)];

La ecuación se escribe como un sistema de dos ecuaciones diferenciales de primer orden (ODE). Estas ecuaciones se evalúan para diferentes valores del parámetro $\mu$. Para una integración más rápida, debería elegir un solver apropiado según el valor de $\mu$.

Para $\mu =1$, cualquiera de los solvers de ODE de MATLAB pueden resolver correctamente la ecuación de van der Pol. El solver ode45 es un ejemplo de ellos. La ecuación se resuelve en el dominio $\left[0,20\right]$ con las condiciones iniciales $\mathit{y}\left(0\right)=2$ y $\frac{\mathrm{dy}}{\mathrm{dt}}{|}_{\mathit{t}=0}=0$.

tspan = [0 20];
y0 = [2; 0];
Mu = 1;
ode = @(t,y) vanderpoldemo(t,y,Mu);
[t,y] = ode45(ode, tspan, y0);

% Plot solution
plot(t,y(:,1))
xlabel('t')
ylabel('solution y')
title('van der Pol Equation, \mu = 1')

Para magnitudes mayores de $\mu$, el problema se vuelve rígido. Esta etiqueta es para problemas que resisten los intentos de ser evaluados con técnicas ordinarias. En su lugar, para una integración rápida son necesarios métodos numéricos especiales. Las funciones ode15s, ode23s, ode23t y ode23tb pueden resolver correctamente problemas rígidos.

Esta solución para la ecuación de van der Pol para $\mu =1000$ utiliza ode15s con las mismas condiciones iniciales. Debe ampliar el lapso de tiempo drásticamente a $\left[0,3000\right]$ para poder ver el movimiento periódico de la solución.

tspan = [0, 3000];
y0 = [2; 0];
Mu = 1000;
ode = @(t,y) vanderpoldemo(t,y,Mu);
[t,y] = ode15s(ode, tspan, y0);

plot(t,y(:,1))
title('van der Pol Equation, \mu = 1000')
axis([0 3000 -3 3])
xlabel('t')
ylabel('solution y')

Problemas de valores de límites

bvp4c y bvp5c resuelven problemas de valores de límites de las ecuaciones diferenciales ordinarias.

La función de ejemplo twoode tiene una ecuación diferencial escrita como un sistema de dos ecuaciones ODE de primer orden. La ecuación diferencial es

type twoode

function dydx = twoode(x,y)
%TWOODE  Evaluate the differential equations for TWOBVP.
%
%   See also TWOBC, TWOBVP.

%   Lawrence F. Shampine and Jacek Kierzenka
%   Copyright 1984-2014 The MathWorks, Inc.

dydx = [ y(2); -abs(y(1)) ];

La función twobc tiene las condiciones de límite para el problema: $\mathit{y}\left(0\right)=0$ y $\mathit{y}\left(4\right)=-2$.

type twobc

function res = twobc(ya,yb)
%TWOBC  Evaluate the residual in the boundary conditions for TWOBVP.
%
%   See also TWOODE, TWOBVP.

%   Lawrence F. Shampine and Jacek Kierzenka
%   Copyright 1984-2014 The MathWorks, Inc.

res = [ ya(1); yb(1) + 2 ];

Antes de llamar a bvp4c, debe proporcionar una suposición para la solución que quiera que se represente en una red. El solver entonces adapta la red a medida que refina la solución.

La función bvpinit reúne la suposición inicial en una forma que puede pasar al solver bvp4c. Para una red de [0 1 2 3 4] y una suposición constante de $\mathit{y}\left(\mathit{x}\right)=1$ y $\mathrm{y\text{'}}\left(\mathit{x}\right)=0$, la llamada a bvpinit es:

solinit = bvpinit([0 1 2 3 4],[1; 0]);

Con esta suposición inicial, puede resolver el problema con bvp4c. Evalúe la solución que devuelve bvp4c en algunos puntos con deval, y luego represente los valores resultantes.

sol = bvp4c(@twoode, @twobc, solinit);

xint = linspace(0, 4, 50);
yint = deval(sol, xint);
plot(xint, yint(1,:));
xlabel('x')
ylabel('solution y')
hold on

Este problema particular de valores de límite tiene exactamente dos soluciones. Puede obtener la otra solución si cambia la suposición inicial a $\mathit{y}\left(\mathit{x}\right)=-1$ y $\mathrm{y\text{'}}\left(\mathit{x}\right)=0$.

solinit = bvpinit([0 1 2 3 4],[-1; 0]);
sol = bvp4c(@twoode,@twobc,solinit);

xint = linspace(0,4,50);
yint = deval(sol,xint);
plot(xint,yint(1,:));
legend('Solution 1','Solution 2')
hold off

Ecuaciones diferenciales con retardo

dde23, ddesd y ddensd resuelven ecuaciones diferenciales con retardo con varios retardos. Los ejemplos ddex1, ddex2, ddex3, ddex4 y ddex5 componen un minitutorial sobre cómo usar estos solvers.

El ejemplo ddex1 muestra cómo resolver un sistema de ecuaciones diferenciales

Puede representar estas ecuaciones con la función anónima

ddex1fun = @(t,y,Z) [Z(1,1); Z(1,1)+Z(2,2); y(2)];

El historial del problema (para $\mathit{t}\le 0$) es constante:

Puede representar el historial como un vector de unos.

ddex1hist = ones(3,1);

Un vector de dos elementos representa los retardos del sistema de ecuaciones.

lags = [1 0.2];

Pase la función, los retardos, el historial de soluciones y el intervalo de integración $\left[0,5\right]$ al solver como entradas. El solver produce una solución continua sobre todo el intervalo de integración que sea susceptible de representar.

sol = dde23(ddex1fun, lags, ddex1hist, [0 5]);

plot(sol.x,sol.y);
title({'An example of Wille and Baker', 'DDE with Constant Delays'});
xlabel('time t');
ylabel('solution y');
legend('y_1','y_2','y_3','Location','NorthWest');

Ecuaciones diferenciales parciales

pdepe resuelve ecuaciones diferenciales parciales en una variable de espacio y tiempo. Los ejemplos pdex1, pdex2, pdex3, pdex4 y pdex5 componen un minitutorial sobre cómo utilizar pdepe.

Este problema de ejemplo utiliza las funciones pdex1pde, pdex1ic y pdex1bc.

pdex1pde define la ecuación diferencial

type pdex1pde

function [c,f,s] = pdex1pde(x,t,u,DuDx)
%PDEX1PDE  Evaluate the differential equations components for the PDEX1 problem.
%
%   See also PDEPE, PDEX1.

%   Lawrence F. Shampine and Jacek Kierzenka
%   Copyright 1984-2014 The MathWorks, Inc.

c = pi^2;
f = DuDx;
s = 0;

pdex1ic establece la condición inicial

type pdex1ic

function u0 = pdex1ic(x)
%PDEX1IC  Evaluate the initial conditions for the problem coded in PDEX1.
%
%   See also PDEPE, PDEX1.

%   Lawrence F. Shampine and Jacek Kierzenka
%   Copyright 1984-2014 The MathWorks, Inc.

u0 = sin(pi*x);

pdex1bc establece las condiciones de límites

type pdex1bc

function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
%PDEX1BC  Evaluate the boundary conditions for the problem coded in PDEX1.
%
%   See also PDEPE, PDEX1.

%   Lawrence F. Shampine and Jacek Kierzenka
%   Copyright 1984-2014 The MathWorks, Inc.

pl = ul;
ql = 0;
pr = pi * exp(-t);
qr = 1;

pdepe exige una discretización espacial x y un vector de tiempo t (en el que le interesa una instantánea de la solución). Resuelva el problema con una red de 20 nodos y solicite la solución en cinco valores de t. Extraiga y represente el primer componente de la solución.

x = linspace(0,1,20);
t = [0 0.5 1 1.5 2];
sol = pdepe(0,@pdex1pde,@pdex1ic,@pdex1bc,x,t);

u1 = sol(:,:,1);

surf(x,t,u1);
xlabel('x');
ylabel('t');
zlabel('u');