Unable to perform assignment because the left and right sides have a different number of elements. ODE45

Question

josh johnson 2022-11-17

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1854298-unable-to-perform-assignment-because-the-left-and-right-sides-have-a-different-number-of-elements

评论： Torsten 2022-11-18

采纳的回答： Cris LaPierre

在 MATLAB Online 中打开

Im try to plot t,y. but i keep getting the effor.

Unable to perform assignment because the

left and right sides have a different

number of elements.

Error in testq3>MDOFODE45 (line 161)

yy(6) = (f(3) - k(9)*y(6) - k(6)*y(5) -

.2*k(6)*y(5) + .2*k(7)*y(4) - .2*k*(9)*y(6)

)/m(9);

Error in

testq3>@(t,y)MDOFODE45(t,y,m,k,Fd,FI,w,phiB,phiAC,phiD)

(line 116)

[t,y] = ode45(@(t,y) MDOFODE45(t, y, m, k,

Fd, FI, w , phiB, phiAC, phiD), tspan,y0) ;

Error in odearguments (line 90)

f0 = feval(ode,t0,y0,args{:}); % ODE15I

sets args{1} to yp0.

Error in ode45 (line 115)

odearguments(FcnHandlesUsed, solver_name,

ode, tspan, y0, options, varargin);

Error in testq3 (line 116)

[t,y] = ode45(@(t,y) MDOFODE45(t, y, m, k,

Fd, FI, w , phiB, phiAC, phiD), tspan,y0) ;

>>

clc;clear;

%---------------------------------variable -------------------

%distances to centre of gravity

I1 = 16; %distance to g side 1

I2 = 14; %distance to g side 1

I3 = 17; %distance to g side 2

I4 = 13; %distance to g side 2

M =9000*1000 ; %mass 9000 tonnes

I=1.35*10^9; %interia

Rout= 1.25 ; %outter radius of legs

Rin = Rout-.02; %inner raidus of legs

E=204*10^9; %youngs modulus

Im=(pi/4)*(Rout^4-Rin^4); %second moment of cylinder

dis1 = 41.42; %distance to pillar D

dis2 = dis1 / 2; %distance to pillar A & C

%------------------------------Case selection ------------------

% i = input('Enter input:');

i=1 % modified to run here

i = 1

switch i

case 1 %variables for sea state 1

w = 1.571; %frequency

T = 4; %Wave period

Wh = 2 ; %wave height

Fd = 3 ; %Max wave Drag (kN)

FI = 60 ; %Max interia wave (kN)

Lo = 24.956; %Deep water wave

case 2 %variables for sea state 2

w = 1.142;

T = 4;

Wh = 2 ;

Fd = 3 ;

FI = 60 ;

Lo = 47.183;

case 3 %variables for sea state 3

w = 0.898;

T = 4;

Wh = 2 ;

Fd = 3 ;

FI = 60 ;

Lo = 76.428;

case 4 %variables for sea state 2

w = 0.698;

T = 4;

Wh = 2 ;

Fd = 3 ;

FI = 60 ;

Lo = 126.341;

end

if (i == 1) %variables for sea state 1

wd = 10:5:30; % wave depth

L = Lo*sqrt(tanh(2*pi*wd/Lo)); % intermediate water waelength

wps = L / T; % wave proagation speed

keq=(3*E*Im)./(wd.^3); %lateral stiffness

elseif( i == 2 ) %variables for sea state 2

wd = 10:5:30;

L = Lo*sqrt(tanh(2*pi*wd/Lo));

wps = L / T;

keq=(3*E*Im)./(wd.^3);

elseif (i == 3) %variables for sea state 3

wd = 20:5:40;

L = Lo*sqrt(tanh(2*pi*wd/Lo));

wps = L / T;

keq=(3*E*Im)./(wd.^3);

elseif(i == 4) %variables for sea state 4

wd = 20:5:40;

L = Lo*sqrt(tanh(2*pi*wd/Lo));

wps = L / T;

keq=(3*E*Im)./(wd.^3);

end

%phi matix space allocation

phiB = zeros(1,5);

phiAC = zeros(1,5);

phiD = zeros(1,5);

for B = 1:5

pillart = 0/wps(B) ; %time to pillar B

pillartAC = dis2/wps(B) ; %time to pillar A & C

pillartD = dis1/wps(B) ; %time to pillar D

phiB(B) = (2*pi*pillart)/T; %phi at pillar B

phiAC(B) = (2*pi*pillartAC)/T; %phi at pillar A & C

phiD(B) = (2*pi*pillartD)/T; %phi at pillar D

end

%------------------------------------------matrices-------------

kmat=[];

for kx = 1:length(L)

keq=(3*E*Im)/(L(kx)^3); %lateral stiffness

kmat=[keq,kmat];

k=[4*keq 0 2*keq*(I1-I2); 0 4*keq 2*keq*(I3-I4); 2*keq*(I1-I2) 2*keq*(I3-I4) 2*keq*(I1^2+I2^2+I3^2+I4^2)];

end

c =0.2 * k;

m =[M 0 0;0 M 0;0 0 I ];

disp(k);

1.0e+09 * 0.0189 0 0.0189 0 0.0189 0.0377 0.0189 0.0377 8.5856

%ODE45 solver

y0=[0 0 0 0 0 0];% inital conditions

tspan=[0 10]; % time span

[t,y] = ode45(@(t,y) MDOFODE45(t, y, m, k, Fd, FI, w , phiB, phiAC, phiD), tspan,y0) ;

%results

% y_ = y(:, 4:end) ;

% ydot_ = y(:, 1:3) ;

figure(1)

plot (t,y(:,1)) % Plot the time

xlabel('Time (s)') % Label the plot

ylabel('Response X_1')

figure(2)

plot (t,y(:,2)) % Plot the time

xlabel('Time (s)') % Label the plot

ylabel('Response X_2')

% ODE Function

function yy = MDOFODE45(t, y, m, k, Fd, FI, w , phiB, phiAC, phiD)

I1 = 16;

I2 = 14;

I3 = 17;

I4 = 13;

FlegB = (3*Fd/4)*(sin(w*t+phiB)) - (3*sin(w*t+phiB)/3)+FI*cos(w*t+phiB); %forces at leg B

FlegAC = (3*Fd/4)*(sin(w*t+phiAC)) - (3*sin(w*t+phiAC)/3)+FI*cos(w*t+phiAC); %forces at leg A & C

FlegD = (3*Fd/4)*(sin(w*t+phiD)) - (3*sin(w*t+phiD)/3)+FI*cos(w*t+phiD); %forces at leg D

FlegBx = FlegB * sind(45); %forces at leg B in x direction

FlegBy = FlegB * cosd(45); %forces at leg B in y direction

FlegACx = FlegAC * sind(45); %forces at leg A & C in x direction

FlegACy = FlegAC * cosd(45); %forces at leg A & C in y direction

FlegDx = FlegD * sind(45); %forces at leg D in x direction

FlegDy = FlegD * cosd(45); %forces at leg D in y direction

%force matrix

f=[FlegACx+FlegBx+FlegACx+FlegDx ; FlegACy+FlegBy+FlegACy+FlegDy ; (FlegACx+FlegBx)*I1-(FlegACx+FlegDx)*I2-(FlegACy+FlegDy)*I3+(FlegBy+FlegACy)*I4];

yy = zeros(6,1);

yy(1) = y(4) ;

yy(2) = y(5) ;

yy(3) = y(6) ;

yy(4) = (f(1) - k(1)*y(1) - k(3)*y(6) - .2*k(1)*y(4) - .2*k(3)*y(6) )/m(1);

yy(5) = (f(2) - k(5)*y(2) - k(6)*y(3) - .2*k(5)*y(5) - .2*k(6)*y(6) )/m(5) ;

yy(6) = (f(3) - k(9)*y(6) - k(6)*y(5) - .2*k(6)*y(5) + .2*k(7)*y(4) - .2*k(9)*y(6) )/m(9);

end

I have also noticed that my k fuction doesnt doesnt full work because the k mat isnt fully statisfied by the keq values. any help would be greatly apprieated

3 个评论
显示 1更早的评论隐藏 1更早的评论

josh johnson 2022-11-18

Much apprecitaed !!

is there a way to get the k to run for each wave depth (wd) then have figures that show the response at each depth ?

Torsten 2022-11-18

Make a loop over k, call ODE45 in the loop for each value of k separately and save the results.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Cris LaPierre 2022-11-17

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1854298-unable-to-perform-assignment-because-the-left-and-right-sides-have-a-different-number-of-elements#answer_1103123

编辑：Cris LaPierre 2022-11-17

在 MATLAB Online 中打开

The problem is that your calculation for yy(6) returns a 3x3 matrtix, and you are assigning it to a 1x1 space. Hence, the left (1x1) and the right (3x3) have a different number of elements.

I think this is your mistake

(f(3) - k(9)*y(6) - k(6)*y(5) - .2*k(6)*y(5) + .2*k(7)*y(4) - .2*k*(9)*y(6) )/m(9)

% ^

I think you want .2*k(9)*y(6) instead

clc;clear;

%---------------------------------variable -------------------

%distances to centre of gravity

I1 = 16; %distance to g side 1

I2 = 14; %distance to g side 1

I3 = 17; %distance to g side 2

I4 = 13; %distance to g side 2

M =9000*1000 ; %mass 9000 tonnes

I=1.35*10^9; %interia

Rout= 1.25 ; %outter radius of legs

Rin = Rout-.02; %inner raidus of legs

E=204*10^9; %youngs modulus

Im=(pi/4)*(Rout^4-Rin^4); %second moment of cylinder

dis1 = 41.42; %distance to pillar D

dis2 = dis1 / 2; %distance to pillar A & C

%------------------------------Case selection ------------------

% i = input('Enter input:');

i=1 % modified to run here

i = 1

switch i

case 1 %variables for sea state 1

w = 1.571; %frequency

T = 4; %Wave period

Wh = 2 ; %wave height

Fd = 3 ; %Max wave Drag (kN)

FI = 60 ; %Max interia wave (kN)

Lo = 24.956; %Deep water wave

case 2 %variables for sea state 2

w = 1.142;

T = 4;

Wh = 2 ;

Fd = 3 ;

FI = 60 ;

Lo = 47.183;

case 3 %variables for sea state 3

w = 0.898;

T = 4;

Wh = 2 ;

Fd = 3 ;

FI = 60 ;

Lo = 76.428;

case 4 %variables for sea state 2

w = 0.698;

T = 4;

Wh = 2 ;

Fd = 3 ;

FI = 60 ;

Lo = 126.341;

end

if (i == 1) %variables for sea state 1

wd = 10:5:30; % wave depth

L = Lo*sqrt(tanh(2*pi*wd/Lo)); % intermediate water waelength

wps = L / T; % wave proagation speed

keq=(3*E*Im)./(wd.^3); %lateral stiffness

elseif( i == 2 ) %variables for sea state 2

wd = 10:5:30;

L = Lo*sqrt(tanh(2*pi*wd/Lo));

wps = L / T;

keq=(3*E*Im)./(wd.^3);

elseif (i == 3) %variables for sea state 3

wd = 20:5:40;

L = Lo*sqrt(tanh(2*pi*wd/Lo));

wps = L / T;

keq=(3*E*Im)./(wd.^3);

elseif(i == 4) %variables for sea state 4

wd = 20:5:40;

L = Lo*sqrt(tanh(2*pi*wd/Lo));

wps = L / T;

keq=(3*E*Im)./(wd.^3);

end

%phi matix space allocation

phiB = zeros(1,5);

phiAC = zeros(1,5);

phiD = zeros(1,5);

for B = 1:5

pillart = 0/wps(B) ; %time to pillar B

pillartAC = dis2/wps(B) ; %time to pillar A & C

pillartD = dis1/wps(B) ; %time to pillar D

phiB(B) = (2*pi*pillart)/T; %phi at pillar B

phiAC(B) = (2*pi*pillartAC)/T; %phi at pillar A & C

phiD(B) = (2*pi*pillartD)/T; %phi at pillar D

end

%------------------------------------------matrices-------------

kmat=[];

for kx = 1:length(L)

keq=(3*E*Im)/(L(kx)^3); %lateral stiffness

kmat=[keq,kmat];

k=[4*keq 0 2*keq*(I1-I2); 0 4*keq 2*keq*(I3-I4); 2*keq*(I1-I2) 2*keq*(I3-I4) 2*keq*(I1^2+I2^2+I3^2+I4^2)];

end

c =0.2 * k;

m =[M 0 0;0 M 0;0 0 I ];

disp(k);

1.0e+09 * 0.0189 0 0.0189 0 0.0189 0.0377 0.0189 0.0377 8.5856

%ODE45 solver

y0=[0 0 0 0 0 0];% inital conditions

tspan=[0 500]; % time span

[t,y] = ode45(@(t,y) MDOFODE45(t, y, m, k, Fd, FI, w , phiB, phiAC, phiD), tspan,y0) ;

%results

% y_ = y(:, 4:end) ;

% ydot_ = y(:, 1:3) ;

figure(1) % ##################### <--------------- modified this line

plot (t,y(:,1)) % Plot the time

xlabel('Time (s)') % Label the plot

ylabel('Response X_1')

figure(2) % ##################### <--------------- modified this line

plot (t,y(:,2)) % Plot the time

xlabel('Time (s)') % Label the plot

ylabel('Response X_2')

% ODE Function

function yy = MDOFODE45(t, y, m, k, Fd, FI, w , phiB, phiAC, phiD)

I1 = 16;

I2 = 14;

I3 = 17;

I4 = 13;

FlegB = (3*Fd/4)*(sin(w*t+phiB)) - (3*sin(w*t+phiB)/3)+FI*cos(w*t+phiB); %forces at leg B

FlegAC = (3*Fd/4)*(sin(w*t+phiAC)) - (3*sin(w*t+phiAC)/3)+FI*cos(w*t+phiAC); %forces at leg A & C

FlegD = (3*Fd/4)*(sin(w*t+phiD)) - (3*sin(w*t+phiD)/3)+FI*cos(w*t+phiD); %forces at leg D

FlegBx = FlegB * sin(45); %forces at leg B in x direction

FlegBy = FlegB * cos(45); %forces at leg B in y direction

FlegACx = FlegAC * sin(45); %forces at leg A & C in x direction

FlegACy = FlegAC * cos(45); %forces at leg A & C in y direction

FlegDx = FlegD * sin(45); %forces at leg D in x direction

FlegDy = FlegD * cos(45); %forces at leg D in y direction

%force matrix

f=[FlegACx+FlegBx+FlegACx+FlegDx ; FlegACy+FlegBy+FlegACy+FlegDy ; (FlegACx+FlegBx)*I1-(FlegACx+FlegDx)*I2-(FlegACy+FlegDy)*I3+(FlegBy+FlegACy)*I4];

yy = zeros(6,1);

yy(1) = y(4) ;

yy(2) = y(5) ;

yy(3) = y(6) ;

yy(4) = (f(1) - k(1)*y(1) - k(3)*y(6) - .2*k(1)*y(4) - .2*k(3)*y(6) )/m(1);

yy(5) = (f(2) - k(5)*y(2) - k(6)*y(3) - .2*k(5)*y(5) - .2*k(6)*y(6) )/m(5) ;

yy(6) = (f(3) - k(9)*y(6) - k(6)*y(5) - .2*k(6)*y(5) + .2*k(7)*y(4) - .2*k(9)*y(6) )/m(9);

end

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Torsten 2022-11-17

在 MATLAB Online 中打开

And most probably

FlegBx = FlegB * sind(45);                %forces at leg B in x direction
FlegBy = FlegB * cosd(45);                %forces at leg B in y direction
FlegACx = FlegAC * sind(45);              %forces at leg A & C in x direction
FlegACy = FlegAC * cosd(45);              %forces at leg A & C in y direction
FlegDx = FlegD * sind(45);                %forces at leg D in x direction
FlegDy = FlegD * cosd(45);                %forces at leg D in y direction

instead of

FlegBx = FlegB * sin(45);                %forces at leg B in x direction
FlegBy = FlegB * cos(45);                %forces at leg B in y direction
FlegACx = FlegAC * sin(45);              %forces at leg A & C in x direction
FlegACy = FlegAC * cos(45);              %forces at leg A & C in y direction
FlegDx = FlegD * sin(45);                %forces at leg D in x direction
FlegDy = FlegD * cos(45);                %forces at leg D in y direction