高等光学例题解析

高等光学----第三章例题讲解clearcloseallbesselj1=@(x)besselj(1,x);forn=1:10z(n)=fzero(besselj1,[(n-1)n]*pi);endx=0:pi/100:10*pi;y=besselj(1,x);plot(z,zeros(1,10),'o',x,y,'-')line([010*pi],[00],'color','black')axis([010*pi-0.51.0])xlabel('U')ylabel('J_1(U)')[(1:n)'z']MATLAB程序：例3.1求解一阶贝塞尔函数前10个零点并标注，做函数曲线1.匿名函数直接调用格式besselj1=@(x)besselj(1,x);2.for循环体Forindex=[matrix](定义索引矩阵)

commandstobeexecuted（执行命令）End3.fzero查找一元连续函数的零点x=fzero(f,x0)%查找f函数在x0附近的零点4.zeros返回一个零矩阵zeros(n,m)返回一个n行m列的零矩阵5.line画线函数line([X1X2],[Y1Y2],S)，就是在点A(X1,Y1)和点B(X2Y2)之间画一条直线，S为其它属性（颜色，线的粗细等）。6.axis（x1x2y1y2）设置坐标轴范围知识点讲解：结果显示：一阶贝塞尔函数零点数值解：一阶贝塞尔函数曲线及其标注零点：functionxp=FindZeroOfBesselj(m,n)xp=zeros(n,1);x0=m+2.5;xp(1)=fzero(@(x)besselj(m,x),x0);

if(n>1.5)xp(2)=fzero(@(x)besselj(m,x),xp(1)+pi);end

if(n>2.5)forii=3:nxp(ii)=fzero(@(x)besselj(m,x),2*xp(ii-1)-xp(ii-2));endendx=linspace(0,ceil(xp(end)),1000);y=besselj(m,x);plot(xp,zeros(1,n),'o',x,y,'-')%`作出m阶贝赛尔函数曲线并标记其零点`line([0ceil(xp(end))],[00],'color','black')axis([0ceil(xp(end))-0.51.0])xlabel('U')ylabel('J(U)')title(['m='num2str(m)',n='num2str(n)])求解m阶贝塞尔函数的前n个根MATLAB程序1.Function用来定义一个函数，Xp是输出变量，(m,n)是输入变量，FindZeroOfBesselj是函数名，在MATLAB的命令窗口中调用FindZeroOfBesselj(m,n)就可得到对应解。2.根据表3.5知道m阶贝塞尔函数的第1个根在哪个数值附近，再用fzero求精确的贝塞尔函数的根，根据求的前n个根不同，所以需要确定在不同的x值附近求精确的贝塞尔函数的根3.

Linespace（x1,x2,N）用于产生x1,x2之间的N点行矩阵4.Ceil函数是在朝正无穷的方向取整知识点讲解：结果显示：10阶贝塞尔函数曲线及其前15个标注零点：10阶贝塞尔函数零点数值解：例3.2MATLAB程序：clearcloseallticVmax=10;N=100;forj=1:NV(j)=j/N*Vmax;Vtemp=V(j);Utemp=NaN;i=0;while(isnan(Utemp)&&i<N+1)

init=Vtemp*(i+1)/N-eps;

tryUtemp=fzero(@(Utemp)...besselj(0,Utemp)/(Utemp*besselj(1,Utemp))-...besselk(0,sqrt(Vtemp^2-Utemp^2))/...(sqrt(Vtemp^2-Utemp^2)*besselk(1,sqrt(Vtemp^2-Utemp^2))),init);catchendi=i+1;

end

U(j)=Utemp;end求解零LP01模的特征方程，做V-U、V-W的关系曲线W=sqrt(V.^2-U.^2);Ymax=ceil(max([U,W]));figuresubplot(1,2,1)plot(V,U,'r');axisequalaxis([0Vmax0Ymax])xlabel('V')ylabel('U')title('LP_{01}V-U')MATLAB程序：subplot(1,2,2)plot(V,W,'r');axisequalaxis([0Vmax0Ymax])xlabel('V')ylabel('W')title('LP_{01}V-W')toc知识点：1.NaN是NotaNunmber的意思，相当于将空值赋值给Utemp2.Isnan函数返回一个与A相同维数的数组，若A的元素为NaN（非数值），在对应位置上返回逻辑1（真），否则返回逻辑0（假）3.&&逻辑运算符，相当于与运算4.eps(a)是|a|与大于|a|的最小的浮点数之间的距离，距离越小表示精度越高。init=Vtemp*(i+1)/N-eps，

此处用-eps表示要将后面的小数点去掉5.try的作用是让Matlab尝试执行一些语句,执行过程中如果出错,则执行catch部分的语句try

尝试执行语句块catch

出错后执行语句块LP模的本征方程：数理知识点：V、U、W之间的关系：结果显示：程序运行的时间结果：V-U和V-W曲线LP01模对应的V、U、W的数值解例3.3根据3.2得到的结果，拟合区间的V和W线性关系式，并给出拟合得到数值最大的误差MATLAB程序:NN=15:24;x=V(NN);y=W(NN);p=polyfit(x,y,1);f=polyval(p,x);maxerr=max(y-f);figureplot(x,y,'o',x,f,'-')xlabel('V')ylabel('W')1.polyfit函数用作线性拟合polyfit（x,y,1)其中x表示源数据点对应的横坐标，y表示源数据点对应的纵坐标，1表示拟合阶数2.polyval函数用作求拟合得到的数值polyval返回n次多项式p在x处的值polyval(p,x)3.max表示返回数组中的最大值知识点：结果显示：V和W线性拟合曲线：拟合系数和最大误差：例3.4在MATLAB中作图，给出LP01模在V分别是0.8、1.6、2.4时，电场分量E相对归一化直径Ra的归一化曲线MATLAB程序：clearcloseallV=[0.80001.60002.4000];U=[0.79741.36701.6453];W=[0.06400.83151.7473];Ra1=-1:0.01:1;Ra2=[-5:0.01:-1];Ra3=[1:0.01:5];E1=zeros(length(V),length(Ra1));E2=zeros(length(V),length(Ra2));E3=zeros(length(V),length(Ra3));fori=1:length(V)E1(i,:)=besselj(0,U(i)*Ra1);E2(i,:)=besselj(0,U(i)).*besselk(0,W(i).*abs(Ra2))./besselk(0,W(i));E3(i,:)=besselj(0,U(i)).*besselk(0,W(i).*abs(Ra3))./besselk(0,W(i));endR=[Ra2Ra1Ra3];E=[E2E1E3];plot(R,E)xlabel('R_a=r/a')ylabel('E')holdonplot([-1-1],[01],'b--',[11],[0,1],'b--')知识点：1.Length求数组长度函数Length（x）求x行或者列中最大长度2.E1(I,:)

MATLAB中"："的使用非常灵活，简单的说就是代表矩阵中的全行或全列元素.E1矩阵的第i行的所有元素。数理知识点：1.LP01模在纤芯区和包层区的归一化电场分布结果显示：电场E相对归一化直径Ra的归一化曲线V=0.8V=1.6V=2.4例3.5clearcloseallV=2.4000;U=1.6453;W=1.7473;Npoint=21;Rx=2;Ry=2;x=linspace(-Rx,Rx,Npoint);y=linspace(-Ry,Ry,Npoint);X=meshgrid(x,y);Y=meshgrid(y,x);Y=Y';R=sqrt(X.^2+Y.^2);E1=besselj(0,U*R);E2=besselj(0,U).*besselk(0,W.*R)./besselk(0,W);E=E1;pos=find(R>=1);E(pos)=E2(pos);imagesc(x,y,E);colormap(gray),colorbarxlabel('x'),ylabel('y'),title('E')figurecontour(x,y,E),colorbaraxisequalxlabel('x'),ylabel('y'),title('E')figure[px,py]=gradient(E);contour(X,Y,E),holdonquiver(X,Y,px,py),holdoffaxisequalMATLAB程序：LP01模在V=2.4时，归一化电场分量E相对归一化直径Ra的二维分布图知识点：Meshgrid数据集成绘图meshgrid（x，y）的作用是分别产生以向量x为行，向量y为列的两个大小相同的矩阵2.Find函数用于返回所需要元素的所在位置find(R>=1)返回R>=1的数值位置3.imagesc(x,y,A)x,y分别为二维向量，Matlab会在[x1,x2]*[y1,y2]范围内染色。4.用colormap来定义图像显示用的颜色查找表，比如用colormap(pink)，可以把黑白图像显示成带粉红色的图像，很有趣的5.[Fx,Fy]=gradient(F)，其中Fx为其水平方向上的梯度，Fy为其垂直方向上的梯度，Fx的第一列元素为原矩阵第二列与第一列元素之差，Fx的第二列元素为原矩阵第三列与第一列元素之差除以2，以此类推：Fx(i,j)=(F(i,j+1)-F(i,j-1))/2。最后一列则为最后两列之差。同理，可以得到Fy。6.contour是等高线绘制函数contour(X,Y,Z),contour(X,Y,Z,n),contour(X,Y,Z,v)绘制Z的等高线。X,Y限制在x,y轴上的范围。7.

quiver(x,y,u,v)该函数使用箭头来直观的显示