Plotting a multivariable function
1 view (last 30 days)
Show older comments
Akshay Pratap Singh
on 20 Feb 2019
Commented: Akshay Pratap Singh
on 20 Feb 2019
I wrote a code for plotting a mutivariable function but getting error like "Error using fplot (line 136)
Invalid parameter '0 ...'.
Error in SFDBMDNLK (line 83)
fplot(SF,x,SF1,x).
How can I resolve it?
code:
clear all
clc
format longEng
syms y1 y2 x
phi=(pi/180)*39;
delta=(pi/180)*26;
gma=18.4;
h=4;
h1=1.91;
h2=0.088;
L=h+h1+h2;
q=0;
beta=1;
alfa=1;
Ra1=-1;
Ra2=-(alfa*(y2))^0.5;
Rp1=3*(beta*(1-y1))^0.5;
Rp2=3*(alfa*(y2))^0.5;
delma1=0.5*(1-Ra1)*delta;
delma2=-0.5*(1-Ra2)*delta;
delmp1=0.5*(Rp1-1)*delta;
delmp2=0.5*(Rp2-1)*delta;
ka1=(2*cos(phi)^2)/(cos(phi)^2*(1+Ra1)+cos(delma1)*(1-Ra1)*(1+sqrt((sin(phi+delma1)*sin(phi))/cos(delma1)))^2);
ka2=(2*cos(phi)^2)/(cos(phi)^2*(1+Ra2)+cos(delma2)*(1-Ra2)*(1+sqrt((sin(phi+delma2)*sin(phi))/cos(delma2)))^2);
kp1=1+0.5*(Rp1-1)*((cos(phi)^2/(cos(delmp1)*(-sqrt((sin(phi+delmp1)*sin(phi))/cos(delmp1))+1)^2))-1);
kp2=1+0.5*(Rp2-1)*((cos(phi)^2/(cos(delmp2)*(-sqrt((sin(phi+delmp2)*sin(phi))/cos(delmp2))+1)^2))-1);
fup1=matlabFunction(kp1*y1*cos(delmp1));
Final_result_p1=gma*(x-h)^2*integral(fup1,0,1);
M21=matlabFunction(kp1*cos(delmp1)*y1);
Final_result_m21=gma*h1*(x-h)^2*integral(M21,0,1);
M22=matlabFunction(kp1*cos(delmp1)*y1^2);
Final_result_m22=gma*(x-h)^3*integral(M22,0,1);
Final_result_m2=Final_result_m21+Final_result_m22;
Hfup1=matlabFunction(kp1*y1*cos(delmp1));
HFinal_result_p1=gma*h1^2*integral(Hfup1,0,1);
T31=matlabFunction(kp2*cos(delmp2));
HFinal_result_T31=gma*(h+h1)*(x-h-h1)*integral(T31,0,1);
T32=matlabFunction(ka2*cos(delma2));
HFinal_result_T32=gma*h1*(x-h-h1)*integral(T32,0,1);
T33=matlabFunction(kp2*cos(delmp2)*y2-ka2*cos(delma2)*y2);
HFinal_result_T33=gma*(x-h-h1)^2*integral(T33,0,1);
T3=HFinal_result_T31+HFinal_result_T32+HFinal_result_T33;
M21=matlabFunction(kp1*cos(delmp1)*y1);
Result_M21=gma*h1*(x-h)^2*integral(M21,0,1);
M22=matlabFunction(kp1*cos(delmp1)*y1^2);
Result_M22=gma*(x-h)^3*integral(M22,0,1);
M23=matlabFunction(kp1*y1*cos(delmp1));
Result_M23=gma*(x-h-h1)*(x-h)^2*integral(M23,0,1);
M31=matlabFunction(kp2*cos(delmp2)*y2);
Result_M31=gma*(h+h1)*0.5*(x-h-h1)^2*integral(M31,0,1);
M32=matlabFunction(ka2*cos(delma2)*y2);
Result_M32=gma*h1*0.5*(x-h-h1)^2*integral(M32,0,1);
M3=Result_M31-Result_M32;
M4=matlabFunction((kp2*cos(delmp2)-ka2*cos(delma2))*y2*(1-y2));
Result_M4=gma*(x-h-h1)^3*integral(M4,0,1);
MT1=-0.5*ka1*gma*x^2*cos(delma1);
MM1=-(1/6)*ka1*gma*x^3*cos(delma1);
MT2=-0.5*ka1*gma*x^2*cos(delma1);
MM2=-(1/6)*ka1*gma*x^3*cos(delma1);
i=0;
for x=0:0.02:L
i=i+1;
if(x<h)
SF(i)=MT1;
SF1(i)=0;
BM(i)=MM1;
BM1(i)=0;
elseif(x>=h && x<(h+h1))
SF(i)=MT2+Final_result_p1;
SF1(i)=0;
BM(i)=MM2+Final_result_m2;
BM1(i)=0;
else
SF(i)=-0.5*ka1*gma*(h+h1)^2*cos(delma1) + HFinal_result_p1 - T3;
SF1(i)=0;
BM(i)=-0.5*ka1*gma*(h+h1)^2*(((h+h1)/3)+(x-h-h1))*cos(delma1)+ M3 - Result_M4; %0.5*kp1*gma*h1^2*((h1/3)+(x-h-h1))*cos(delmp1)-0.5*gma*(x-h-h1)^2*(kp2*(h+h1)*cos(delmp2)-ka2*h1*cos(delma2))-(1/6)*(kp2*cos(delmp2)-ka2*cos(delma2))*(x-h-h1)^3;
BM1(i)=0;
end
end
x=0:0.02:L;
subplot(2,1,1);
fplot(SF,x,SF1,x)
xlabel('Length of the beam in m')
ylabel('Shear Force in KN')
title('Shear force diagram')
col_header={'x',SF};
xlswrite('data.xlsx',[x(:),SF(:)],'Sheet1','A2');
xlswrite('data.xlsx',col_header,'Sheet1','A1');
subplot(2,1,2)
fplot(BM,x,BM1,x)
xlabel('Length of the beam in m')
ylabel('Bending Moment in KN-m')
title('Bending Moment diagram')
col_header={'x',BM};
xlswrite('data.xlsx',[x(:),BM(:)],'Sheet1','A2');
xlswrite('data.xlsx',col_header,'Sheet1','A1');
1 Comment
Accepted Answer
More Answers (0)
See Also
Categories
Find more on MuPAD 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!