Question

Use Matlab and provide the code

Consider the equation for the forced oscillation of a damped system with hardening spring given by (cf. Problem 14.27). Solve this equation numerically in MatlabB with F( given by the step function excitation and ramp excitation, m 1, c= 0.5, k= 1, and μ= 0.01, 0.1, 1. Compare with the results obtained when μ-0.

Answer 1

clc;
clear all;
close all;

m =1; c=0.5;k=1; mu = [0.001,0.1,1,0];
F1 = @(t) heaviside(t); F2 = @(t) t*heaviside(t);

% let x = X1; x' = X2;

syms x(t)

for i = 1:numel(mu)
%assumed init conditions x' = 0, x = 0
x0 = 0; xp0 = 0;
fprintf(['\n assumed initial coditions x'' = ' num2str(xp0) ', x = ' num2str(x0)])

%assumed integration time range [0 30]
a = 0; b = 30;
fprintf(['\n assumed time (integration) range a = ' num2str(a) ', b = ' num2str(b)])
  
fprintf(['\n mu = ' num2str(mu(i))])
  
[V1] = odeToVectorField(m*diff(x,2)+ c*diff(x)+k*x + mu(i)*x^3 == F1(t))
[V2] = odeToVectorField(m*diff(x,2)+ c*diff(x)+k*x + mu(i)*x^3 == F2(t))

M1 = matlabFunction(V1,'vars',{'t','Y'})
M2 = matlabFunction(V2,'vars',{'t','Y'})

solutionF1 = ode45(M1,[a b],[x0 xp0]) ;
solutionF2 = ode45(M2,[a b],[x0 xp0]) ;

figure(1)
fplot(@(x) deval(solutionF1,x,1),[a b])
hold on
grid on
xlabel('t')
ylabel('x(t)')
title(['mx'''' + cx'' + kx + mux^3 = u(t)'])
legend1Info{i} = strcat('mu = ',num2str(mu(i)));
figure(2)
fplot(@(x) deval(solutionF2,x,1),[a b])
hold on
grid on
xlabel('t')
ylabel('x(t)')
title(['mx'''' + cx'' + kx + mux^3 = t u(t)'])
legend2Info{i} = strcat('mu = ',num2str(mu(i)));
end

figure(1),legend(legend1Info);
figure(2),legend(legend2Info);

assumed initial coditions x' 0, x0 assumed time (integration) range a0, b30 Y [2] Y[2]/2 heaviside (t)Y1 V2 = Y2] Y[2]/2 t heaviside (t) - Y1 e(t, Y) [Y (2);heaviside (t)-Y (1)-Y (2(1.0./2.0)] 12 e (t, Y) [Y (2) t. *heaviside (t)-Y (1)-Y (2).*(1.0./2.0) 1