must return a column vector

1 view (last 30 days)
ryan
ryan on 11 Dec 2014
Commented: Star Strider on 11 Dec 2014
I know this is a dimension issue, but i just can't see it at the moment. could someone help me out.
n=12;
x=zeros(1,2*n);
ct=zeros(n,n);
l=3; % in m
t=6;
%%%%%%%%%mass matrix %%%%%%%%%%%%%%
g=5000;
m=[ g,0,0,0,0,0,0,0,0,0,0,0;
0,g,0,0,0,0,0,0,0,0,0,0;
0,0,g,0,0,0,0,0,0,0,0,0;
0,0,0,g,0,0,0,0,0,0,0,0;
0,0,0,0,g,0,0,0,0,0,0,0;
0,0,0,0,0,g,0,0,0,0,0,0;
0,0,0,0,0,0,g,0,0,0,0,0;
0,0,0,0,0,0,0,g,0,0,0,0;
0,0,0,0,0,0,0,0,g,0,0,0;
0,0,0,0,0,0,0,0,0,g,0,0;
0,0,0,0,0,0,0,0,0,0,g,0;
0,0,0,0,0,0,0,0,0,0,0,g];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%Damping ratio matrix%%%%%%%
h=0.1;
zeta=[ h,0,0,0,0,0,0,0,0,0,0,0;
0,h,0,0,0,0,0,0,0,0,0,0;
0,0,h,0,0,0,0,0,0,0,0,0;
0,0,0,h,0,0,0,0,0,0,0,0;
0,0,0,0,h,0,0,0,0,0,0,0;
0,0,0,0,0,h,0,0,0,0,0,0;
0,0,0,0,0,0,h,0,0,0,0,0;
0,0,0,0,0,0,0,h,0,0,0,0;
0,0,0,0,0,0,0,0,h,0,0,0;
0,0,0,0,0,0,0,0,0,h,0,0;
0,0,0,0,0,0,0,0,0,0,h,0;
0,0,0,0,0,0,0,0,0,0,0,h];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%Stiffness matrix%%%%%%%%%%
b=4*((3*100)/l); %kN/m
k=[ 2*b,-b,0,0,0,0,0,0,0,0,0,0;
-b,2*b,-b,0,0,0,0,0,0,0,0,0;
0,-b,2*b,-b,0,0,0,0,0,0,0,0;
0,0,-b,2*b,-b,0,0,0,0,0,0,0;
0,0,0,-b,2*b,-b,0,0,0,0,0,0;
0,0,0,0,-b,2*b,-b,0,0,0,0,0;
0,0,0,0,0,-b,2*b,-b,0,0,0,0;
0,0,0,0,0,0,-b,2*b,-b,0,0,0;
0,0,0,0,0,0,0,-b,2*b,-b,0,0;
0,0,0,0,0,0,0,0,-b,2*b,-b,0;
0,0,0,0,0,0,0,0,0,-b,2*b,-b;
0,0,0,0,0,0,0,0,0,0,-b,2*b];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
a=m^(-0.5);
kt=a*k*a;
[p,q]=eig(kt);
omega=sqrt(q);
for zt=1:n
ct(zt,zt)=zeta(zt)*2*omega(zt,zt);
end
af=a^(-1);
c=af*ct*af;
v0=[0;0;0;0;0;0;0;0;0;0;0;0];%initial velocities
x0=[0;0;0;0;0;0;0;0;0;0;0;0];%initial displacement
ti=0; %initial time
tf=100; %final time
init=[v0,x0]; %initial velocity and displacement initialize
[t,ny]=ode45(@yd,[ti,tf],init);
for kl=1:n
plot(t,ny(:,kl+n),'b');
title(['Response x',num2str(kl)]);
xlabel('Time, seconds');
ylabel('Displacement, m');
grid
figure(gcf)
pause
end
here is the file it references
if true
function yd=yd(t,x)
yd(1)=(-((800)*x(13)+(68.1858)*x(1)*(-400)*x(14))+(0))/(5000);
yd(2)=(-((-400)*x(13)*(800)*x(14)*(-400)*x(15))+(0))/(5000);
yd(3)=(-((-400)*x(14)*(800)*x(15)*(-400)*x(16))+(0))/(5000);
yd(4)=(-((-400)*x(15)*(800)*x(16)*(-400)*x(17))+(0))/(5000);
yd(5)=(-((-400)*x(16)*(800)*x(17)*(-400)*x(18))+(0))/(5000);
yd(6)=(-((-400)*x(17)*(800)*x(18)*(-400)*x(19))+(0))/(5000);
yd(7)=(-((-400)*x(18)*(800)*x(19)*(-400)*x(20))+(0))/(5000);
yd(8)=(-((-400)*x(19)*(800)*x(20)*(-400)*x(21))+(0))/(5000);
yd(9)=(-((-400)*x(20)*(800)*x(21)*(-400)*x(22))+(0))/(5000);
yd(10)=(-((-400)*x(21)*(800)*x(22)*(-400)*x(23))+(0))/(5000);
yd(11)=(-((-400)*x(22)*(800)*x(23)*(-400)*x(24))+(0))/(5000);
yd(12)=(-((-400)*x(23)*(800)*x(24))+(5000*sin(3*t)))/(5000);
yd(13)=x(1);
yd(14)=x(2);
yd(15)=x(3);
yd(16)=x(4);
yd(17)=x(5);
yd(18)=x(6);
yd(19)=x(7);
yd(20)=x(8);
yd(21)=x(9);
yd(22)=x(10);
yd(23)=x(11);
yd(24)=x(12);
end

Accepted Answer

Star Strider
Star Strider on 11 Dec 2014
I didn’t run your code, but see if making the first statement after the function statement a preallocation of ‘yd’ to be a column vector will help:
yd = zeros(24,1);
Also, it’s not good practice to name your function the same as the variables within it or the variable it returns. Change the function name to ‘ydf’ or something that is meaningful in the context of your code. Also change the name of the .m file containing it and of course the reference in ode45.

More Answers (1)

ryan
ryan on 11 Dec 2014
yup that got rid of all errors. thanks.

Products

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!