% 二分法求零点
clear
k=160/(11.8*8.85*12);
w0=sqrt(0.534/k);
l=50*w0;
H=1;
w=@(y0,x)sqrt(0.534/k)-0.183*H*(l*x-y0);
m=@(y0,x)exp((k/0.026)*(w0^2-w(y0,x)^2));
f=@(y0)integral(@(x) m(y0,x),0,1, 'ArrayValued', true)-1; % intergral 积分区间只能
是[0,1]，所以要先换元
if H>=0
y0=l+1;
else
y0=-1;
end
sumold=f(y0);
dy=0.1;toly=1e-5;
while (dy>toly)&&(abs(sumold)>1e-5)
if H>=0
y0=y0-dy;
else
y0=y0+dy;
end
sumnew=f(y0);
if sumold*sumnew<=0
if H>=0
y0=y0+dy;
else
y0=y0-dy;
end
