How to fit data to a forced oscillation model?

I am trying to fit data to the differential equation y(t)'' + ay'(t) + by(t) = c cos (2t ) with lsqcurvefit but I don't get a good fit. I think the problem is in the initial values, but I don't know how to estimate them.

function Y=model(theta,t)
    
    y0 = [theta(1); theta(2)];    
    [T_out,Y_out]=ode45(@DifEq,t,y0);
    Y=Y_out;
    
    % y1 =y ; y2 = y'
    % y1' = y2
    % y2' = -a*y2 -b*y1 +c*cos(omega*t)
    function dY=DifEq(t,y)
        dYdt=zeros(2,1);
        dYdt=[y(2);-theta(3).*y(2)-theta(4).*y(1)+theta(5).*cos(2.*t)];
        dY=dYdt;
    end

end

DATA = importdata('Data.txt',' ');
t = DATA(:,1);
p = DATA(:,2);
v = DATA(:,3);

theta0=[1, 1, 0.1, 10, 1000];

lb=[-inf, -inf, 0, 0, 0 ];
ub=[inf, inf, inf, inf, inf];

[theta,resnorm,res,exitflag,output,lambda,J]=lsqcurvefit(@model,theta0,t,[p,v],lb,ub);

tv = linspace(min(t), max(t), 10*size(t,1))';

Yfit= model(theta, tv);
Pfit=Yfit(:,1);
Vfit=Yfit(:,2);

figure;
subplot(2,1,1);
plot(t, p)
hold on
hlp = plot(tv, Pfit);
hold off
grid
xlabel('Time')
ylabel('Position')

subplot(2,1,2);
plot(t, v)
hold on
hlp = plot(tv, Vfit);
hold off
grid
xlabel('Time')
ylabel('Velocity')

[Mathieu NOE]

hello

your data are quite noisy , or said differently, the second order DE with a sin / cos excitation is not likely to match your experimental data , so I doubt lsqcurvefit will do any miracle anyway

Answer 1

[Star Strider]

I encourage you to experiment with the file I have attached.