% Oppgave 5.2

% Definisjon av diverse parametere

m = 0.2533;
k = 10.0;

omega_n = sqrt(k/m)

p_0 = 10;
t_d = 0.6;
omega = pi/t_d
t_total = 5.0;
delta_t = 0.1;


% Definisjon av antall diskrete tidspunkter
% 'round' betyr at resultater rundes av til n?rmeste heltall

ant_pkt = round(t_total/delta_t) + 1


% Initialisering av konstanter som inng?r i metoden fra seksjon 5.2

A = cos(omega_n*delta_t);
B = sin(omega_n*delta_t)/omega_n;
C = (1 - cos(omega_n*delta_t))/(2*k);
D = (1 - cos(omega_n*delta_t))/(2*k);

A_mark = -omega_n*sin(omega_n*delta_t);
B_mark = cos(omega_n*delta_t);
C_mark = sin(omega_n*delta_t)*omega_n/(2*k);
D_mark = sin(omega_n*delta_t)*omega_n/(2*k);


% Deklarasjon vektorer (eller  n x 1-matriser) for diskrete tidspunkter,
% lastverdier, samt verdier for u og u_dot.
% 'zeros' betyr at vektoren fylles med 0-er.

t = (0:delta_t:t_total);
load = zeros(1,ant_pkt);

u = zeros(1,ant_pkt);
u_dot = zeros(1,ant_pkt);


% L?kke for beregning av lastvektoren

for i = 1:ant_pkt

  if t(i) > t_d
    load(i) = 0;
  else
    load(i) = p_0*sin(omega*t(i));
  end

end
  

% L?sning av svingelikningen vha metoden fra seksjon 5.2
% t = 0 => u(1) = 0, u_dot(1) = 0; initialverdier - 
% oppfylt av initialiseringen ovenfor

% L?kke for beregning av l?sningen

for i = 1:(ant_pkt-1)

  u(i+1) = A*u(i) + B*u_dot(i) + C*load(i) + D*load((i+1));
  u_dot(i+1) = A_mark*u(i) + B_mark*u_dot(i) + C_mark*load(i) + D_mark*load(i+1);

end
  

% Plotting av l?sningen

plot(t,u)
xlabel('tid')
ylabel('forskyvning')
title('Svingelikning med impulslast og uten dempning')