1、5.1 创建函数例:x=cos(at)+b y=|x|+c创建一个comxy函数计算,文件名为comxy.mfunctionx,y=comxy(a,t,b,c)x=cos(a*t)+b y=abs(x)+c 大写字母开头5.1.1 函数文件第五章 函数15.1 创建函数例:x=c o s(a t)+b y=|x|+c2functionx,y=comxy(t,a,b,c)%变量名可不一致%x=sin(at)+bx=cos(a*t)+b;y=abs(x)+c;n=3;b=linspace(1,1.4,n);for k=1:3u,v=comxy(0:pi/4:pi,sqrt(1.8/(1+k)3),
2、b(k),1/0.85)end调用程序2 f u n c t i o n x,y =c o m x y(t,a,b,c)334functionx,y=comxy(t)%x=sin(at)+bglobal A B Cx=cos(A*t)+B;y=abs(x)+C;global A B Cn=3;C=1/.85;b=linspace(1,1.4,n);for k=1:n A=sqrt(1.8/(1+k)3);B=b(k);m,n=comxy(0:pi/4:pi)end全局变量的用法,大写,空格(不用逗号)4 f u n c t i o n x,y =c o m x y(t)g l o b a l5
3、functionx,y=comxy(t,a,b,c)%x=sin(at)+bx=cos(a*t)+b;y=abs(x)+c;clear alln=3;b=linspace(1,1.4,n);t=0:pi/4:pi;p=zeros(n,length(t);q=zeros(n,length(t);for k=1:3 p(k,:),q(k,:)=comxy(t,sqrt(1.8/(1+k)3),b(k),1/.85)end5 f u n c t i o n x,y =c o m x y(t,a,b,c)65.2.2 inlinefox=inline(x.2.*cos(a*x)-b,x,a,b)b=f
4、ox(pi/3 pi/3.5,4,1)5.2.3 子函数6 5.2.2 i n l i n e f o x=i n l i n e(x.2.7785.3 用户自定义函数、函数句柄和feval函数句柄是一种引用函数的方法。函数句柄由和函数名构成。feval就是把已知的数据或符号带入到一个定义好的函数句柄中。clcclear allsyms t%定义符号变量f=(x,y)x3+y2k1=feval(f,2,t)%k2=f(2,t)k3=feval(f,1,2)%k4=f(1,2)8 5.3 用户自定义函数、函数句柄和f e v a l 函数句柄是一种95.4 以数组变量为输入参数的MATLAB函数
5、5.4.1 多项式拟合-polyfit/polyval函数polyfit()采用最小二乘法对给定数据进行多项式拟合,其具体使用方法如下:c=polyfit(x,y,n),采用n次多项式p来拟合数据x和y。y(x)=c1xn+c2xn-1+cnx+cn+1 (5.1)求得c后,可以用下列语句求出式5.1y=polyval(c,xnew)9 5.4 以数组变量为输入参数的MA T L A B 函数5.4.1 101 011function nd=datand=50,0.13;70,0.092;90,0.072;110,0.057;130,0.046;150,0.037;170,0.028;190,
6、0.020;210,0.015;230,0.010;250,0.007;ncs=data;c=polyfit(ncs(:,1),ncs(:,2),4);r=input(enter r(0r0.2);su=input(enter su(50su p=1 0 3 12-7%f(x)=x4+3x2+12x-7p=1 0 3 12 -7 roots(p)ans=0.7876+2.4351i 0.7876-2.4351i -2.0872 0.5121 2 2 在MA T L A B 7 语言里,多项式由一个行向量表示,设为p23如果已知某个多项式的根,那么,使用poly函数,可以很轻松地产生其对应的多项
7、式。r=roots(1,-10,35,-50,24)c=poly(r)使用conv函数对多项式进行乘法运算。格式为c=conv(a,b),其中a和b为两个多项式的系数向量,c为相乘所生成的多项式的系数向量。a=1 2 3 4;b=5 6 7 8;c=conv(a,b)2 3 如果已知某个多项式的根,那么,使用p o l y 函数,可以很轻248.1.2 数值积分的实现方法1变步长辛普生法基于变步长辛普生法,MATLAB给出了quad函数来求定积分。该函数的调用格式为:I,n=quad(fname,a,b,tol,trace)其中fname是被积函数名。a和b分别是定积分的下限和上限。tol用来
8、控制积分精度,缺省时取tol=0.001。trace控制是否展现积分过程,若取非0则展现积分过程,取0则不展现,缺省时取trace=0。返回参数I即定积分值,n为被积函数的调用次数。2 4 8.1.2 数值积分的实现方法25求定积分。(1)建立被积函数文件fesin.m。function f=fesin(x)f=exp(-0.5*x).*sin(x+pi/6);(2)调用数值积分函数quad求定积分。S,n=quad(fesin,0,3*pi)%n可省略S=0.9008n=772 5 求定积分。26二重定积分的数值求解使用MATLAB提供的dblquad函数可以求二重定积分的数值解。该函数的调
9、用格式为:I=dblquad(f,a,b,c,d,tol,trace)该函数求f(x,y)在a,bc,d区域上的二重定积分。参数tol,trace的用法与函数quad完全相同。2 6 二重定积分的数值求解27计算二重定积分(1)建立一个函数文件fxy.m:function f=fxy(x,y)global Ki;Ki=Ki+1;%ki用于统计被积函数的调用次数f=exp(-x.2/2).*sin(x.2+y);(2)调用dblquad函数求解。global Ki;Ki=0;I=dblquad(fxy,-2,2,-1,1)Ki2 7 计算二重定积分28常微分方程常微分方程(Ordinary di
10、fferential equations,ODE)n初值问题初值问题-给出初始值给出初始值n边值问题边值问题-给出边界条件给出边界条件与初值常微分方程解算有关的指令与初值常微分方程解算有关的指令ode23 ode45 ode113 ode23tode15s ode23sode23tb2 8 常微分方程(O r d i n a r y d i f f e r e n t i a l 29n常用调用格式:常用调用格式:t,y=ode45(odefun,tspan,y0)t,y=ode45(odefun,tspan,y0)u 参数说明:参数说明:uodefun odefun 表示表示f(t,y)f(
11、t,y)的的函数句柄或函数句柄或inlineinline函数函数,t t是标量,是标量,y y是是标量或向量;标量或向量;utspan tspan 如果是二维向量如果是二维向量t0,tft0,tf,表示自变量初值,表示自变量初值t0t0和和tftf;如;如果是高维向量果是高维向量t0,t1,t0,t1,tn,tn,则表示输出结点列向量;,则表示输出结点列向量;uy0 y0 表示初值向量表示初值向量y0y0;ut t 表示结点列向量表示结点列向量(t0,t1,(t0,t1,tn)T;,tn)T;uy y 表示数值解矩阵,每一列对应表示数值解矩阵,每一列对应y y的一个分量;的一个分量;n完整调用
12、格式:完整调用格式:t,y=ode45(odefun,tspan,y0t,y=ode45(odefun,tspan,y0,options,p1,p2,options,p1,p2,)uoptions options 为计算参数设置(如精度要求等),默认为计算参数设置(如精度要求等),默认表示;表示;up1p1,p2p2,附加传递参数,附加传递参数,odefunodefun表示为表示为f(t,y,p1,p2,f(t,y,p1,p2,)2 9 常用调用格式:t,y =o d e 4 5(o d e f u n,t s30一一.解解ODE的基本机理的基本机理:2.把高阶方程转换成一阶微分方程组把高阶方
13、程转换成一阶微分方程组1.列出微分方程列出微分方程初始条件初始条件令令)0()0(210yyYyyyy21,),(),(),(2121YtfYtfYtfyyY(5.2)(5.3)(5.4),(tyyfy 00)0(,)0(yyyy3 0 一.解O D E 的基本机理:2.把高阶方程转换成一阶微分31例:著名的例:著名的Van der Pol方程方程0)1(2yyyy 令令 yyyy21,降为一阶降为一阶1221221)1(yyyyyyY21yyY初始条件初始条件2010210)0()0(yyyyY3 1 例:著名的V a n d e r P o l 方程令 降为一阶初始条323.根据式根据式(
14、5.3)编写计算导数的编写计算导数的M函数文件函数文件-ODE文件文件把把t,Y作为输入宗量,把作为输入宗量,把 作为输出宗量作为输出宗量Y%M function file name:dYdt.m function Yd=f(t,Y)Yd=f(t,Y)的展开式的展开式例例Van der Pol方程方程%M function file name:dYdt.m function Yd=f(t,Y)Yd=zeros(size(Y);Yd(1)=Y(2);Yd(2)=-(Y(1).2-1)*Y(2)-Y(1);3 2 3.根据式(5.3)编写计算导数的M函数文件-O D E 文334.使编写好的使编写
15、好的ODE函数文件和初值函数文件和初值 供微分方程解算指令供微分方程解算指令(solver)调用)调用0YSolver解算指令的使用格式解算指令的使用格式t,Y=solver(ODE函数文件名函数文件名,t0,tN,Y0,tol);ode45)()()()()()(222112110201tytytytytytyY输出宗量形式输出宗量形式)()1(:,1tyY)()2(:,2tyY说明:说明:t0:初始时刻;:初始时刻;tN:终点时刻:终点时刻Y0:初值;:初值;tol:计算精度:计算精度3 3 4.使编写好的O D E 函数文件和初值 供微分方34例题例题1:著名的著名的Van der Po
16、l方程方程0)1(2yyyy%主程序主程序 (程序名:程序名:VanderPol _ex1.m)t0=0;tN=20;tol=1e-6;Y0=0.25;0.0;t,Y=ode45(f,t0,tN,Y0,tol);subplot(121),plot(t,Y)subplot(122),plot(Y(:,1),Y(:,2)解法:采用解法:采用ODE命令命令3 4 例题1:著名的V a n d e r P o l 方程%主程序 35Van der Pol方程方程%子程序子程序 (程序名:程序名:dYdt.m)function Ydot=dYdt(t,Y)Ydot=Y(2);-Y(2)*(Y(1)2-1
17、)-Y(1);或写为或写为function Ydot=dYdt(t,Y)Ydot=zeros(size(Y);Ydot(1)=Y(2);Ydot(2)=-Y(2)*(Y(1).2-1)-Y(1);3 5 V a n d e r P o l 方程%子程序 (程序名:363 637n常微分方程一阶方程组边值问题常微分方程一阶方程组边值问题MATLABMATLAB标准形式:标准形式:n调用格式调用格式uSinit=bvpinit(tinit,yinit)Sinit=bvpinit(tinit,yinit)由在粗略结点由在粗略结点tinittinit的预估解的预估解yinityinit生成粗略解网格生
18、成粗略解网格sinitsinituSol=bvp4c(odefun,bcfun,sinit)Sol=bvp4c(odefun,bcfun,sinit)odefun odefun是微分方程组函数,是微分方程组函数,bcfunbcfun为边值条件函数,为边值条件函数,sol.xsol.x为求为求解结点,解结点,sol.ysol.y是是y(t)y(t)的数值解的数值解uSx=deval(sol,ti)Sx=deval(sol,ti)计算由计算由bvp4cbvp4c得到的解在得到的解在titi的值的值btabyaygtytfty,0)(),(),(,()(3 7 常微分方程一阶方程组边值问题MA T
19、L A B 标准形式:38例题:求解边值问题例题:求解边值问题解:首先改写成方程组解:首先改写成方程组 边界条件为边界条件为2)4(,0)0(,0 zzzz)1()2(),2()1(yyyy02)1(,0)1(ybya3 8 例题:求解边值问题39求解用求解用M M函数函数eg6_5fun.meg6_5fun.m%M%M函数函数eg6_5fun.meg6_5fun.mclear;close;clear;close;sinit=bvpinit(0:4,1;0)sinit=bvpinit(0:4,1;0)odefun=inline(odefun=inline(y(2);-abs(y(1)y(2);
20、-abs(y(1),t t,y y););bcfun=inline(bcfun=inline(ya(1);yb(1)+2ya(1);yb(1)+2,yaya,ybyb););sol=bvp4c(odefun,bcfun,sinit)sol=bvp4c(odefun,bcfun,sinit)t=linspace(0,4,101);t=linspace(0,4,101);y=deval(sol,t);y=deval(sol,t);plot(t,y(1,:),sol.x,sol.y(1,:),plot(t,y(1,:),sol.x,sol.y(1,:),o o,sinit.x,sinit.x,sin
21、it.y(1,:),sinit.y(1,:),s s)legend(legend(解曲线解曲线,解点解点,粗略解粗略解)3 9 求解用M函数e g 6 _ 5 f u n.m40最优化是求最优解,也就是在某个区间内有条件约束或者无条件约束地找到函数的最大值或者最小值。MATLAB使用数字方法求函数的最小值。使用迭代算法,也就是有些步骤要重复许多次。现在,假设要求函数f在某个区间内的最小值xmin。迭代方法需要一个初始估计值x0。从x0开始找到一个更接近xmin的新值x1,这个值的好坏取决于使用的数学方法。直到找到有足够精度的近似值xi才停止迭代,也就是绝对值|xminxi|足够小。这里提到了标
22、准MATLAB系统的两个最优化命令,fmin命令可以求单变量函数的最小值;fmins命令可以求多变量函数的最小值,同时它还要求有一个初始向量。在新版本中,fmin和fmins分别被fminbnd和fminsearch取代。Matlab中没有求函数f的最大值的命令,但可以通过求其相反函数h=f 的最小值间接求得。4 0 最优化是求最优解,也就是在某个区间内有条件约束或者无条件41函数格式:x=fminbnd(fcn,x1,x2,options)求函数在区间(x1,x2)内取最小值时的x值,采用黄金分割法和二次插值。fcn是目标函数名。如果没有局部最小值,则返回区间内的最小x值。向量options
23、为控制参数,如options(1)=1,显示中间结果;options(2)表示得到的结果x的精度,缺省为104,options可省略。x=fmins(fcn,x0,options)求函数fcn的最小值。由用户自己给出一个初始估计向量x0,优化参数options可省略。4 1 函数格式:42例:在区间 0,2 内求函数f(x)=x3-2x-5的最小值。f=inline(x.3-2*x-5);x=fminbnd(f,0,2)x=0.8165f(x)ans=-6.08874 2 例:在区间 0,2 内求函数f(x)=x 3-2 x-5符号数学工具箱是操作和解决符号表达式的符号数学工具箱(函数)集合,
24、有复合、简化、微分、积分以及求解代数方程和微分方程的工具。另外还有一些用于线性代数的工具,求解逆、行列式、正则型式的精确结果,找出符号矩阵的特征值而无由数值计算引入的误差。符号数学工具箱中的工具是建立在功能强大的称作Maple软件的基础上。它最初是由加拿大的滑铁卢(Waterloo)大学开发的。当要求MATLAB进行符号运算时,它就请求Maple去计算并将结果返回到MATLAB命令窗口。因此,在MATLAB中的符号运算是MATLAB处理数字的自然扩展。符号数学工具箱是操作和解决符号表达式的符号数学工具箱(函数)44定义符号变量定义符号变量参与符号运算的对象可以是符号变量、符号表达式或符号矩阵。
25、符号变量要先定义,后引用。1.sym函数:每次定义一个符号变量。x=sym(x)%创建变量x。x=sym(x,real)%创建实变量x。x=sym(x,unreal)去掉x的实数属性。a=sym(alpha)%创建变量a,打印为alpha。2.syms函数sysms arg1 arg2 argN。一次可定义多个符号变量,使用起来更方便。注意各变量间用空格隔开。4 4 定义符号变量45符号表达式符号表达式符号表达式的书写格式与数值表达式一样,由符号变量、函数、运算符组成,如:syms a b c xf=a*x2+b*x+c 也可以用f=sym(a*x2+b*x+c)把符号表达式ax2+bx+c赋
26、给f,但此时没有定义变量a b c x。有两种方法创建符号函数,一是用符号表达式,如上面f=a*x2+b*x+c 定义了符号函数f,从而可以用diff等符号函数对f进行运算。如diff(f)将返回结果2*a*x+b。4 5 符号表达式46极限极限函数 limit格式 limit(F,x,a)%计算符号表达式F=F(x)当xa时的极限值。limit(F,a)%用命令findsym(F)确定F中的自变量,设为变量x,再计算F的极限值,当xa时。limit(F)%用命令findsym(F)确定F中的自变量,设为变量x,再计算F的极限值,当x0时。limit(F,x,a,right)或limit(F,
27、x,a,left)%计算符号函数F的单侧极限:左极限xa-或右极限xa+。例:syms x a t h;Then limit(sin(x)/x)=1limit(1/x,x,0,right)=inflimit(sin(x+h)-sin(x)/h,h,0)=cos(x)4 6 极限47导数(包括偏导数)导数(包括偏导数)函数 diff格式 diff(S,v)、diff(S,sym(v)%对表达式S中指定符号变量v计算S的1阶导数。diff(S)%对表达式S中的符号变量v计算S的1阶导数,其中v=findsym(S)。diff(S,n)%对表达式S中的符号变量v计算S的n阶导数,其中v=findsy
28、m(S)。diff(S,v,n)%对表达式S中指定的符号变量v计算S的n阶导数。例:syms x tdiff(sin(x2)%结果为 2*cos(x2)*xdiff(t6,t,6)%结果为 returns 7204 7 导数(包括偏导数)48符号函数的积分符号函数的积分函数 int格式 R=int(S,v)%对符号表达式S中指定的符号变量v计算不定积分。注意的是,表达式R只是函数S的一个原函数,后面没有带任意常数C。R=int(S)%对符号表达式S中的符号变量v计算不定积分,其中v=findsym(S)。R=int(S,v,a,b)%对表达式s中指定的符号变量v计算从a到b的定积分R=int(S,a,b)%对符号表达式s中的符号变量v计算从a到b的定积分,其中v=findsym(S)。例:syms x nint(xn)%或int(xn,x)、int(xn,x),结果均为 x(n+1)/(n+1)int(sin(2*x),0,pi/2)%结果为 14 8 符号函数的积分494 9505 0515 1