1、1,MATLAB与数学建模,2008年7月,2,第一章 微积分、线性代数的基础知识 及其在MATLAB中的实现,1.1 数学建模中常用的微积分知识在MATLAB中的实现,MATLAB中的基本符号运算函数： x=sym(x) 创建一个符号变量x syms x y z 创建多个符号变量x，y，z r=collect(S,v) 合并同类项，S是符号表达式，v是变量或表达 式，r是合并同类项后的结果,3,例1. 按不同的方式合并同幂项EXPR= sym(x2+x*exp(-t)+1)*(x+exp(-t); expr1=collect(EXPR)expr2=collect(EXPR,exp(-t),4
2、. factor(S) 符号计算的因式分解，S是待分解的符号多项式,例2. factor指令的使用 （1）除x外不含其他自由变量的情况 (2) 含其他自由变量的情况之一 syms a x; f2=x2-a2;f1=x4-5*x3+5*x2+5*x-6; factor(f2)factor(f1) （3）对正整数的质数分解factor(1025),4,5. expand(S) 对符号多项式或函数S进行展开,例：syms x y;expand(x+1)3)expand(sin(x+y),6. r=simple(S) 或 r=simplify(S) 对符号表达式S进行化简,例3. 简化 （1）运用si
3、mplify简化syms x;f=(1/x3+6/x2+12/x+8)(1/3);sfy1=simplify(f), sfy1=simplify(sfy1) （2）运用simple简化g1=simple(f), g2=simple(g1),说明：simple给出的简化式比simplify给出的更短,5,例4. 简化 ff=cos(x)+sqrt(-sin(x)2);ssfy1=simplify(ff), ssfy2=simplify(ssfy1)gg1=simple(ff), gg2=simple(gg1),7. subs(S,old,new) 把符号变量中的变量old用new代替，new可以
4、是一个符号，也可以是具体的数,例5. 用简单的算例演示subs的置换规则 （1）产生符号函数 (4) 双精度数值置换syms a x; f3=subs(f,a,x,2, pi/3)f=a*sin(x)+5; (2) 符号变量置换 （5）数值数组置换之一f1=subs(f,sin(x),sym(y) f4=subs(subs(f,a,2), x,0: pi/6:pi)(3) 符号常数置换 （6）数值数组置换之二f2=subs(f,a,x,2,sym(pi/3) f5=subs(f,a,x,0:6,0: pi/6:pi),6,1.1.1 导数、极值和积分、Taylor公式及MATLAB中的实现,1
5、. 极限运算,注意：在左右极限不相等或左右极限有一个不存在时，MATLAB的默认状态是求右极限。,7,例8. 求极限 与极限syms x;y1=(1+4*x)(1/x);y2=(exp(x)-1) /x;limit(y1,x,0)limit(y2,x,0),例9. 求极限syms x;y=sqrt(x)-2(-1/x);limit(y,x,0, right),2. 求导运算 （1）一元函数的求导diff(f) diff(f,n),例10. 求函数 的二阶导数syms x;f=3*x3+5*x+1;diff(f,2),例11. 设 ，求syms x;y=3*x2-2*x+1;B=diff(y),
6、x=1;eval(B),8,（2）多元函数的偏导数diff(f,xi) diff(f,xi,n),例12. 求 关于x的偏导数。syms x y;z=x2*sin(2*y);B=diff(z,x),3. 积分运算 （1）一元函数的不定积分int(f) 求函数f对默认变量的不定积分，用于函数只有一个变量的情况int(f,v) 求符号函数f对变量v的不定积分,例13. 计算syms x ; y=1/(sin(x)2*cos(x)2);int(y)pretty(int(y),例14. 计算syms x z; B=int(x/(1+z2),z),9,（2）一元函数的定积分int(f,x,a,b) 用微
7、积分基本公式计算定积分,例15. 求syms x ;y=(x2+sin(x)/(1+x2);int(y,x,-1,1),（3）多重积分运算int(int(f,y),x) 计算不定积分int(int(f,y,c,d),x,a,b) 计算不定积分,例16. 计算syms x y ;int(int(x2+y2+1,y,x,x+1),x,0,1),10,4. 函数的Taylor展开taylor(f) 将函数f展开成默认变量的6阶麦克劳林(Maclaurin)公式taylor(f,n) 将函数f展开成默认变量的n阶麦克劳林(Maclaurin)公式taylor(f,n,v,a) 将函数f(v)在v=a处
8、展开成n阶Taylor公式,例17. 将函数 展开为x的6阶麦克劳林(Maclaurin)公式syms x ;f=x*atan(x)-log(sqrt(1+x2);taylor(f),例18. 将函数 展开为关于(x-2)的最高次为4的幂级数syms x ;f=1/x2;taylor(f,4,x,2);pretty(taylor(f,4,x,2),11,1.1.2 数值微分与数值积分在MATLAB中的实现,数值微分数值微分是用离散的方法近似地计算函数y=f(x)在某点x=a处的导数值，通常仅当函数以离散数值形式给出时才有这种必要。diff(x),2. 用数值方法计算定积分(1)复合梯形公式 t
9、rapz(x,y) (2)复合辛普生公式 quad(fun,a,b,tol,trace)quadl(fun,a,b,tol,trace),例19. 用复合梯形公式和复合辛普生公式求syms x ; x=2:0.1:5; y=log(x)./(x.2); t=trapz(x,y);ff=inline(log(x)./(x.2) , x ); q=quad(ff,2,5);disp(blanks(3) 梯形法求积分 blanks(3) 辛普生法求积分),t,q,12,(3)用数值方法计算二重积分 dblquad(fun,xmin,xmax,ymin,ymax)计算二重积分,例20. 计算 ，其中D
10、是由y=1, x=4, x=0, y=0所围成的区域 Q = dblquad(inline(x*y), 0, 4, 0, 1),1.1.3 线性方程和非线性方程在MATLAB中的各种求解方法,求多项式方程的根roots(p)solve(p),例21. 求方程 的所有根p=1 -4 9 -10r=roots(p)或 s1=sym(x3-4*x2+9*x-10);solve(s1),13,2. 求超越方程的根 例22. 求方程 的根，其中p，r为常数ff=sym(p*sin(x)=r)solve(ff,x),solve命令还可以用来求解方程组 例23. 求方程组 的根x,y= solve(x+y=
11、1,x-11*y=5,x,y),1.1.4 MATLAB中求和及求极值方法,求和 (1) 向量或矩阵求和sum(x) 求向量x的和或者是矩阵每一列向量的和,例24. a=1:5; A=1 2 3;2 3 4;7 8 9;sum(a)sum(A),14,(2) 级数求和symsum(s,v,a,b) 对表达式s的符号变量v从v=a到v=b进行求和,例25. 求syms k n;f=k3;symsum(f,k,0,n-1),例26. 求syms x k;symsum(xk/sym(k!),k,0,inf),2. 求函数的极值点 (1) 求一元函数的极值问题fminbnd(fun,x1,x2) 在区
12、间x1,x2内求函数fun的极小值点,例27. f = inline(sin(x)+3);x = fminbnd(f,2,5),(2) 求多元函数的极值问题fminsearch(fun,x0) 用单纯形法求多元函数fun在x0附近的极值点fminunc(fun,x0) 用拟牛顿法求多元函数fun在x0附近的极值点,15,1.1.5 函数插值与曲线的拟合,函数插值 (1) 一维插值Y1=interp1(x,y,X1, method)功能：根据已知的数据(x,y)，用method方法进行插值，然后计算X1对应的函数值Y1.,例28. 对 ， 用11个节点作三种插值，并比较其结果.x0=-5:0.5
13、:5; y0=1./(1+x0.2);x=-5:0.2:5;y1=interp1(x0,y0,x,linear);y2=interp1(x0,y0,x,spline);y3=interp1(x0,y0,x,nearst);subplot(2,2,1),plot(x0,y0,r-p),title(y=1/(1+x2);subplot(2,2,2),plot(x0,y0,r-,x,y1),title(linear);subplot(2,2,3),plot(x0,y0,r-,x,y2),title(spline);subplot(2,2,4),plot(x0,y0,r-,x,y3),title(ne
14、arst);,16,(1) 二维插值Z1=interp2(x,y,z,X1,Y1, method)功能：根据已知的数据(x,y,z)，用method方法进行插值，然后计算(X1,Y1)对应的函数值Z1.,2. 曲线的拟合 (1) 多项式拟合p,s=polyfit(x,y,n)功能：对于已知的数据组x,y进行多项式拟合，拟合的多项式的阶数为n，其中p为多项式的系数矩阵，s为预测误差估计值的矩阵。,例29. x取0至1之间的数，间隔为0.1；y为2.3，2.5， 2.1，2.5，3.2，3.6，3.0，3.1，4.1，5.1，3.8，分别用二次、三次、和七次拟合曲线来拟合这组数据，试观察三组拟合曲
15、线哪个效果最好。clf;x=0:0.1:1;y=2.3,2.5, 2.1,2.5,3.2,3.6,3.0,3.1,4.1,5.1,3.8;p2=polyfit(x,y,2);p3=polyfit(x,y,3);p7=polyfit(x,y,7);,17,disp(二次拟合曲线),p2 disp(三次拟合曲线),p3 disp(七次拟合曲线),p7 x1=0:0.1:1; y2=polyval(p2,x1); y3=polyval(p3,x1); y7=polyval(p7,x1); plot(x,y,rp,x1,y2,-,x1,y3,k-,x1,y7) legend(拟合点, 二次拟合, 三次
16、拟合, 七次拟合),(2) 非线性最小二乘拟合leastsq(f,x0) f是M函数文件,例30. 用表中的一组数据拟合 中系数r，k，并给出图像。,18,建立函数文件ct.mfunction y=ct(x)t=0.25 0.5 1 1.5 2 3 4 6 8c=19.21 18.15 15.36 14.10 12.98 9.32 7.45 5.24 3.01y=c-x(1)*exp(-x(2).*t) (2) 建立命令文件x0=10,0.5x=leastsq(ct,x0)tt=0:0.2:8yy=x(1).*exp(-x(2).*tt)plot(tt,yy,-rp),19,1.2 数学建模中
17、常用的线性代数基础知识在MATLAB中的实现,1.2.1 MATLAB中向量和矩阵的基本运算,设A、B为两个矩阵 A+B，A-B k*A A*B AB 左除 A-1B，A必须为方阵 A/B 右除 AB-1 ，B必须为方阵 det(A) 求|A|，A必须为方阵 inv(A)或A-1 An A或transpose(A) rank(A) rref(A) 矩阵行变化化简，求矩阵A阶梯形的行最简形式,20,例31. 计算,例32. 求矩 阵 的转置、逆和行列式。 syms a b c d;A=a b;c d;A inv(A) det(A),例33. 求矩 阵 的秩与行最简形。,1.2.2矩阵的变换与分解
18、及其在MATLAB中的实现,矩阵的对角元素 (1) 函数diag可以将一个矩阵的对角元素提取出来diag(A) 由矩阵A的对角线元素得到一个列向量,例34. A=pascal(3)diag(A),21,(2) 用该函数来产生第k阶对角线上的元素diag(A,k) 其中，k=0表示主对角线；k0表示在主对角线以上； k0表示在主对角线以下。,例35. a=1 2 3A=diag(a,0), B=diag(a,1), C=diag(a,-1),(3) 函数blkdiag可以根据输入创建一个分块对角矩阵S=blkdiag(a,b,c,) 功能：根据输入的a,b,c等参数构造一个分块对角矩阵,例36.
19、 a=1;b=1 2;3 4;c=9; d=1 -1;2 -2;3 -3S=blkdiag(a,b,c,d),2. 矩阵的分解 矩阵的奇异值分解 s=svd(A)U,S,V=svd(A)U,S,V=svd(A,0) (2) 矩阵的LU分解 L,U=lu(A) (2) 矩阵的QR分解 Q,R=qr(A) (2) 矩阵的Cholesky分解 C=chol(A),22,1.2.3 MATLAB中矩阵特征值和特征向量的求解方法,矩阵的特征值与特征向量d=eig(A) 返回方阵A的全部特征值组成的特征向量dV,D= eig(A) 返回矩阵A的特征值D与特征向量矩阵V，满足AV=VDpoly(A) 求矩阵
20、A的特征多项式,例37. 求矩 阵 的特征多项式、特征值和特征向量。 A=4 6 0;-3 -5 0;-3 -6 1P=poly(A)V,D= eig(A),注：求矩阵A的特征根 D=roots(poly(A) 。,2. 矩阵的相似对角化,例38. 化上述矩 阵A为对角阵。 D=V-1AV,23,例39. 求一个正交变换将二次型化为标准型A=1 1 0 -1;1 1 -1 0;0 -1 1 1;-1 0 1 1P,D= eig(A)syms x1 x2 x3 x4;X=x1;x2;x3;x4Y=P*X,1.2.4 线性方程组的直接求解法在MATLAB中的实现,齐次线性方程组的直接求解法,例40
21、. 用基础解系表示齐次方程组 的通解。A=1 1 1 1 1;3 2 1 1 -3;0 1 2 2 6;5 4 3 3 -1; B=null(A,r);syms k1 k2 k3; x=k1*B(:,1)+k2*B(:,2)+k3*B(:,3),24,2. 非齐次线性方程组的直接求解法 （1） 求逆法对AX=B，若|A|不等于0，则解由inv(A)*B获得 （2） 左除法与右除法,例41. 求非齐次线性方程组 的通解。A=2 3 1 4;1 -2 4 -5;3 8 -2 13; 4 -1 9 -6; rref(A),3. 符号方程组的求解方法 （1） 线性方程组AX=B的符号解X=lineso
22、lve(A,B) 只给出特解,例42. 求 的解。syms a b c;A=a 0 0;0 b 0; B=1;cX=linesolve(A,B),25,（2） 非线性方程组的符号解x1,x2,x3,=solve(e1,e2,e3,),例43 解e1=sym(a+b+x=y);e2=sym(2*a*x-b*y=-1);e3=sym(2*(a+b)=x+y);e4=sym(a*y+b*x=4);a,b,x,y=solve(e1,e2,e3,e4),26,第二章 常微分方程的符号解,r=dsolve(eqn1,eqn2,.,var)r=dsolve(eqn1,eqn2,. ,var) 功能：求微分方
23、程或微分方程组的通解,r=dsolve(eqn1,eqn2,., cond1, cond2,.,var)r=dsolve(eqn1,eqn2,., cond1,cond2,., var) 功能：求微分方程(组)满足初始条件的特解,在Matlab中，约定D1表示一次微分，D2表示二次微分，依此类推，Dn表示n次微分，符号Dy相当于Dy/Dt，默认自变量为t。,例1 求解一阶微分方程 及x=0时的特解。dsolve(Dy=1+y2, x) dsolve(Dy=1+y2, y(0)=1 ,x),例2 求解二阶微分方程dsolve(x2*D2y+x*Dy+(x2-1/2)*y=0,y(pi/2)=2,
24、Dy(pi/2)=-2/pi,x),27,第三章 概率论基础知识及其在MATLAB中的实现,3.1 随机事件及其概率,rand(N) 返回一个NN的随机矩阵rand(N,M) 返回一个NM的随机矩阵rand(P1,P2,Pn) 返回一个P1P2Pn的随机矩阵round(X) 对向量或矩阵X的每个分量四舍五入取整fix(x) 对x朝零方向取整floor(x) 对x朝负无穷大方向取整ceil(x) 对x朝正无穷大方向取整,二项分布binopdf(X,N,P) 计算二项分布的密度函数binocdf(X,N,P) 计算二项分布的累积分布函数binoinv(X,N,P) 计算二项分布的逆累积分布函数bi
25、nornd(N,P,m,n) 产生服从二项分布的mn阶随机矩阵binostat(N,P) 求二项分布的数学期望与方差,3.2 随机变量的分布及其数学特征,28,2. 泊松分布poisspdf(X,LMD)poisscdf(X,LMD)poissinv(Y,LMD)poissrnd(LMD,m,n)poisstat(LMD),3. 超几何分布hygepdf(M,n,k,N)hygecdf(M,n,k,N)hygeinv(P,n,k,N)hygernd(n,k,N,mr,mc)hygestat(n,k,N),4. 均匀分布unifpdf(X,A,B)unifcdf(X,A,B)unifinv(P,
26、A,B)unifrnd(A,B,m,n)unifstat(A,B),5. 指数分布exppdf(X,L)expcdf(X,L)expinv(P,L)exprnd(X,L,m,n)expstat(L),6. 正态分布normpdf(X,M,C)normcdf(X,M,C)norminv(P,M,C)normrnd(M,C,m,n)normstat(M,C),7. 分布chi2pdf(X,N)chi2cdf(X,N)chi2inv(P,N)chi2rnd(N,m,n)chi2stat(N),29,8. t分布tpdf(X,N)tcdf(X,N)tinv(P,N)trnd(N,m,n)tstat(N
27、),9. F分布fpdf(X,M,N)fcdf(X,M,N)finv(P,M,N)frnd(M,N,m,n)fstat(M,N),3.3 样本的数字特征,1. 集中趋势（位置）度量,geomean(X) 样本的几何平均值harmmean(X) 样本的调和平均值mean(X) 样本的算术平均值median(X) 样本的中位数trimmean(X,P) 修正的样本均值，计算时剔除样本数据中最高P和最低P后的均值,注：上述命令中，若X为向量，则返回该向量的各度量值；若X为矩阵，则返回矩阵各列的度量值组成的行向量。,30,2. 离中趋势（散布）度量,vqr(X) 样本的内四分位数间距mad(X) 样本
28、的平均绝对差rang(X) 样本的极差var(X) 样本的方差std(X) 样本的标准差cov(X)或cov(X,Y) 样本的协方差,3. 中心矩 moment(X,k),4. 峰度 kurtosis(X)偏度 skewness(X),5. 相关系数 corrcoef(X),31,第四章 线性规划及其在MATLAB中的实现,X=lp(c,A,b)X=lp(c,A,b,v1)X=lp(c,A,b,v1,v2)X=lp(c,A,b,v1,v2,x0)X=lp(c,A,b,v1,v2,x0,ne)X=lp(c,A,b,v1,v2,x0,ne,dis)X,lag=lp(c,A,b,)X,lag,how=lp(c,A,b,),例1 求解两种方法：(1) c=1,1 (2) c=1,1 a=1,-1 a=1,-1;-1,0 b=1 b=1,0v1=0 X=lp(c,a,b)X=lp(c,a,b,v1),
