How to solve damped forced vibration analysis using ode45 ?
18 views (last 30 days)
Show older comments
SOMYA RANJAN PATRO
on 23 Jul 2021
Commented: Star Strider
on 23 Jul 2021
I am trying to solve equation of motion for a damped forced vibration analysis for a SDOF system. My input Force is function of time having random white noise. I am using ode45 to solve the problem. Here is my code.
% Main code
clear variables
close all
clc
t = 0:20; % Time duration
x = [1; 6]; % Initial conditions
[t, x] = ode45(@eqm,t,x); % ode45 command
Below is the function eqm
% Function
function dxdt = eqm(t, x)
F = 1:20 % For simplicity I have taken force as an array
m = 1; % mass
c = 1; % damping
k = 1; % stiffness
dxdt = [x(2); F/m - c/m*x(2) - k*x(1)];
end
After I run this code it shows error "Error using vertcat. Dimensions of arrays being concatenated are not consistent". I checked online but most of the examples contains only periodic functions of force such as F = sin(wt) but I didn't find anything for random load. Any help will be greatly appreciated.
0 Comments
Accepted Answer
Star Strider
on 23 Jul 2021
One problem is with the ‘dxdt’ assignment. MATLAB interprets the spaces as delimiters, so the second row has 3 columns, as MATLAB sees it. Put parentheses around the second row terms (to preserve readability) and after providing ‘F’ as an argument rather than internally as a vector, and a few other tweaks, it works —
t = 0:20; % Time duration
x0 = [1; 6]; % Initial conditions
F = 1:20; % For simplicity I have taken force as an array
for k = 1:numel(F)
[t, x] = ode45(@(t,x)eqm(t,x,F(k)),t,x0); % ode45 command
tv{k} = t;
xv{k} = x;
end
figure
for k = 1:numel(F)
subplot(5,4,k)
plot(tv{k}, xv{k})
title(sprintf('k = %2d',k))
grid
end
% Function
function dxdt = eqm(t, x, F)
m = 1; % mass
c = 1; % damping
k = 1; % stiffness
dxdt = [x(2); (F/m - c/m*x(2) - k*x(1))];
end
.
4 Comments
More Answers (0)
See Also
Categories
Find more on Numerical Integration and 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!