A never used variable's index out of bounds! From the function quadl? weird, dont know why.
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
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 Kommentare
Akzeptierte Antwort
Walter Roberson
am 7 Okt. 2011
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 Kommentare
Weitere Antworten (0)
Siehe auch
Kategorien
Mehr zu Matrix Indexing finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!