Error ---> Brace indexing is not supported for variables of this type.
20 views (last 30 days)
Show older comments
Please any help would be much appreciated. I cant get rid of this error message
clc
clear
g = 9.81;
L = 0.1;
m = 0.5;
zeta = 2;
x0 = 0.25;
v0 = 0;
tspan = [0 10];
step_sizes = [0.1, 0.05, 0.025];
num_steps = numel(step_sizes);
x_euler = cell(1, num_steps);
x_midpoint = cell(1, num_steps);
times_euler = zeros(1, num_steps);
times_midpoint = zeros(1, num_steps);
ode_func = @(t, y) [y(2); (-zeta * y(2) - (m * g * L * y(1))) / m];
for i = 1:num_steps
tic;
[t_euler, x_euler{i}] = euler(ode_func, tspan, [x0; v0], step_sizes(i));
times_euler(i) = toc;
end
for i = 1:num_steps
tic;
[t_midpoint, x_midpoint{i}] = midpoint(ode_func, tspan, [x0; v0], step_sizes(i));
times_midpoint(i) = toc;
disp(['Length of t_midpoint{' num2str(i) '}: ' num2str(length(t_midpoint))]);
end
figure;
hold on;
for i = 1:num_steps
if numel(t_euler) == numel(t_midpoint{i})
plot(t_euler, x_euler{i}(:, 1), 'DisplayName', ['Euler (h=' num2str(step_sizes(i)) ')']);
plot(t_midpoint{i}, x_midpoint{i}(:, 1), 'DisplayName', ['Midpoint (h=' num2str(step_sizes(i)) ')']);
else
disp("Error: Unequal number of time points between Euler and Midpoint methods.");
break;
end
disp(['Size of t_midpoint{' num2str(i) '}: ' num2str(size(t_midpoint{i}))]);
disp(['Size of x_midpoint{' num2str(i) '}: ' num2str(size(x_midpoint{i}))]);
end
xlabel('Time');
ylabel('Position');
title('Comparison of Euler vs Midpoint Method');
legend('Location', 'best');
hold off;
disp('Computation times for Euler Method:');
disp(times_euler);
disp('Computation times for Midpoint Method:');
disp(times_midpoint);
function [t, y] = euler(ode_func, tspan, y0, h)
t = tspan(1):h:tspan(2);
n = numel(t);
y = zeros(n, numel(y0));
y(1, :) = y0;
for i = 2:n
t_i = t(i);
y(i, :) = y(i-1, :) + h * ode_func(t_i, y(i-1, :)')';
end
end
function [t, y] = midpoint(ode_func, tspan, y0, h)
t = tspan(1):h:tspan(2);
n = numel(t);
y = zeros(n, numel(y0));
y(1, :) = y0;
for i = 2:n
t_i = t(i);
k1 = ode_func(t_i, y(i-1, :)');
k2 = ode_func(t_i + h/2, (y(i-1, :) + (h/2) * k1)');
y(i, :) = y(i-1, :) + h * k2';
end
end
0 Comments
Answers (1)
Torsten
on 16 Apr 2024 at 18:20
clc
clear
g = 9.81;
L = 0.1;
m = 0.5;
zeta = 2;
x0 = 0.25;
v0 = 0;
tspan = [0 10];
step_sizes = [0.1, 0.05, 0.025];
num_steps = numel(step_sizes);
x_euler = cell(1, num_steps);
x_midpoint = cell(1, num_steps);
times_euler = zeros(1, num_steps);
times_midpoint = zeros(1, num_steps);
ode_func = @(t, y) [y(2); (-zeta * y(2) - (m * g * L * y(1))) / m];
for i = 1:num_steps
tic;
[t_euler{i}, x_euler{i}] = euler(ode_func, tspan, [x0; v0], step_sizes(i));
times_euler(i) = toc;
end
for i = 1:num_steps
tic;
[t_midpoint{i}, x_midpoint{i}] = midpoint(ode_func, tspan, [x0; v0], step_sizes(i));
times_midpoint(i) = toc;
disp(['Length of t_midpoint{' num2str(i) '}: ' num2str(numel(t_midpoint{i}))]);
end
figure;
hold on;
for i = 1:num_steps
if numel(t_euler{i}) == numel(t_midpoint{i})
plot(t_euler{i}, x_euler{i}(:, 1), 'DisplayName', ['Euler (h=' num2str(step_sizes(i)) ')']);
plot(t_midpoint{i}, x_midpoint{i}(:, 1), 'DisplayName', ['Midpoint (h=' num2str(step_sizes(i)) ')']);
else
disp("Error: Unequal number of time points between Euler and Midpoint methods.");
break;
end
disp(['Size of t_midpoint{' num2str(i) '}: ' num2str(size(t_midpoint{i}))]);
disp(['Size of x_midpoint{' num2str(i) '}: ' num2str(size(x_midpoint{i}))]);
end
xlabel('Time');
ylabel('Position');
title('Comparison of Euler vs Midpoint Method');
legend('Location', 'best');
hold off;
disp('Computation times for Euler Method:');
disp(times_euler);
disp('Computation times for Midpoint Method:');
disp(times_midpoint);
function [t, y] = euler(ode_func, tspan, y0, h)
t = tspan(1):h:tspan(2);
n = numel(t);
y = zeros(n, numel(y0));
y(1, :) = y0;
for i = 2:n
t_i = t(i);
y(i, :) = y(i-1, :) + h * ode_func(t_i, y(i-1, :)')';
end
end
function [t, y] = midpoint(ode_func, tspan, y0, h)
t = tspan(1):h:tspan(2);
n = numel(t);
y = zeros(n, numel(y0));
y(1, :) = y0;
for i = 2:n
t_i = t(i);
k1 = ode_func(t_i, y(i-1, :)');
k2 = ode_func(t_i + h/2, (y(i-1, :) + (h/2) * k1)');
y(i, :) = y(i-1, :) + h * k2';
end
end
0 Comments
See Also
Categories
Find more on Spreadsheets in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!