How can I obtain the derivative of a vector? Displacement-Velocity but what about acceleration? I am using the space state formulation

Dear Sir/Madam,

I am using space state formulation to obtain the displacement and velocity for a SDOF system under free vibration. The code is below. Could you kindly let me know how to estimate the acceleration, please?

Thank you!

clc;  %command window
clear all;%clear workspace
m = 80;  %Defining mass of SDOF system
k = 10000;%Defining stiffness of SDOF system
global m k 

dt = .02; %Defining Time Step:
t = 0:dt:4;% Defining time vector.
y0 = [0.085 0.1]; %initial vel and disp [vel disp] %velocity. 
[tsol, yvectorl] = ode45('statespace',[0:dt:4],y0); order.


%%function statespace
function dy = testode1(t,y)
global m k
%dy(1) = -k*y(2)/m
%dy(2) = y(1)
dy=[-k*y(2)/m;y(1)];
%% fin function state space

%%plot
figure()
plot(t,yvector(:,2)) %Disp.(2)
figure()
plot(t,yvector(:,1)) %Vel. (1)

Answer 1

[Mitchell Thurston]

(Altered the code provided so that it ran properly)

clc;  %command window
clear;%clear workspace
close all; % close figures
m = 80;  %Defining mass of SDOF system
k = 10000;%Defining stiffness of SDOF system
dt = .02; %Defining Time Step:
t = 0:dt:4;% Defining time vector.
y0 = [0.085; 0.1]; %initial vel and disp [vel disp] %velocity. 
[tsol, yvectorl] = ode45(@testode1,t,y0); 

%%plot
figure()
plot(tsol,yvectorl(:,2)), title("Displacement") %Disp.(2)
figure()
plot(tsol,yvectorl(:,1)), title("Velocity") %Vel. (1)

% THIS NEEDS TO GO AT THE END OF THE FILE
%%function statespace
function dy = testode1(t,y)
global m k
dy=[-k*y(2)/m;y(1)];
end

Because you know the state space, you can find the acceleration based on the state:

accel_vector = yvectorl(:,2).*(-k/m);

but for it to apply more generally, you can just use the diff function

accel_vector = diff(yvectorl(:,1));

[Walter Roberson]

diff() assumes uniform time steps, and so would not normally be usable for ode45() results. It would be in this case because a time vector with equal steps was passed in. On the other hand, diff() would scale everything incorrectly becuase it does not know the time interval. And you would end up with one fewer entry than the original.

gradient() passing in the time and values will take the time step into account, and will use more accurate central differences except at the first and last point.

[Fedilberto Gonzalez]

Mr. Thurston, thank you very much for such a detailed answer!

Answer 2

[Walter Roberson]

estimated_acceleration = gradient(t, yvector(:,1));

By the way, please see http://www.mathworks.com/help/matlab/math/parameterizing-functions.html for information on how to stop using global.

[Fedilberto Gonzalez]

Mr. Roberson, thank you very much for your reply and help!

[Joe]

Hello, this is a little dated but I have a question regarding the solution and statespace formulation. Why is it that when the displacement from above solution is differentiated wrt time, the velocity vector derived is not same as the velocity vector from the state space solution?

That is, diff(yvector(:,2))./diff(t) is not equal to yvector(:1). Thanks in advance