theta=0;
R=[];
for ws=0:0.5:400;
y=dscf_foo(1,ws/3500,2.7,1,0.25,theta,3);
q=trapz(ws,y);
R=[R,q];
end
function d=dscf_foo(u1,alpha1,rho1,rho2,v,theta,gama)
r=1; %半径
b1=sqrt((2-2*v)/(1-2*v))*alpha1;
alpha2=alpha1*gama;
x1=alpha1*r;
x2=alpha2*r;
y1=b1*r;
N=rho2/rho1;
ds=0;
for n=0:10
if n==1
e=1;
else
e=2;
end
E1131=(n^2+n-0.5.*y1.^2).*besselh(n,x1)-x1.*besselh(n-1,x1);
E1111=(n^2+n-0.5.*y1.^2).*besselj(n,x1)-x1.*besselj(n-1,x1);
E1231=n*(-(n+1)*besselh(n,y1)+y1.*besselh(n-1,y1));
E1112=(N/2)*y1.^2.*besselj(n,x2);
E4131=-n*(-(n+1).*besselh(n,x1)+x1.*besselh(n-1,x1));
E4111=-n*(-(n+1).*besselj(n,x1)+x1.*besselj(n-1,x1));
E4231=-(n^2+n-0.5*y1.^2).*besselh(n,y1)+y1.*besselh(n-1,y1);
E7131=x1.*besselh(n-1,x1)-n*besselh(n,x1);
E7112=x2.*besselj(n-1,x2)-n*besselj(n,x2);
E7231=n*besselh(n,y1);
E7111=x1.*besselj(n-1,x1)-n*besselj(n,x1);
E2111=-(n^2+n-x1.^2+0.5*y1.^2).*besselj(n,x1)+x1.*besselj(n-1,x1);
E2131=-(n^2+n-x1.^2+0.5*y1.^2).*besselh(n,x1)+x1.*besselh(n-1,x1);
E2231=n*((n+1)*besselh(n,y1)-y1.*besselh(n-1,y1));
A=[E1131 E1231 E1112;E4131 E4231 0;E7131 E7231 -E7112];
B1=[E1111;E4111;E7111];
B=-e*1i^n.*B1;
X=lsqminnorm(A,B);
An=X(1,:);
Bn=X(2,:);
% Cn=X(3,:);
D=2*u1/r^2*(e*1i^n.*E2111+An.*E2131+Bn.*E2231).*cos(n.*theta);
ds=ds+D;
end
d=ds/(-u1.*b1^2);
end
错误使用 trapz (line 47) 维度参数必须是处于索引范围内的正整数标量。
问问大家,这个问题怎么解决呢,上面是报错的程序,下面是调用的function