Using MatLAB to solve 2nd Order DDE, weird output graph but not sure why?
    3 views (last 30 days)
  
       Show older comments
    
Hi, I have an assignment for Uni to solve the 2nd Order ODE:
m*x''(t) + c*x'(t) + kx = P(t)
given initial conditions, x(t=0) = 1.0, dx/dt(t=0) = 0.0, and the values of parameters as m=1 kg, c=0.1 N•s/m, k=42 N/m , P(t)=0 N, and between 0: t : 10, and dt=0.2s,
I have a code which gives me a graph that is incorrect, however I am very new to using MATLAB and I am probably missing something obvious in my code, so any help would be very appreciated.
Code is attatched below!
Thanks!
clc;  %Clear Windows%
clear all;
f_1= @(t,V,x) -0.1*V-42*x;   %Define Functions%
f_2= @(t,V,x) V;
h=0.2;      %Define Time Step and Time Range%
t=0:h:10;
x = zeros(1,length(t));     %Create Matricies for loop%
V = zeros(1,length(t));     
x(1)=1; %Input initial conditions%
V(1)=0;
for i=1:50  %Computational loop%
      k_1_1 = h*f_1(  t(i),        x(i),            V(i));
      k_1_2 = h*f_2(  t(i),        x(i),            V(i));
      k_2_1 = h*f_1(  t(i)+h/2  ,  x(i)+k_1_1/2  ,  V(i)+k_1_2/2  );
      k_2_2 = h*f_2(  t(i)+h/2  ,  x(i)+k_1_1/2  ,  V(i)+k_1_2/2  );
      k_3_1 = h*f_1(  t(i)+h/2  ,  x(i)+k_2_1/2  ,  V(i)+k_2_2/2  );
      k_3_2 = h*f_2(  t(i)+h/2  ,  x(i)+k_2_1/2  ,  V(i)+k_2_2/2  );
      k_4_1 = h*f_1(  t(i)+h    ,  x(i)+k_3_1    ,  V(i)+k_3_2    );
      k_4_2 = h*f_2(  t(i)+h    ,  x(i)+k_3_1    ,  V(i)+k_3_2    ); 
      xn(i+1)= x(i) + (1/6)*(k_1_1+2*k_2_1+2*k_3_1+k_4_1)*h;
      x(i+1) = xn(i+1)+ (1/6)*(k_1_1+2*k_2_1+2*k_3_1+k_4_1)*h;
end
%Exact Solution%
%%%%%%%%%%%%%%%%%%%%%
A=1;
B=0.0077154;
t = 0:0.2:10;
y=A.*exp(-0.05.*t) .* cos(6.48055.*t)+B.*exp(-0.05.*t) .* sin(6.48055.*t);
%%%%%%%%%%%%%%%%%%%%%
%Graphs%
figure1 = figure;
% Create axes
axes1 = axes('Parent',figure1,'FontSize',20,'FontName','Times New Roman');
box(axes1,'on');
grid(axes1,'on');
hold(axes1,'all');
axis([axes1],[0 1 0 11]);
% Create xlabel
xlabel('t','FontSize',20,'FontName','Times New Roman');
% Create ylabel
ylabel('x','FontSize',20,'FontName','Times New Roman');
% Create multiple lines using matrix input to plot
plot1 = plot(x,t,'Parent',axes1,'LineWidth',2);
plot2 = plot(y,t,'Parent',axes1,'LineWidth',2);
%Label%
set(plot1(1),'Marker','diamond','DisplayName','RK4');
set(plot2(1),'Marker','o','Color',[1 0 0],'DisplayName','Exact');
% Create legend
legend(axes1,'show');
0 Comments
Answers (0)
See Also
Categories
				Find more on Ordinary Differential Equations 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!