精通MATLAB科学计算(第3版)(王正林)09-3r.pdf
《精通MATLAB科学计算(第3版)(王正林)09-3r.pdf》由会员分享,可在线阅读,更多相关《精通MATLAB科学计算(第3版)(王正林)09-3r.pdf(39页珍藏版)》请在得力文库 - 分享文档赚钱的网站上搜索。
1、第 章 线 性 方 程 组 求 解 在 自 然 科 学 和 工 程 技 术 中,很 多 问 题 可 以 归 结 为 求 解 线 性 方 程 组。采 用 MATLAB,不 仅 可 以 利 用 其 提 供 的 相 关 函 数 直 接 解 决 一 些 简 单 的 线 性 方 程 组,而 且 可 以 通 过 简 洁 的 编 程 来 解 决 一 些 复 杂 的 线 性 方 程 组。通 过 本 章 的 介 绍,读 者 既 能 应 用 MATLAB中 相 应 的 函 数 求 解 线 性 方 程 组,又 能 通 过 编 程,灵 活 使 用 迭 代 法 和 其 他 的 特 殊 解 法 来 求 解 线 性 方 程
2、 组。9.1 求 逆 法 MATLAB中 求 解 线 性 方 程 最 直 接 的 方 法 是 矩 阵 求 逆 法,它 适 用 于 系 数 矩 阵 的 数 据 无 规 律 且 系 数 矩 阵 的 阶 数 比 较 小 的 情 况。如 果 系 数 矩 阵 的 阶 数 太 大 的 话,系 数 矩 阵 求 逆 就 需 要 花 费 很 长 时 间。在 MATLAB中 用 求 逆 法 求 解,线 性 方 程 可 以 直 接 用 左 除 法“”求 解,也 可 以 用 求 逆 函 数 inv()求 解。【例 9-1】左 除 法 和 求 逆 法 求 解 线 性 方 程 组 实 例。用 左 除 法 和 求 逆 法
3、求 解 如 下 线 性 方 程 组,X+2%2+3x3=1-x+3工 2+7工 3=49x i+3x3=7解:用 左 除 法 求 解,在 MATLAB命 令 窗 口 中 输 入 求 解 程 序:A=l 2 3;-1 3 7;9 0 3 b=l 4 7 z;x=A b输 出 计 算 结 果 为:X=0.3 3 3 3-1.6 6 6 71.3 3 3 3用 逆 矩 阵 法 求 解,在 MATLAB命 令 窗 口 中 输 入 求 解 程 序:第 9 章 线 性 方 程 组 求 解 A=l 2 3;-1 3 7;9 0 3;b=l 4 7 z;x=in v(A)*b输 出 计 算 结 果 为:X=0
4、.3 3 3 3-1.6 6 6 71.3 3 3 3由 计 算 结 果 可 知,方 程 组 的 解 为 的 4 2户 3=03333,-1.6667,1.3333。9.2 分 解 法 矩 阵 分 解 法 是 指 根 据 一 定 的 原 理 用 某 种 算 法 将 系 数 矩 阵 分 解 成 若 干 个 矩 阵 的 代 数 运 算,常 用 的 分 解 是 乘 积 分 解,而 分 解 后 的 新 矩 阵 一 般 是 某 种 特 殊 矩 阵。常 用 的 分 解 有 L U分 解、Q R分 解、Cholesky分 解、Schur分 解、Hessenberg分 解 和 奇 异 分 解 等。9.2.1
5、LU分 解 法 矩 阵 的 L U分 解 也 叫 三 角 分 解,还 有 的 书 称 之 为 Doolittle分 解,它 是 将 矩 阵 分 解 为 一 个 单 位 下 三 角 矩 阵 与 上 三 角 矩 阵 的 乘 积。只 要 矩 阵 非 奇 异,这 种 分 解 总 是 可 以 进 行 的。MATLAB提 供 了 lu函 数 来 求 矩 阵 A 的 L U分 解,常 用 格 式 为:也,S=lu(Z):产 生 一 个 上 三 角 矩 阵 U 和 一 个 变 换 形 式 的 下 三 角 矩 阵 Z(进 行 了 行 交 换),使 得 A=LU;,(/,P=lu(/1):产 生 一 个 上 三
6、角 矩 阵 U 和 一 个 下 三 角 矩 阵 乙 及 置 换 矩 阵 P,使 得 PA=LU;分 解 以 后,方 程 组=的 解 可 写 为 x=LA(ZA)或 x=(A(ZAP*)。【例 9-21 L U分 解 法 求 解 线 性 方 程 组 实 例。用 L U分 解 法 求 解 以 下 线 性 方 程 组,1.5xi+3x2-0.83+414=42xi+9%3+10%4=0 lx+4.8x2 0.6x3+=11 4%|+12.3%2 4M+5X 4=2 解:用 L U分 解 法 求 解,在 MATLAB命 令 窗 口 中 输 入 求 解 程 序:A=1.5 3-0.8 4;2 0 9 1
7、 0;-7 4.8-0.6 1;1 4 1 2.3-4 5;b=4 0 1-2;Lz U=lu(A);1 187精 通 M ATLAB科 学 计 算(第 2 忡-x=U(L b)输 出 计 算 结 果 为:x=-0.3 7 2 1-0.8 2 9 1-1.5 3 3 61.4 5 4 7由 计 算 结 果 可 知,方 程 组 的 解 为 即 X2,X3,X4=-0.3721,-0.8291,-1.5336,1.4547。9.2.2 QR分 解 法 矩 阵 的 Q R 分 解 就 是 把 矩 阵 分 解 为 一 个 正 交 矩 阵 和 一 个 上 三 角 矩 阵 的 乘 积。MATLAB中 对
8、矩 阵 A 进 行 Q R 分 解 的 函 数 调 用 格 式 如 下:Q,K=qr(/):产 生 一 个 正 交 矩 阵。和 一 个 上 三 角 矩 阵 A,使 得;Q,K,E=qr(/):产 生 一 个 正 交 矩 阵 Q、一 个 上 三 角 矩 阵 R 以 及 一 个 置 换 矩 阵 E,使 得 AE=QR;实 现 分 解 以 后,方 程 组 4r=6 的 解 可 写 为 x=K(Q功)或 x=E(K(QS)。【例 9-31 Q R 分 解 法 求 解 线 性 方 程 组 实 例。用 Q R 分 解 法 求 解 以 下 线 性 方 程 组,x+0.5x2+0.3333%3+0.25%4=
9、10.5%1+0.3333%2+0.25x3+0.2%4=20.3333x1+0.25x2+0.2%3+0.1667x4=20.25xi+0.2必+0.1667x3+0.1429x4=1。解:用 Q R 分 解 法 求 解,在 M A T L A B 命 令 窗 口 中 输 入 求 解 程 序:A=1.0 0 0 00.5 0 0 00.3 3 3 30.2 5 0 00.5 0 0 00.3 3 3 30.2 5 0 00.2 0 0 00.3 3 3 30.2 5 0 00.2 0 0 00.1 6 6 70.2 5 0 0;0.2 0 0 0;0.1 6 6 7;0.1 4 2 9;b=
10、l 2 2 1,;Q,R=q r(A)x=R(Q b)输 出。、A 的 计 算 结 果 为:Q=-0.8 3 8 1 0.5 2 2 6-0.1 5 4 0-0.0 2 6 3-0.4 1 9 1-0.4 4 1 7 0.7 2 7 8 0.3 1 5 7-0.2 7 9 4-0.5 2 8 8-0.1 3 9 5-0.7 8 9 2-0.2 0 9 5-0.5 0 2 1-0.6 5 3 6 0.5 2 6 1-1.1 9 3 2-0.6 7 0 5-0.4 7 4 9-0.3 6 9 80-0.1 1 8 5-0.1 2 5 7-0.1 1 7 5188 第 9 章 线 性 方 程 组 求
11、 解 0 0-0.0062-0.00960 0 0 0.0002输 出 解 的 计 算 结 果 为:x=1.0e+003*0.1160-1.44003.6000-2.3800由 计 算 结 果 可 知,方 程 组 的 解 为 的 广 2,孙 孙=口 的,-1440,3600,-2380。9.2.3 Cholesky 分 解 法 当 系 数 矩 阵 N 正 定 且 对 称 时,它 能 分 解 为 以 下 的 形 式:其 中 K 为 上 三 角 矩 阵,肝 为 R 的 转 置,是 下 三 角 矩 阵,这 种 分 解 称 为 Cholesky分 解。MATLAB中 与 Cholesky分 解 有 关
12、 的 函 数 如 下:R=chol(Z):对 进 行 Cholesky 分 解,使 得 A=RXR;K,p=chol(4):对 A 进 行 Cholesky 分 解,使 得 A=RrRa如 果 4 对 称 正 定,则 p=0;否 则,p 为 一 正 整 数。如 果 4 满 秩,/?是 为 2-1 阶 的 上 三 角 矩 阵,且 有 相/?=4 1:仍-1),l:(p-l);实 现 分 解 后,方 程 组 的 解 写 成 上 大(如 财 的 形 式。【例 9-4 Cholesky分 解 法 求 解 线 性 方 程 组 实 例。用 Cholesky分 解 法 求 解 以 下 线 性 方 程 组,9
13、xi-36切+3013=1*-36xi+1922 一 180工 3=130%1-180%2+180%3=1 0解:用 Cholesky分 解 法 求 解,在 MATLAB命 令 窗 口 中 输 入 求 解 程 序:A=9-36 30;-36 192-180;30-180 180);b=ones(3,1);R=chol(A)x=R(Rzb)输 出/?的 计 算 结 果 为:R=3.0000-12.0000 10.00000 6.9282-8.66030 0 2.2361输 出 解 的 计 算 结 果 为:189精 通 MATLAB科 学 计 算(第 2 版 i-X=1.8 3 3 31.0 8
14、3 30.7 8 3 3由 计 算 结 果 可 知,方 程 组 的 解 为 莺/2,必=1.8333,1.0833,0.7833。起 合 在 M A TLA B中,当 N 正 定 但 不 对 称 时,用 K=chol(N)也 能 得 出 一 个 矩 阵,MATLAB并 不 报 错,但 是 只 要 验 证 一 下 Rt*R就 会 看 出 R*R 并 不 等 于 A,因 此 用 公 式 x=R(*M)就 会 得 出 错 误 的 结 果。190 第 9 章 线 性 方 程 组 求 解 9.2.4 其 他 分 解 法 MATLAB中 还 提 供 了 矩 阵 的 奇 异 值 分 解、Hessenberg
15、分 解 和 Schur分 解 等 函 数,可 用 来 求 解 线 性 方 程 组,如 表 9-1所 示。表 9-1 其 他 的 矩 阵 分 解 函 数 函 数 及 常 用 语 法 说 明 U SM=svd(j)奇 异 值 分 解,使 得 P.H=hess()Hessenberg分 解.使 得 A=P*H*/,其 中 P 为 酉 矩 阵 U,T=schur(/l)Schui分 解,使 得 力。,其 中 U 为 正 交 矩 阵,T为 Schur矩 阵【例 9-5】奇 异 值 分 解 法 求 解 线 性 方 程 组 实 例。用 奇 异 值 分 解 法 求 解 以 下 线 性 方 程 组,X)+0.5
16、x2+0.3333%3+0.25x4+0.2x5=10.5%1+0.3333x2+0.25%3+0.2x4+0.1667x5=0 UzSrV=svd(A)b=l 0 1 0 1 zr x=V*in v(S)*UI*b输 出 分 解 的 I/、S、H 的 计 算 结 果 为:U=-0.7 6 7 9 0.6 0 1 9-0.2 1 4 2 0.0 4 7 2 0.0 0 6 2-0.4 4 5 8-0.2 7 5 9 0.7 2 4 1-0.4 3 2 7-0.1 1 6 7-0.3 2 1 6-0.4 2 4 9 0.1 2 0 5 0.6 6 7 4 0.5 0 6 2-0.2 5 3 4-
17、0.4 4 3 9-0.3 0 9 6 0.2 3 3 0-0.7 6 7 2-0.2 0 9 8-0.4 2 9 0-0.5 6 5 2-0.5 5 7 6 0.3 7 6 2S=1.5 6 7 1 0 0 0 00 0.2 0 8 5 0 0 00 0 0.0 1 1 4 0 00 0 0 0.00 0 3 00 0 0 0 0.,0000V=-0.7 6 7 9 0.6 0 1 9-0.2 1 4 2 0.0 4 7 2 0.0 0 6 2 191精 通 MATLAB科 学 计 算(第 2 蚪-0.4458-0.3216-0.2534-0.2098-0.2759-0.4249-0.443
18、9-0.42900.72410.1205-0.3096-0.5652-0.43270.66740.2330-0.5576-0.11670.5062-0.76720.3762输 出 的 计 算 结 果 为:x=1.0e+004*0.0865-1.24674.4659-5.84042.5426由 计 算 结 果 可 知,方 程 组 的 解 为,x1,X2,x3,X4,x5=104*0.0865,-1.2467,4.4659,-5.8404,2.5426。选 W 奇 异 值 分 解 很 有 用,将 系 数 矩 阵 进 行 奇 异 值 分 解 以 后,A 的 奇 异 值 都 在 S的 对 角 线 上,
19、这 样 就 可 以 粗 略 估 计 N 的 条 件 数,看 看 N 是 否 是 病 态 矩 阵,以 此 决 定 相 应 的 求 解 方 法。如 果 系 数 矩 阵 对 称,那 么 奇 异 值 分 解 后 的。和/都 是 正 交 矩 阵,而 S 是 对 角 矩 阵,这 样 求 解 方 程 就 十 分 方 便。【例 9-6 Hessenberg分 解 法 求 解 线 性 方 程 组 实 例。用 Hessenberg分 解 法 求 解 以 下 线 性 方 程 组,闲+工 2+工 3+工 4=1X+2X2+3汹+414=4X+3X 2+6x3+1014=7X+4X 2+1013+204=-2 o解:用
20、 Hessenberg分 解 法 求 解,在 MATLAB命 令 窗 口 中 输 入 求 解 程 序:A=1 1 1 1;1 2 3 4;1 3 6 10;1 4 10 20;P,H=hess(A)b=l 4 7-21;x=P*inv(H)*P1*b输 出 分 解 的 P、的 计 算 结 果 为:P=0.7754-0.6247-0.0925 0-0.6092-0.7015-0.3698 0192 第 9 章 线 性 方 程 组 求 解 0.166200.34310-0.92450 1.0,00000.2147-0.2654 0 0-0.2654 1.0844 1.0802 00 1.0802
21、7.7009-10.81670 0-10.8167 20.0000输 出 的 计 算 结 果 为:10.0000-33.000036.0000-12.0000由 计 算 结 果 可 知,方 程 组 的 解 为,x i,x2,X 3=10,-3 3,3 6,-1 2 o起 合 Hessenberg矩 阵 指 的 是 第 一 子 对 角 线 以 下 的 元 素 为 0 的 矩 阵。Hessenberg分 解 法 分 解 后 产 生 正 交 矩 阵 和 Hessenberg矩 阵,正 交 矩 阵 的 逆 就 是 其 转 置;如 果 系 数 矩 阵 对 称 的 话,产 生 的 Hessenberg矩
22、阵 就 是 三 对 角 矩 阵,这 种 矩 阵 的 求 逆 速 度 也 很 快,因 此 用 这 种 方 法 求 解 线 性 方 程 组 也 很 快。【例 9-7】S chur分 解 法 求 解 线 性 方 程 组 实 例。S chur分 解 法 求 解 以 下 线 性 方 程 组,汨+工 2+工 3+工 4=1X+2X2+3冷+414=4X 1+3x2+6x3+10 x4=7X+42+10工 3+20%4=-2 o解:用 S chur分 解 法 求 解,在 M A T L A B命 令 窗 口 中 输 入 求 解 程 序:A=1 1 1 1;1 2 3 4;1 3 6 10;1 4 10 20
23、;b=l 4 7-2 UzT=schur(A)x=U*inv(T)*Uf*b输 出 分 解 的 U、7 的 计 算 结 果 为:U=0.3087-0.7873 0.5304 0.0602-0.7231 0.1632 0.6403 0.20120.5946 0.5321 0.3918 0.4581-0.1684-0.2654-0.3939 0.8638T=0.0380 0 0 00 0.4538 0 01 1939.3精 通 M ATLAB科 学 计 算(第 2 版 i-0 0 2.2034 00 0 0 26.3047输 出 的 计 算 结 果 为:x=10.0000-33.000036.00
24、00-12.0000由 计 算 结 果 可 知,方 程 组 的 解 为 XI,X2,X3,X4=10,-33,36,-12。通 多 Schur矩 阵 是 上 三 角 矩 阵,且 其 对 角 元 素 为 被 分 解 矩 阵 的 特 征 值。Schur分 解 法 在 分 解 后 产 生 正 交 矩 阵 和 Schur矩 阵,如 果 系 数 矩 阵 对 称 的 话,产 生 的 Schur矩 阵 就 是 对 角 阵,显 然 对 角 阵 的 逆 矩 阵 非 常 容 易 得 到,因 此 用 这 种 方 法 求 解 线 性 方 程 组 也 很 快。迭 代 法 在 实 际 应 用 中(例 如 在 运 筹 学、
25、图 论 等 领 域 中),往 往 会 出 现 这 样 一 种 情 况:系 数 矩 阵 阶 数 很 高(例 如 IC T),但 系 数 矩 阵 含 零 元 素 相 对 较 多,如 果 用 前 面 的 几 种 分 解 法 去 求 解 的 话,非 零 元 素 反 而 增 多 了,这 就 有 点 得 不 偿 失 了。本 节 介 绍 另 一 类 常 用 的 求 解 方 法 迭 代 法。迭 代 法 是 将 求 一 组 解 转 换 为 求 一 个 近 似 解 序 列 的 过 程,并 用 最 终 的 近 似 解 来 逼 近 真 实 解。迭 代 法 需 要 考 虑 以 下 3 个 重 要 的 问 题。(1)迭
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 精通 MATLAB 科学 计算 王正林 09
限制150内