matlab标量的矩阵次幂,matlab计算时报错错误使用 ^ 用于对矩阵求幂的维度不正确。请检查并确保矩阵为方阵并且幂为标量。要执行按元素矩阵求幂，请使用 '.^'...

matlab报错：错误使用  ^ 用于对矩阵求幂的维度不正确。请检查并确保矩阵为方阵并且幂为标量。要执行按元素矩阵求幂，请使用 '.^'

源代码是：

if t1==2 %x轴为Db

for i=1:100

Db=0.01+0.0002*i;

syms y x        %y=ubv  x=cdv

if rhob>rhol

%默认第一个条件

[x, y]=solve(x==24/129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu,y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

if Re>0.1 && Re<1000      %第二个

syms y x

[x, y]=solve(x==24/Re*(1+0.14*Re^0.7),y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif Re>=1000 && Re<=35000    %第三个

syms y x

[x, y]=solve(x==0.445,y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif  Re>=35000 && Re<2.9*10^5   %第四个

syms y x

[x, y]=solve(x==24/Re*(1+0.1806*Re^0.6459)+0.4251/(1+6880.95*Re^(-1)),y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif Re>=2.9*10^5 && Re<=3.85*10^5      %第五个

syms y x

[x, y]=solve(x==-4.16*10^(-6)*Re+1.67,y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif  Re>=3.85*10^5 && Re<=1*10^6   %第六个

syms y x

[x, y]=solve(x==7.45*10^(-8)*Re+0.0385,y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif  Re>1*10^6                     %第七个

syms y x

[x, y]=solve(x==0.19-8*10^4/Re,y==sqrt(2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

end

else

%默认第一个条件

[x, y]=solve(x==24/129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu,y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

if Re>0.1 && Re<1000      %第二个

syms y x

[x, y]=solve(x==24/Re*(1+0.14*Re^0.7),y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif Re>=1000 && Re<=35000    %第三个

syms y x

[x, y]=solve(x==0.445,y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif  Re>=35000 && Re<2.9*10^5   %第四个

syms y x

[x, y]=solve(x==24/Re*(1+0.1806*Re^0.6459)+0.4251/(1+6880.95*Re^(-1)),y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif Re>=2.9*10^5 && Re<=3.85*10^5      %第五个

syms y x

[x, y]=solve(x==-4.16*10^(-6)*Re+1.67,y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif  Re>=3.85*10^5 && Re<=1*10^6   %第六个

syms y x

[x, y]=solve(x==7.45*10^(-8)*Re+0.0385,y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

elseif  Re>1*10^6                     %第七个

syms y x

[x, y]=solve(x==0.19-8*10^4/Re,y==sqrt(-2*m*g*(1-rhol/rhob)/x/A/rhol)+ulv);

x=vpa(x(1),5);

y=vpa(y(1),5);

Re=129*abs(y-ulv)*196.848*39.3701*Db*rhol/1000/mu;

end

end

Nre(i)=double(Re);

ubv(i)=double(y);

Cdv(i)=double(x);

db(i,1)=Db;

FI(i)=(0.4536*0.3048)*3.78*10^(-4)*rhob/1000*(39.3701)^3*Db^3*(196.848*ubv)^2/39.3701/Dc; %N  %惯性力

Fd(i)=(0.4536*0.3048)*1.18*10^(-5)*Cdv*rhol/1000*(39.3701*Db)^2*Vperf^2;  %拖拽力 N

Fu(i)=0.3927*Cdv*rhol*Vl^2*(Db^2-Db^2*theta(i)/180+(Dp/pi)*(Db^2-Dp^2)^(1/2)); %脱落力N

Fh(i)=Dp/sqrt(Db^2-Dp^2)*(1.76*10^(-4)*Dp^2*rhol*Q^2*(1.062/n^2/Dp^4/alpha^2-1/Dc^4));

view=strcat('已经到计算Db=',num2str(Db));

disp(view);

end

axes(handles.axes1) ;

plot(db,FI,'b','linewidth',3);

hold on

plot(db,Fd,'r','linewidth',3);

hold on

plot(db,Fu,'g','linewidth',3);

hold on

plot(db,Fh,'y','linewidth',3);

legend('惯性力','拖拽力','脱落力','保持力','Location','Northeast');

grid on

xlabel('封堵球直径(m)');

ylabel('合力 N');

end