Implementing implicit method for system of ODEs

18 visualizaciones (últimos 30 días)
Romio
Romio el 25 de Jun. de 2022
Editada: Torsten el 26 de Jun. de 2022
I am trying to implement implicit Runge Kutta for solving a system of ODEs. Currently I was able to do it for scalar ODE (y'=-y), but I got stuck on how to modify it so that it work for a system of ODEs. Something like y1' = y2, y2' = -y1.
Everythin before the loop is fine if we replace the function handle f by the following function
function dydt = f(t,y)
dydt(1) = y(2);
dydt(2) = -y(1);
end
and modify the initial value to be something like
y0 = [1;0]
But How do we modify the process of solving the nonlinear system inside the loop? Also I feel I have hard-coded the method for this specific 2 stage method. Is there any way to generalize the loop for any implicit runge kutta method with size s by s?
Here is the code
n = 50; % number of time steps
t0 = 0; % initial time
T = 10; % final time
dt = (T-t0)/n; % size of timestep
y0 = 1; % initial data
m = length(y0); % dim. of problem
% intitalize time and solution vector
t = zeros(n+1,1);
y = zeros(n+1,m);
% initial conditions
t(1,1) = t0;
y(1,:) = y0(:);
% RHS of ODE
f = @(t,y) -1 * y;
% Runge Kutta coefficients
A = [1/2 + 1/(2*sqrt(3)) , 0;
-1/sqrt(3), 1/2 + 1/(2*sqrt(3))];
b = [1/2;1/2];
c = [1/2 + 1/(2*sqrt(3));1/2 - 1/(2*sqrt(3))];
s = length(c);
for i = 1:n
XI = @(xi) [xi(1) - (y(i,:) + dt*A(1,1)*f(t(i)+c(1)*dt,xi(1))...
+ dt*A(1,2)*f(t(i)+c(2)*dt,xi(2)));
xi(2) - (y(i,:) + dt*A(2,1)*f(t(i)+c(1)*dt,xi(1))...
+ dt*A(2,2)*f(t(i)+c(2)*dt,xi(2)))];
xi = fsolve(XI,zeros(s,1));
t(i+1,1) = t(i,1) + dt;
y(i+1,:) = y(i,:) + dt * (b(1)*f(t(i)+c(1)*dt,xi(1))+...
b(2)*f(t(i)+c(2)*dt,xi(2)));
end
plot(t,y)
xlabel('t')
ylabel('y(t)')

Respuestas (1)

Torsten
Torsten el 26 de Jun. de 2022
Editada: Torsten el 26 de Jun. de 2022
n = 50; % number of time steps
t0 = 0; % initial time
T = 2*pi; % final time
dt = (T-t0)/n; % size of timestep
y0 = [1;0]; % initial data
m = length(y0); % dim. of problem
% intitalize time and solution vector
t = zeros(1,n+1);
y = zeros(m,n+1);
% initial conditions
t(1) = t0;
y(:,1) = y0;
% RHS of ODE
f = @(t,y) [y(2);-y(1)];
% Runge Kutta coefficients
A = [1/2 + 1/(2*sqrt(3)) , 0;
-1/sqrt(3), 1/2 + 1/(2*sqrt(3))];
b = [1/2;1/2];
c = [1/2 + 1/(2*sqrt(3));1/2 - 1/(2*sqrt(3))];
s = length(c);
xiguess = zeros(m*s,1);
for i = 1:n
XI = @(xi) [xi(1:m) - f(t(i)+c(1)*dt,y(:,i)+dt*(A(1,1)*xi(1:m)+A(1,2)*xi(m+1:m*s)));
xi(m+1:m*s) - f(t(i)+c(2)*dt,y(:,i)+dt*(A(2,1)*xi(1:m)+A(2,2)*xi(m+1:m*s)))];
xi = fsolve(XI,xiguess);
t(i+1) = t(i) + dt;
y(:,i+1) = y(:,i) + dt * (b(1)*xi(1:m) + b(2)*xi(m+1:m*s));
xiguess = xi;
end
Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient. Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient.
plot(t,y)
xlabel('t')
ylabel('y(t)')

Categorías

Más información sobre Ordinary Differential Equations en Help Center y File Exchange.

Productos


Versión

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by