Plotting for loop values

Hello! I am trying to plot values using a for loop but everytime I run the code I get a figure with nothing on it except for axis values. I am not sure how to fix this.

EA=30000;
alpha=30;
beta=30;
gamma=180-(alpha+beta);
a=20;
n=6;
L_ab=a;
L_bc=(a*sind(alpha))/(sind(beta));
L_ac=(a*sind(alpha+beta))/(sind(beta));

for h = 0.1:0.1:1
    syms r phi
    U=((n*EA)/2)*(L_ab*((2*(r/L_ab)*sin(pi/n)-1)^2)+L_bc*(((sind(beta)/sind(alpha))*sqrt(((h/L_ab)^2)...
    -(2*((r/L_ab)^2)*cos(phi))+2*((r/L_ab)^2)))-1)^2+L_ac*((sind(beta)/sind(alpha+beta))*sqrt(((h/L_ab)^2)...
    -(2*((r/L_ab)^2)*cos(phi+(2*pi/n)))+(2*((r/L_ab)^2)))-1)^2);
    assume(phi>=0 & phi<=2*pi)
    assumeAlso(r<=2)
    assumeAlso(r,'positive')
    eqns=[diff(U,r)==0,diff(U,phi)==0]
    sol=vpasolve(eqns,[r,phi])
    sol1=sol.r
    sol2=sol.phi
    plot(h,sol1,'r'); hold on
    plot(h,sol2,'b')
end

Answer 1

[David Hill]

Try indexing

h=.1:.1:1;
for H = 1:length(h)
    syms r phi
    U=((n*EA)/2)*(L_ab*((2*(r/L_ab)*sin(pi/n)-1)^2)+L_bc*(((sind(beta)/sind(alpha))*sqrt(((h(H)/L_ab)^2)...
    -(2*((r/L_ab)^2)*cos(phi))+2*((r/L_ab)^2)))-1)^2+L_ac*((sind(beta)/sind(alpha+beta))*sqrt(((h(H)/L_ab)^2)...
    -(2*((r/L_ab)^2)*cos(phi+(2*pi/n)))+(2*((r/L_ab)^2)))-1)^2);
    assume(phi>=0 & phi<=2*pi)
    assumeAlso(r<=2)
    assumeAlso(r,'positive')
    eqns=[diff(U,r)==0,diff(U,phi)==0]
    sol=vpasolve(eqns,[r,phi])
    sol1(H)=sol.r
    sol2(H)=sol.phi
end
plot(h,sol1,'r',h,sol2,'b'); 

[Ashmika Gupta]

It did not work as I got this error:

More equations than variables is only supported for polynomial systems.

[David Hill]

Are you getting a solution?

sol=vpasolve(eqns,[r,phi]);

Is the above line passing?

[Ashmika Gupta]

No, it is not. I am not sure why as h is defined above.