Error running for loop

Hello,

I am working with the following MATLAB code, it works properly if I want to calculate the values at some fixed points, but if I run the same code using for loop (for n=0:4), it reports an error (Error using vertcat

Dimensions of arrays being concatenated are not consistent.)

How could I resolve this issue.

I would be very grateful for your kind support.

Here's my MATLAB code:

syms P
l=100; a=50; b=20;  h0=50; E=200; h1=50;%n=1.5;
for n=1:0.5:4
I0=b*h0^3/12; I1=b*h1^3/12;
L0=sqrt(P/(E*I0-P*n));%L0=simplify(L00);
L1=sqrt(P/(E*I1-P*n));%L1=simplify(L11);
v=0.3;
s=0.5; 
f=1.93-3.07*s+14.53*s^2-25.11*s^3+25.8*s^4;
k=(E*I1)/(6*pi*h1*f*(1-v^2)); 
a11=cos(L1*l) - cos(L1*a) - L1*a*sin(L1*l) + L1*l*sin(L1*l);
a12=sin(L1*l) - sin(L1*a) + L1*a*cos(L1*l) - L1*l*cos(L1*l);
a13=cos(L0*a) - 1;
a14=sin(L0*a) - L0*a;
a21=-(L1*sin(L1*a) - L1*sin(L1*l)) - L1^2*cos(L1*a)*((P*n - E*I1)/k);
a22=(L1*cos(L1*a) - L1*cos(L1*l)) - L1^2*sin(L1*a)*((P*n - E*I1)/k);
a23=L0*sin(L0*a);
a24=(L0 - L0*cos(L0*a));
a31=-L1^2*cos(L1*a)*(P*n - E*I1);
a32=-L1^2*sin(L1*a)*(P*n - E*I1);
a33=L0^2*cos(L0*a)*(P*n - E*I0);
a34=L0^2*sin(L0*a)*(P*n - E*I0);
a41=L1^3*sin(L1*a)*(P*n - E*I1);
a42=-L1^3*cos(L1*a)*(P*n - E*I1);
a43=-L0^3*sin(L0*a)*(P*n - E*I0);
a44=L0^3*cos(L0*a)*(P*n - E*I0);
A=[a11 a12 a13 a14;a21 a22 a23 a24;a31 a32 a33 a34;a41 a42 a43 a44];
delta1=det(A);
delta2=subs(delta1);
delta=simplify(delta2);
F=matlabFunction(delta);
%FF=fplot(F)
%P=vpasolve(delta,P,1e3)
%P=feval(F,1e3 )% sol=fzero(F,1e3)
P=[];
for j = 0.01:10
    PP = fzero(F, j);
    P = [P; PP];
end
P_cr = min(P(P > 0)) 
D=plot(n,P_cr,'k.','linewidth',1.5)
 %axis([0 4 0 1e6])
DD=double(D)
hold on
end

Answer 1

[Cameron B]

if true
  syms P
l=100; a=50; b=20;  h0=50; E=200; h1=50;%n=1.5;
for n=1:0.02:4
I0=b*h0^3/12; I1=b*h1^3/12;
L0=sqrt(P/(E*I0-P*n));%L0=simplify(L00);
L1=sqrt(P/(E*I1-P*n));%L1=simplify(L11);
v=0.3;
s=0.5; 
f=1.93-3.07*s+14.53*s^2-25.11*s^3+25.8*s^4;
k=(E*I1)/(6*pi*h1*f*(1-v^2)); 
a11=cos(L1*l) - cos(L1*a) - L1*a*sin(L1*l) + L1*l*sin(L1*l);
a12=sin(L1*l) - sin(L1*a) + L1*a*cos(L1*l) - L1*l*cos(L1*l);
a13=cos(L0*a) - 1;
a14=sin(L0*a) - L0*a;
a21=-(L1*sin(L1*a) - L1*sin(L1*l)) - L1^2*cos(L1*a)*((P*n - E*I1)/k);
a22=(L1*cos(L1*a) - L1*cos(L1*l)) - L1^2*sin(L1*a)*((P*n - E*I1)/k);
a23=L0*sin(L0*a);
a24=(L0 - L0*cos(L0*a));
a31=-L1^2*cos(L1*a)*(P*n - E*I1);
a32=-L1^2*sin(L1*a)*(P*n - E*I1);
a33=L0^2*cos(L0*a)*(P*n - E*I0);
a34=L0^2*sin(L0*a)*(P*n - E*I0);
a41=L1^3*sin(L1*a)*(P*n - E*I1);
a42=-L1^3*cos(L1*a)*(P*n - E*I1);
a43=-L0^3*sin(L0*a)*(P*n - E*I0);
a44=L0^3*cos(L0*a)*(P*n - E*I0);
A=[a11 a12 a13 a14;a21 a22 a23 a24;a31 a32 a33 a34;a41 a42 a43 a44];
delta1=det(A);
delta2=subs(delta1);
delta=simplify(delta2);
F=matlabFunction(delta);
%FF=fplot(F)
%P=vpasolve(delta,P,1e3)
P1=[];
for j = 0.01:10
  PP = fzero(F, j);
  P1 = [P1; PP];
end
P_cr = min(P1(P1 > 0));
D=plot(n,P_cr,'k.','linewidth',1.5);
%axis([0 4 0 1e6])
DD=double(D);
hold on
end
end

[Cameron B]

Delete the “if true” at the beginning and the final “end”. I think it’s a formatting thing when I do this on my phone.

[HINA]

Yes, it works for me now. Thank you so much for your time and support. I have one more question regarding this code, that would be great if you could give any suggestion for that.

Since I am interested to get the smallest positive value of P, my question is that the syntax that i used to get smallest positive value is correct or not? as after running this file still I am not able to get the desired solutions.

I am talking about the following part:

P1=[];

for j = 0.01:10

PP = fzero(F, j);

P1 = [P1; PP];

end

P_cr = min(P1(P1 > 0));

[Cameron B]

Yea I was looking at the output and it seems as if something is off. I think the way you went about finding the lowest positive value is ok, but I’m not sure about the other code. I know this is related to fracture toughness calculations for 3 point bends, but something does seem to be off a bit. However I can’t be sure. Try this:

for j = 0.001:0.001:10 PP = fzero(F, j); P1 = [P1; PP]; end

Or you can try something else for the j range.

Answer 2

[Colo]

Hi,

I copied your code and ran it. First loop works fine because P is 1x1sym.

in line 38 you change the size of P to 10x1: P=[P;PP]

Thus, a11 becomes 10x10 and a41 10x1 during the following iteration. So the dimensions cannot match in A.

The solution is to keep P as a 1x1sym.

You can either add a line at the beginning of the loop declaring P again as syms:

syms P

or use another variable instead of deleting P in line 36 and the following nested for loop.

I hope this helped :)

[HINA]

Thank you so much for your answer, but Can you please help me how to keep P as 1x1sym. As I am not an expert in MATLAB. I 'll be very thankful to you for that.