《第6章MATLAB解方程与最优化问题求解》由会员分享，可在线阅读，更多相关《第6章MATLAB解方程与最优化问题求解（20页珍藏版）》请在金锄头文库上搜索。
1、第 6 章 MATLAB 解方程与最优化问题求解 在科学研究和工程应用的许多领域，很多问题都常常归结为解方程，包括线性 方程、非线性方程和常微分方程。研究方程的解析解固然能使人们更好地掌握问题 的规律。但是，在很多情况下无法求出其解析解，这时数值解法就是一个十分重要 的手段。MATLAB 为解决这类问题提供了极大的方便。 另外，在日常生活和实际工作中，人们对于同一个问题往往会提出多个解决方 案，并通过各方面的论证从中提取最佳方案。最优化方法就是专门研究如何从多个 方案中科学、合理地提取出最佳方案的科学。 本章介绍有关解方程和最优化问题的 MATLAB 实现方法。 MATLAB 线性方程组求解
2、MATLAB 非线性方程数值求解 MATLAB常微分方程初值问题的数值解法 MATLAB 最优化问题求解 6.1 线性方程组求解 将包含 n 个未知数，由 n 个方程构成的线性方程组表示为 11 11221nn121 12222nn2n1 1n22nnnna xa xa xb a xaxaxba xaxaxb 其矩阵表示形式为 Ax=b 其中 11121n1121222n22n1n2nnnnaaaxb aaaxbAxbaaaxb ，。， 第 6 章 MATLAB 解方程与最优化问题求解 137 在 MATLAB 中，关于线性方程组的解法一般可以分为两类：一类是直接解法，就是在没 有舍入误差的情
3、况下，通过有限步的矩阵初等运算来求得方程组的解；另一类是迭代解法， 就是先给定一个解的初始值，然后按照一定的迭代算法进行逐步逼近，求出更精确的近似解。 6.1.1 直接解法 1利用左除运算符的直接解法 线性方程组的直接解法大多基于高斯消元法、主元素消元法、平方根法和追赶法等。在 MATLAB 中，这些算法已经被编成了现成的库函数或运算符，因此，只需调用相应的函数或 运算符即可完成线性方程组的求解。其中，最简单的方法就是使用左除运算符“” 。程序会 自动根据输入的系数矩阵判断选用哪种方法进行求解。 对于线性方程组 Ax=b，可以利用左除运算符“”求解： x=Ab 当系数矩阵 A 为 NN 的方阵
4、时，MATLAB 会自行用高斯消元法求解线性方程组。若右 端项 b 为 N1 的列向量，则 x=Ab 可获得方程组的数值解 x（N1 的列向量） ；若右端项 b 为 NM 的矩阵，则 x=Ab 可同时获得与系数矩阵 A 相同的 M 个线性方程组的数值解 x（为 NM 的矩阵） ，即 x(:,j)=Ab(:,j)，j=1,2,M。注意，如果矩阵 A 是奇异的或接近奇异的，则 MATLAB 会给出警告信息。 当系数矩阵 A 不是方阵时，称这样的方程组为欠定方程组或超定方程组，MATLAB 将会 在最小二乘意义下求解。 例 6-1 用直接解法求解下列线性方程组： 123412423412342xx5
5、xx13 x5x7x92xxx6 x6xx4x0 解：程序如下： A=2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4; b=13,-9,6,0; x=Ab x = -66.5556 25.6667 -18.7778 26.5556 在 2.4.3 节中曾介绍过利用矩阵求逆来解线性方程组， 即 x=A-1b， 其结果与使用左除运算相同。 2利用矩阵的分解求解线性方程组 矩阵分解是指根据一定的原理用某种算法将一个矩阵分解成若干个矩阵的乘积。常见的 矩阵分解有 LU 分解、QR 分解、Cholesky 分解、Schur 分解、Hessenberg 分解及奇异分解等。 这里着
6、重介绍前 3 种常见的分解。通过这些分解方法求解线性方程组的优点是运算速度快， 可以节省存储空间。 MATLAB 程序设计教程（第二版） 138 （1）LU 分解。矩阵的 LU 分解就是将一个矩阵表示为一个交换下三角矩阵和一个上三角 矩阵的乘积形式。线性代数中已经证明，只要方阵 A是非奇异的，LU 分解总是可以进行的。 MATLAB 提供的 lu 函数用于对矩阵进行 LU 分解，其调用格式为： L,U=lu(X)：产生一个上三角矩阵 U 和一个变换形式的下三角矩阵 L（行交换） ，使 之满足 X=LU。注意，这里的矩阵 X 必须是方阵。 L,U,P=lu(X)：产生一个上三角矩阵 U 和一个下
7、三角矩阵 L 以及一个置换矩阵 P，使 之满足 PX=LU。当然矩阵 X 同样必须是方阵。 当使用第种格式时，矩阵 L 往往不是一个下三角矩阵，但可以通过行交换使其成为一 个下三角矩阵。设 111 A543211 则对矩阵 A 进行 LU 分解的命令如下： A=1,-1,1;5,-4,3;2,1,1 A = 1 -1 1 5 -4 3 2 1 1 L,U=lu(A) L = 0.2000 -0.0769 1.0000 1.0000 0 0 0.4000 1.0000 0 U = 5.0000 -4.0000 3.0000 0 2.6000 -0.2000 0 0 0.3846 为检验结果是否正
8、确，输入命令： LU=L*U LU = 1 -1 1 5 -4 3 2 1 1 说明结果是正确的。例中所获得的矩阵 L 并不是一个下三角矩阵，但经过各行互换后， 即可获得一个下三角矩阵。 利用第种格式对矩阵 A 进行 LU 分解的命令如下： L,U,P=lu(A) L = 1.0000 0 0 0.4000 1.0000 0 0.2000 -0.0769 1.0000 U = 5.0000 -4.0000 3.0000 第 6 章 MATLAB 解方程与最优化问题求解 139 0 2.6000 -0.2000 0 0 0.3846 P = 0 1 0 0 0 1 1 0 0 LU=L*U %这
9、种分解其乘积不为 A LU = 5 -4 3 2 1 1 1 -1 1 inv(P)*L*U %考虑矩阵 P 后其乘积等于 A ans = 1 -1 1 5 -4 3 2 1 1 实现 LU 分解后，线性方程组 Ax=b 的解 x=U(Lb)或 x=U(LPb)，这样可以大大提高运 算速度。 例 6-2 用 LU 分解求解例 6-1 中的线性方程组。 命令如下： A=2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4; b=13,-9,6,0; L,U=lu(A); x=U(Lb) x = -66.5556 25.6667 -18.7778 26.5556 或采用 LU
10、分解的第种格式，命令如下： L,U,P=lu(A); x=U(LP*b) 将得到与上面同样的结果。 （2）QR 分解。对矩阵 X 进行 QR 分解，就是把 X 分解为一个正交矩阵 Q 和一个上三角 矩阵 R 的乘积形式。QR 分解只能对方阵进行。MATLAB 的函数 qr 可用于对矩阵进行 QR 分 解，其调用格式为： Q,R=qr(X)：产生一个正交矩阵 Q 和一个上三角矩阵 R，使之满足 X=QR。 Q,R,E=qr(X)：产生一个正交矩阵 Q、一个上三角矩阵 R 及一个置换矩阵 E，使之满 足 XE=QR。 设 111 A5432710 则对矩阵 A 进行 QR 分解的命令如下： MAT
11、LAB 程序设计教程（第二版） 140 A=1,-1,1;5,-4,3;2,7,10; Q,R=qr(A) Q = -0.1826 -0.0956 -0.9785 -0.9129 -0.3532 0.2048 -0.3651 0.9307 -0.0228 R = -5.4772 1.2780 -6.5727 0 8.0229 8.1517 0 0 -0.5917 为检验结果是否正确，输入命令： QR=Q*R QR = 1.0000 -1.0000 1.0000 5.0000 -4.0000 3.0000 2.0000 7.0000 10.0000 说明结果是正确的。利用第种格式对矩阵 A 进行
12、 QR 分解： Q,R,E=qr(A) Q = -0.0953 -0.2514 -0.9632 -0.2860 -0.9199 0.2684 -0.9535 0.3011 0.0158 R = -10.4881 -5.4347 -3.4325 0 6.0385 -4.2485 0 0 0.4105 E = 0 0 1 0 1 0 1 0 0 Q*R/E %验证 A=Q*R*inv(E) ans = 1.0000 -1.0000 1.0000 5.0000 -4.0000 3.0000 2.0000 7.0000 10.0000 实现 QR 分解后，线性方程组 Ax=b 的解 x=R(Qb)或
13、x=E(R(Qb)。 例 6-3 用 QR 分解求解例 6-1 中的线性方程组。 命令如下： A=2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4; b=13,-9,6,0; Q,R=qr(A); x=R(Qb) x = -66.5556 25.6667 第 6 章 MATLAB 解方程与最优化问题求解 141 -18.7778 26.5556 或采用 QR 分解的第种格式，命令如下： Q,R,E=qr(A); x=E*(R(Qb) 将得到与上面同样的结果。 （3）Cholesky 分解。如果矩阵 X 是对称正定的，则 Cholesky 分解将矩阵 X 分解成一个 下三
14、角矩阵和一个上三角矩阵的乘积。 设上三角矩阵为 R， 则下三角矩阵为其转置， 即 X=RR。 MATLAB 函数 chol(X)用于对矩阵 X 进行 Cholesky 分解，其调用格式为： R=chol(X)：产生一个上三角矩阵 R，使 RR=X。若 X 为非对称正定，则输出一个 出错信息。 R,p=chol(X)：这个命令格式将不输出出错信息。若 X 为对称正定的，则 p=0，R 与 上述格式得到的结果相同； 否则p为一个正整数。 如果X为满秩矩阵， 则R为一个阶数为q=p-1 的上三角矩阵，且满足 RR=X(1:q,1:q)。 设 211 A121113 则对矩阵 A 进行 Cholesky 分解的命令如下： A=2,1,1;1,2,-1;1,-1,3; R=chol(A) R = 1.4142 0.7071 0.7071
