0

I am currently working on Matlab code to solve a second-order differential equation. From there, I convert the equation into a system of two first-order differential equations. I am unsure how solve the system of equations with the initial values provided below using Euler's method first and then using 2nd order Runge-Kutta method.

In my code below, I tried transposing the initial values array y0, but I still am receiving the following error:

%%
Error using indexing
Invalid argument at position 2. Symbolic function expected 1 input arguments but received 2.

Error in matlab_sample0 (line 49)
    y_num(:,i + 1) = y_num(:,i) + (h * f(t(i), 
 y_num(:,i)));
%%

clear;
clc;
clf;

syms f(t)

%y" = g(t) - 2y' - 2y

% input values
h = pi/50;
t0 = 0;
tfinal = 2*pi;
t = 0:h:tfinal;
n = round(tfinal/h);

% Euler's method
y0 = [1; 0];
y_num = zeros(2, n + 1);
y_num(:,1) = y0;
for i = 1:n
    y_num(:,i + 1) = y_num(:,i) + (h * f(t(i), y_num(:,i)));
end

%exact solution
exact_sol = @(t,x) [x2; -2*x2-2*x1];

y_exact = zeros(1, n + 1);
for i = 1:(n + 1)
    y_exact(i) = exact_sol(t(i));
end

%plot figure
figure;
plot(t, y_num, '-k', t, y_exact, '--r')
xlabel('t'); ylabel('y'); 
title("Euler's method Comparison between numerical and exact solutions");
  • Please post your code as text (highlighted with code button) and then ask *specific* questions about it. There are several examples of multi-variable RK integration on MATLAB Answers site as well as the doc for ode45( ). – James Tursa Apr 10 '23 at 00:02
  • @JamesTursa Hi! Sorry for the lack of clarity. I've updated my post with the code I'm currently working on and the code I referenced when discussing my previous project. Unfortunately, I'm not really sure if I can be more specific because I feel pretty lost about the entire project so I'm not even sure what questions to ask beyond what I wrote in my post. – candyms801 Apr 10 '23 at 00:32

1 Answers1

0

You can use the same loop as in the scalar version. Only you need to take care of rows and columns. If, as before, the result is formatted as list of columns, now y_num = zeros(2,n+1), then you need to address the full column as state vector, else you just get a data place somewhere else,

    y_num(:, i + 1) = y_num(:, i) + (h * f(t(i), y_num(:, i))); 

and similar in the 2nd order method.

Lutz Lehmann
  • 25,219
  • 2
  • 22
  • 51
  • Thank you for the clarification. I tried to edit the code and run it, but I encountered an error stating that the left and right sides of my y_num(1) = y0 line do not have the same amount of elements. I'm not sure where I went wrong! – candyms801 Apr 10 '23 at 17:38
  • This requires the same modification, `y_num(:, 1)=y0`. `y_num(1)` is only the first scalar data element in the structure (it is a design "error" of matlab that the flat linear data array behind the matrix is accessible at the same level as the matrix elements). – Lutz Lehmann Apr 10 '23 at 17:42
  • Thank you so much! I think I only have one more question and it's regarding inputting the initial values. I set y0 = [1 0] but now I am getting an error in my for loop because only one initial value was expected. Is there a way for me to use both initial values? – candyms801 Apr 10 '23 at 17:56
  • This is a row vector, the states are column vectors. So either transpose it `[1 0]'` or use `[1; 0]` – Lutz Lehmann Apr 10 '23 at 18:45
  • @candyms801 Can you post your current code? Both initial conditions and propagation equations need to be using column vectors everywhere for the states. – James Tursa Apr 10 '23 at 18:47
  • @LutzLehmann I've updated my post with the section of code I'm currently working on. I opted to use the latter suggestion for formatting the initial values but I still have an error :( – candyms801 Apr 10 '23 at 19:09
  • Remove the `syms f(t)` line and add a proper definition for `f(t,y)`. In your former code you had this for the ode45 integration. – Lutz Lehmann Apr 10 '23 at 19:15