Not enough input agruments
4 views (last 30 days)
Show older comments
% Set the circuit parameters and initial conditions
res = 10 * 9;
C = 0.0001/9;
L = 0.01 * 4;
V0 = 0;
I = 1;
omega = 2 * pi * 0;
I0 = [1; 0];
% Solve the system using ode45()
[t_ode45, I_ode45] = ode45(@(t, I) Question3(t, I, res, L, C, V0, omega), [0, 40 * L / res], I0);
% Set the time step and the number of time steps
h = 0.001;
N = ceil(40 * L / res / h);
% Use the forward Euler method to solve the system
t_euler = zeros(N, 1);
I_euler = zeros(N, 2);
I_euler(1, :) = I0;
for n = 1:N-1
t_euler(n+1) = t_euler(n) + h;
I_euler(n+1, :) = I_euler(n, :) + h * Question3(t_euler(n), I_euler(n, :), res, L, C, V0, omega)';
end
% Use the backward Euler method to solve the system
t_beuler = zeros(N, 1);
I_beuler = zeros(N, 2);
I_beuler(1, :) = I0;
for n = 1:N-1
t_beuler(n+1) = t_beuler(n) + h;
f = @(I) I - I_beuler(n, :) - h * Question3(t_beuler(n+1), I, res, L, C, V0, omega)';
I_beuler(n+1, :) = fsolve(f, I_beuler(n, :));
end
% Plot the results
figure;
plot(t_ode45, I_ode45(:, 1), 'b-', 'LineWidth', 1.5);
hold on;
plot(t_euler, I_euler(:, 1), 'r--', 'LineWidth', 1.5);
plot(t_beuler, I_beuler(:, 1), 'm-.', 'LineWidth', 1.5);
xlabel('Time (s)');
ylabel('Current (A)');
title('RLC Circuit in Series');
legend('ode45', 'Forward Euler', 'Backward Euler');
% Define the function for the RLC circuit in series
function dIdt = Question3(t, I, res, L, C, V0, omega)
dIdt = [I(2); (V0 * omega * cos(omega * t) - res * (I(2) - I(1)) / C) / L];
end
1 Comment
Answers (1)
Chunru
on 28 Apr 2023
Function should be placed at the end of script or separate file:
% Set the circuit parameters and initial conditions
res = 10 * 9;
C = 0.0001/9;
L = 0.01 * 4;
V0 = 0;
I = 1;
omega = 2 * pi * 0;
I0 = [1; 0];
% Solve the system using ode45()
[t_ode45, I_ode45] = ode45(@(t, I) Question3(t, I, res, L, C, V0, omega), [0, 40 * L / res], I0);
% Set the time step and the number of time steps
h = 0.001;
N = ceil(40 * L / res / h);
% Use the forward Euler method to solve the system
t_euler = zeros(N, 1);
I_euler = zeros(N, 2);
I_euler(1, :) = I0;
for n = 1:N-1
t_euler(n+1) = t_euler(n) + h;
I_euler(n+1, :) = I_euler(n, :) + h * Question3(t_euler(n), I_euler(n, :), res, L, C, V0, omega)';
end
% Use the backward Euler method to solve the system
t_beuler = zeros(N, 1);
I_beuler = zeros(N, 2);
I_beuler(1, :) = I0;
for n = 1:N-1
t_beuler(n+1) = t_beuler(n) + h;
f = @(I) I - I_beuler(n, :) - h * Question3(t_beuler(n+1), I, res, L, C, V0, omega)';
I_beuler(n+1, :) = fsolve(f, I_beuler(n, :));
end
% Plot the results
figure;
plot(t_ode45, I_ode45(:, 1), 'b-', 'LineWidth', 1.5);
hold on;
plot(t_euler, I_euler(:, 1), 'r--', 'LineWidth', 1.5);
plot(t_beuler, I_beuler(:, 1), 'm-.', 'LineWidth', 1.5);
xlabel('Time (s)');
ylabel('Current (A)');
title('RLC Circuit in Series');
legend('ode45', 'Forward Euler', 'Backward Euler');
% Function should be placed at the end of script or separate file
function dIdt = Question3(t, I, res, L, C, V0, omega)
dIdt = [I(2); (V0 * omega * cos(omega * t) - res * (I(2) - I(1)) / C) / L];
end
0 Comments
See Also
Categories
Find more on Creating and Concatenating Matrices 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!