A never used variable's index out of bounds! From the function quadl? weird, dont know why.
1 回表示 (過去 30 日間)
古いコメントを表示
I run this code and it gave me error message, and it seems that this problem comes from the function quadl, because the variable y never appear in my code, y is from the function quadl.I set some breakpoint,still couldnt find why.If anyone can help fix this problem , i appreciate it very much. anyone can help find what the problem is. I am really annoyed by this problem.thank you so much!
here is the error:
"??? Attempted to access y(13); index out of bounds because numel(y)=1.
Error in ==> quadl at 72
if ~isfinite(y(13))
Error in ==> G_yi at 20
f(i)=quadl(@sumin,Friction1,Friction2,[],[],i);
Error in ==> quadl at 64
y = feval(f,x,varargin{:}); y = y(:).';
Error in ==> double_int at 11
SS=quadl(@G_yi,Protrusion1,Protrusion2,[],[],Friction1,Friction2);
Error in ==> main at 4
PR=double_int(0,Di,0,Di) "
here is the code, one main code and 3 functions.
Di=4e-3;
PR=double_int(0,Di,0,Di)
function SS=double_int(innlow,innhi,outlow,outhi)
Protrusion1=outlow;Protrusion2=outhi;Friction1=innlow;Friction2=innhi;
SS=quadl(@G_yi,Protrusion1,Protrusion2,[],[],Friction1,Friction2);
function f=G_yi(Protrusion,Friction1,Friction2)
Protrusion=Protrusion(:);n=length(Protrusion);
save G_yi.mat;
f=zeros(n,1);
for i=1:n
f(i)=quadl(@sumin,Friction1,Friction2,[],[],i);
end
f=f(:);
function fun=sumin(Friction,i)
load G_yi.mat;
MeanShearStress=5 ;
Density=1e3;
Ustar=sqrt(MeanShearStress/Density);
N=15;
D=zeros(N,1);
p=zeros(N,1);
for k=1:N-1;
D(1)=0.5e-3;
p(1)=1/N;
D(k+1)=D(k)+0.5e-3;
p(k+1)=p(k);
end;
KinematicViscosity=1.004e-6;
D50=4e-3;
Roughness=2*D50;
RoughnessRenolds=Ustar*Roughness/KinematicViscosity;
if RoughnessRenolds>100
Intensity=Ustar*2.14;
Skewness=0.43;
Flatness=2.88;
else
Intensity=Ustar*(-0.187*log(RoughnessRenolds)+2.93);
Skewness=0.102*log(RoughnessRenolds);
Flatness=0.136*log(RoughnessRenolds)+2.30;
end
CoefficientC=-0.993*log(RoughnessRenolds)+12.36;
SpecificWeightofSand=1.8836e4;
SpecificWeightofWater=9.789e3 ;
Di=4e-3;
Thickness=1.5*D50;
Y1=0.25*Thickness;
Y2(i)=0.25*Thickness+Protrusion(i);
Sum=zeros(N,1);
for m=1:N-1
Kapa=0.4;
ff1=@(t) (Ustar*CoefficientC.*t/Thickness).*sqrt((0.5*Di)^2-(t-Protrusion(i)-Y1+0.5*Di).^2);
ff2=@(t) sqrt((0.5*Di)^2-(t-Protrusion(i)-Y1+0.5*Di).^2);
MeanBedVelocity=zeros(n,1);
MM=zeros(n,1);
MM(i)= quadl(ff2,Y1,Y2(i));
if Y2(i)<=Thickness
MeanBedVelocity(i)=quadl(ff1,Y1,Y2(i))./(MM(i));
else
MeanBedVelocity(i)=(quadl(ff1,Y1,Thickness)+quadl(ff1,Thickness,Y2(i)))./(quadl(ff2,Y1,Y2(i)));
end
Yb=zeros(n,1);
if MeanBedVelocity(i)<=Ustar*CoefficientC
Yb(i)=(MeanBedVelocity(i)*Thickness)/(Ustar*CoefficientC);
else
Yb(i)=Thickness*exp(Kapa*(MeanBedVelocity(i)/Ustar-CoefficientC));
end;
ParticleRenolds=zeros(n,1);
ParticleRenolds(i)=MeanBedVelocity(i).*Protrusion(i)/KinematicViscosity;
Cd=zeros(n,1);
if ParticleRenolds(i)<=1754
Cd(i)=(24./ParticleRenolds(i)).*(1+0.15*ParticleRenolds(i).^0.687);
else
Cd(i)=0.36;
end;
Cl=zeros(n,1);
Cl(i)=Cd(i);
h1=zeros(n,1);
h1(i)=Yb(i)-Y1-Protrusion(i)+0.5*Di;
h2=Di*(Friction+0.5*D(m)-0.5*Di)/(Di+D(m));
Ld=h1(i)+h2;
Ll=sqrt((0.5*Di)^2-h2.^2);
Lw=Ll;
Pei=zeros(N,1);
Phi=zeros(N,1);
for j=2:N;
Pei(1)=p(1)*Di/(Di+D(1));
Pei(j)=p(j)*Di/(Di+D(j))+ Pei(j-1);
Phi(j)=1-Pei(j);
end;
HidingFactor=(Pei(N)/Phi(N))^0.6;
EffectiveShearStress=HidingFactor*MeanShearStress;
DimensionlessEffectiveShearStress= EffectiveShearStress/((SpecificWeightofSand-SpecificWeightofWater)*Di);
ff3=@(t) sqrt((0.5*Di)^2-(t-Protrusion(i)-Y1+0.5*Di).^2);
A=zeros(n,1);
A(i)=quadl(ff3,Y1,Y2(i));
Br=Ustar*sqrt((2*Lw*pi*Di^2)./((Cd(i).*Ld+Cl(i).*Ll)*6.*A(i)*DimensionlessEffectiveShearStress)); % Rolling Threshold
Bl=Ustar*sqrt((2*pi*Di^2)/(Cl(i)*6.*A(i)*DimensionlessEffectiveShearStress)); %Lifting Threshold
syms Ub ;
Pr=zeros(n,1);
PD=zeros(n,1);
PD(i)=int((exp(-((Ub-MeanBedVelocity(i))/Intensity).^2/2)/(Intensity*sqrt(2*pi))).*(1+(Skewness/factorial(3))*(((Ub-MeanBedVelocity(i))/Intensity).^3-3*((Ub-MeanBedVelocity(i))/Intensity))+(Flatness-3)*(((Ub-MeanBedVelocity(i))/Intensity).^4-6*((Ub-MeanBedVelocity(i))/Intensity).^2+3)/factorial(4))/Intensity,-Bl,-Br)...
+int((exp(-((Ub-MeanBedVelocity(i))/Intensity).^2/2)/(Intensity*sqrt(2*pi))).*(1+(Skewness/factorial(3))*(((Ub-MeanBedVelocity(i))/Intensity).^3-3*((Ub-MeanBedVelocity(i))/Intensity))+(Flatness-3)*(((Ub-MeanBedVelocity(i))/Intensity).^4-6*((Ub-MeanBedVelocity(i))/Intensity).^2+3)/factorial(4))/Intensity,Br,Bl);
Pr(i)= double(PD(i));
Sum(1)=p(1)*Pr(i);
Sum(m+1)=p(m)*Pr(i)+Sum(m);
end;
fun=Sum(m+1);
0 件のコメント
採用された回答
Walter Roberson
2011 年 10 月 7 日
The function passed to quad must take exactly one parameter. Your sumin() takes two parameters.
Does your sumin() return a vector with as many results as there were values in the single input parameter?
0 件のコメント
その他の回答 (0 件)
参考
カテゴリ
Help Center および File Exchange で Matrix Indexing についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!