整理matlab16常用计算方法

1、精品文档常用计算方法1.超越方程的求解一超越方程为x(2lnx-3)-100=0求超越方程的解。算法方法一：用迭代算法。将方程改为100x=21n(%)-3其中看是一个初始值，由此计算终值值换成X0的值，重新进行计算；否则程序P1_1abs.m如下。%超越方程的迭代算法clearx0=30;xx=;while1x=100/(2*log(x0)-3);xx=xx,x;iflength(xx)>1000,break,endifabs(x0-x)<1e-4,break,endx0=x;endfigurex。取最大误差为e=10-4,当|x-Xo|>e时，就用x的|x-xo|<

2、e为止。%青除变量沏始值%空向量炕限循环%迭代运算%连接结果皈口果项数太多则退出循环(暗示发散)他精度足够高时退出循环%替换初值%吉束循环喔建图形窗口gridonfs=16;title(超越方程的迭代折线','fontsize'xlabel('itn','fontsize',fs)ylabel('itx','fontsize',fs)plot(xx,'.-','LineWidth'连接对口网格%字体大小,fs)麻题%标签%而签text(length(xx),xx(end),n

3、um2str(xx(end),图示用下标作为自变量画迭代的折线。要迭代19次才能达到精度，超越方程的解为,2,'MarkerSize',12)%画迭彳线.-'表示每个点用.来表示，再用线'fontsize',fs)%显示结果如P0_20_1图所示，当最大误差为10-4时，需27.539。算法方法二：用求零函数和求解函数。将方程改为函数f(x)-2ln(x)-3-100xMATLAB求零函数为fzero,fzero函数的格式之一是x=fzero(f,x0)其中，f表示求解的函数文件，x0是估计值。fzero函数的格式之二是x=fzero(f,x1,x2)精

4、品文档精品文档solve,计算非线性方程和方程组的符号解。,16),16)其中，x1和x2表示零点的范围。另外MATLAB还有求解函数程序P1_2fzero.m如下。力超越方程的求法clearx=10:0.1:100;f=inline('2*log(x)-3-100./x'figureplot(x,f(x),'LineWidth',2)gridonx0=fzero(f,20,30);%x0=fzero(f,20);holdonplot(x0,f(x0),'.')title('超越方程的解','fontsize'xl

5、abel('itx','fontsize',16)ylabel('itf','fontsize',16)text(x0,0,num2str(x0),'fontsize'x0=solve('2*log(x)-3-100./x'plot(double(x0),0,'o')%青除变量%!变量向量%定义内线函数用的是字符窜喔建图形窗口%画曲线对口网格%求方程的零点求方程的零点%呆持图像%画零点%题%标签%而签%标记零点%求超越方程的符号解%再画零点（double是把字符转化成数字）P12.导

6、数的计算正弦函数y=sinx的导数是余弦函数y1=cosx,余弦函数的导数是负的正弦函数，用MATLAB的数值导数和符号导数求正弦函数的一阶和二阶导数，并与其解析解进行比较。程序P2diff.m如下。%正弦函数导数的计算方法cleardx=0.01*2*pi;x=0:dx:2*pi;y=sin(x);f1=diff(y)/dx;f1=f1(1),(f1(1:end-1)+f1(2:end)/2,f1(end);figure%青除变量%间隔%!变量向量须函数%通过差分求导数%求平均值喔建图形窗口精品文档精品文档f2=diff(f1)/dx;f2=f2(1),(f2(1:end-1)+f2(2:e

7、nd)/2,f2(end);d2y_dx2=diff(y,2);df2=subs(d2y_dx2,sx,x);figureplot(x,-sin(x),x,f2,'.',x,df2,'o')gridon%通过差分求导数%求平均值%求二阶符号导数%符号替换数值喔建图形窗口%画二阶导数和差分以及符号导数曲线%!口网格plot(x,cos(x),x,f1,'.')%画一阶导数和数值差分曲线%plot(x,cos(x),x(1:end-1),f1,'.')%数值导数（点）偏左%plot(x,cos(x),x(2:end),f1,'

8、.')%数值导数（点）偏右symssx%定义符号变量y=sin(sx);%建立符号函数dy_dx=diff(y);%求符号导数df1=subs(dy_dx,sx,x);%符号替换数值holdon%呆持图像plot(x,df1,'ro')%画符号导数曲线gridon对口网格legend（'解析导数，数值差分','符号导数',4）洞例title（'正弦函数的一阶导数','FontSize',16）出口标题legend（'解析导数，数值差分，符号导数'，4）%图例title（'正弦函数的二

9、阶导数，'FontSize',16）%J口标题图示（1）如P2a图所示，正弦函数的一阶导数的数值解（点）与解析解（线）符合得很好。（2）如P2b图所示，正弦函数的二阶导数的数值解（点）和符号解（圈）与解析解（线）符合得很好，不过二阶数值导数在端点与精确值有一点偏离。P2a图图P2b3.积分的计算求证：函数y=eaxsinbx的积分为_1axS=-22e（asinbx-bcosbx）Cab其中a=-0.5,b=2。积分下限为0。上限为x,画出定积分的函数曲线。证明利用分部积分得精品文档精品文档ax1ax1axaxS=esinbxdx=sinbxde=esinbx-becosbxd

10、xaa1 axbax1axbaxax=esinbx-cosbxde=esinbx-ecosbxbesinbxdxaaaa.一1b2一即S=e(asinbx-bcosbx)2Saa由此可证不定积分。当x=0时,S应该为零，所以C-ba2b2因此，从0开始的积分为ax，e(asinbx-bcosbxb)利用复数积分的方法更简单。由于eaxibxdx=eax1bxC=eax1bx.C-aiba2b2其中C'表示复常数。根据欧拉公式eix=cosx+isinx,上式两边取虚部即可证明同一结果。上式两边取实部还可证明eaxcosbxdx=-2eax(bsinbxacosbx)Cab算法设被积函数

11、为y=f(x),取间隔为限,取上限为x=nAx,则积分可用求和公式近似表示nS='f(xi)lxi=1积分既能用上式近似计算，也能用积分的解析式计算，还能用数值积分和符号积分计算。程序P3quad.m如下。力数值积分和符号积分方法cleara=-0.5;b=2;dx=0.1;xm=6;%青除变量%旨数的常数%正弦函数的常数%间隔%±限x=0:dx:xm;油变量向量s1=(exp(a*x).*(-b*cos(b*x)+a*sin(b*x)+b)/(aA2+bA2);“积分的解析解y=exp(a*x).*sin(b*x);s2=cumtrapz(y)*dx;figureplot(

12、x,s1,x,s2,'.')gridons='exp(',num2str(a),'*x).*sin('f=inline(s);s3=0;fori=2:length(x)%被积函数%梯形法积分喔建图形窗口%画积分曲线对口网格,num2str(b),'*x)'%积分函数字符串%化为内线函数，才可以被调用(画成)施1个积分值%按自变量循环精品文档精品文档s3=s3,quad(f,0,x(i);endholdonplot(x,s3,"or')symssasbsxss=exp(sa*sx)*sin(sb*sx);sy=in

13、t(ss,sx)ssy=subs(sy,sa,sb,a,b);%连接积分quad又tf积分下限0,上限x(i)%吉束循环%呆持图像%画数值积分曲线%定义符号变量%被积符号函数%tsx进行符号积分%替换常数s4=subs(ssy,sx,x);%替换向量因为sx与sa,sb的长度不一样，不能同时替代plot(x,s4-s4(1),"ko',"MarkerSize',10)%符号积分曲线tit='ityrm=eA",num2str(a),"itAxrmsin",num2str(b),"itx"谁成数学公式：

14、人表示上标title(tit,"rm的积分","FontSize',16)琳题legend("公式法","梯形法","数值法","符号法",4)位置，-1把图例放到外面)图示如P3图所示，梯形法积分(点)与积分的解析解(线)符合得很好，4.微分方程的求解方法(1)求一阶微分方程的解dy_2ydxx1当x=0时，y=2,这是初始条件。用微分方程的数值解和符号解画出函数曲线，并与解析解进行比较。对口图例(4表示右下角，0表示电脑选择最佳4412dy2x22dxx1dxP3图(2)求

15、二阶微分方程的解初始条件为y(0)=1,y'(0)=2。用微分方程的数值解和符号解画出函数曲线和导数的曲线,并与解析解进行比较。解析(1)分离变量得dy2dxyx1积分得lny=2ln(x+1)+C利用初始条件可得C=ln2,因此y=2(x+1)2程序P4_1ode.m如下。%L阶常微分方程的解析解，数值解和符号解clear%青除变量x=linspace(0,2,50);%变量向量y1=2*(x+1).A2;%单析解精品文档精品文档f=inline('2*y/(x+1)');x2,y2=ode45(f,x,2);ys=dsolve('Dy-2*y/(x+1)&#

16、39;定)Dy=dy_dxy3=subs(ys,'x',x);%微分方程右边化为内线函数%求微分方程的数值解（,'y（0）=2','x'）%求微分方程的符号特解ode常微分方程）dsolve（初始条件决figureplot(x,y1,x,y2,'.',x,y3,'o')gridonlegend('解析解','数值解','符号解',4)xlabel('itx','FontSize',16)ylabel('ity','

17、;FontSize',16)title('一阶常微分方程的解','FontSize'%1符号改为向量求数值解喔建图形窗口%画曲线对口网格洞例%黄坐标%纵坐标,16）麻题（圈）与解析解（线）符合得图P42dyFy2xdxx21图示如P4_1图所示，一阶微分方程的数值解（点）和符号解很好。P4_1图解析（2）由于y'=dy/dx,分离变量得积分得lny'-ln(x2+1)=Ci当x=0时，y'=2,所以Ci=ln2,因此y'=2(x2+1)再积分O2Qy=(2x22)dx=-x32xC23当x=0时，y=1,所以C2=1,因此

18、y=-x32x13设y(1)=y,y(2)=dy/dx,可得两个一阶微分方程dy(1)_、2dy(2)_2xy(2)y(2)2ddxdxx1将两个一阶微分方程设计成函数文件，以便求数值解。程序P0_23_2ode.m如下。精品文档精品文档%二阶常微分方程的解析解，数值解和符号解clearx=linspace(0,3,30);y1=1+2*x+2*x.A3/3;dy1=2*x.A2+2;x2,Y=ode45('p4_2fun',x,1,2);认了在(0,0)点附初值y2=Y(：,1);dy2=Y(:,2);ys=dsolve('D2y-2*x*Dy/(xA2+1)'

19、;条件决定)D2y标对y求二阶倒数dy=diff(ys);y3=subs(ys,'x',x);dy3=subs(dy,'x',x);figuresubplot(2,1,1)plot(x,y1,x2,y2,'.',x,y3,'o'：gridon%青除变量%!变量向量%华析解姗析解的导数麻微分方程的数值解1,21,表示第一个初始条件默%取出函数%取出导数1'表示第一列所有元素取出,'y(0)=1','Dy(0)=2','x')%t微分方程的符号特解(由初始'x'指

20、出求得自变量是xlegend('解析解','数值解','符号解',2)xlabel('itx','FontSize',16)ylabel('ity','FontSize',16)title('二阶常微分方程解的函数','FontSize'subplot(2,1,2)plot(x,dy1,x2,dy2,gridon',x,dy3,'o')legend('解析解','数值解','符号解'

21、;,2)title('二阶常微分方程解的导数','FontSize'xlabel('itx','FontSize',16)ylabel('dity/ditx','FontSize',16)%符号导数%等符号改为向量求函数的数值解%等符号改为向量求导数的数值解喔建图形窗口%选子图%画曲线对口网格洞例%黄坐标呦坐标,16)麻题%选子图%画曲线对口网格洞例,16)麻题%黄坐标%纵坐标程序在求微分方程的数值解时将调用函数文件%二阶常微分方程的数值解的函数文件P42fun.mfunctionf=fun(x,y

22、)f=y(2);2*x*y(2)/(xA2+1);一阶导数表达式二阶导数表达式图示如P4_2图所示，二阶微分方程的数值解符合得很好。(包括函数和导数)和符号解与解析解都5.偏导数的计算和等量异号点电荷的电场两个异号点电荷带电量为iQ(Q>0),相距为解析如B5图所示，等量异号点电荷在场点2a,画出电场线和等势线。P(x,y)产生的电势为kQkQri其中，k为静电力常量，1和2是场点P到电荷的距离精品文档B5图精品文档电场强度可根据电势梯度计算其中，劈形算符为在xy平面上，场强只有两个分量ri”(xa)y,r2/(x-a)y；:xEx.:U:xEy.:U两个点电荷在P点产生的电场强度的大小

23、分别为kQE1=2E2kQ场强的两个分量也能根据公式计算riEx=EicosZ-E2cosZkQ(xa)kQ(x-a)3riEy=Esin%。E2sin%kQykQy算法取a为坐标单位，做无量纲处理。则电势可表示为3ri其中，Uo=kQ/a。Uo是Q在原点产生的电势，作为电势的单位。*.ri和r2是约化距离*q772:2*ri=J(x+i)+y,2(2)(4)(6)(7a)(7b)(i*)(2*)Ex其中，x=x/a,y=y/a。x和y是无量纲的坐标或约化坐标。场强的x分量用梯度可表示为a：:(x/a)*-u即Ex=-E0(5a*):x其中，Eo=Uo/a,U*=U/Uo。E0是场强的单位，U

24、:是无量纲的电势。同理可得*LLUEy=-Eo-(5b：)二y两个点电荷的电场强度的两个分量用公式可表示为精品文档精品文档*1*1Ex=E0(H，Ey=E0(ri2*3一*3ri口(7*)将物理量无量纲化之后，只要作纯数值计算就行了。MATLAB的梯度函数gradient可直接计算场强的数值分量，场强的数值解和解析解可相互比较。等势线可根据等值线指令contour绘制，电场线可根据流线指令streamline绘制。程序P0_24gradient.m如下。%等量异号点电荷的电场线和等势线（请在“创建图形窗口”处设置断点，以观察画图过程clearxm=2.5;ym=2;x=linspace(-xm

25、,xm,400);y=linspace(-ym,ym,400);X,Y=meshgrid(x,y);R1=sqrt(X+1).A2+Y.A2);R2=sqrt(X-1).A2+Y.A2);U=1./R1-1./R2;u=-4:0.5:4;figureC=contour(X,Y,U,u,'LineWidth',2);clabel(C,'fontsize',16)gridonholdonplot(-xm;xm,0;0,0;0,-ym;ym)%青除变量苏黄坐标范围%纵坐标范围苏黄坐标向量%纵坐标向量%坐标网点（矩阵）%左边第一个正电荷到场点的距离%右边第二个负电荷到场点的距离%计算电势%等势线的电势向量喔建图形窗口%画等势线并取等势线的坐标%标记等势线的值对口网格%呆持图像%画水平和竖直线plot(-1,0,'o',1,0,'o','MarkerSize',12)%0电荷Ex,Ey=gradient(-U,