MATLAB怎么求解線性方程組? 附完整代碼和案例

線性方程組的直接解法方法很多,包括Gauss消去法、選主元消去法、平方根法和追趕法等等。但是在MATLAB中,可以直接利用“\”或者“/”來解決問題。這兩種方法的內(nèi)部包含非常多的自適應(yīng)算法,比如對超定方程使用最小二乘法;對欠定方程給出誤差范數(shù)最小的一個解;對三對角陣方程組使用追趕法。
一. 直接求解:矩陣除法
對線性方程求解,MATLAB調(diào)用格式如下:
x=A\B
調(diào)用此函數(shù)時,矩陣A、B必須具有相同的行數(shù)。如果矩陣A沒有正確縮放或者接近奇異矩陣,運行代碼時,MATLAB就會顯示警告信息。
對矩陣A可以分成如下三種情況:
- 如果A是標量,那么A\B就等同于A.\B
- 如果A是
方陣,B是n行矩陣,那么A\B就是方程A*X=B的解
- 如果A是
矩陣,而且
,B是m行矩陣,那么A\B返回方程組A*X=B的最小二乘解
例題1
解如下方程:
解:
MATLAB代碼如下:
clc;clear; A=[0.4096 0.1234 0.3678 0.2943;0.2246 0.3872 0.4015 0.1129; 0.3645 0.1920 0.3781 0.0643;0.1784 0.4002 0.2786 0.3927]; b=[0.4043;0.1150;0.4240;-0.2557]; x=A\b; x' %將結(jié)果輸出為行向量
運行結(jié)果:
ans =
- 0.332683779325239
- -1.412140862756104
- 1.602847655449206
- -0.500293786006566
例題2
A為4階的幻方矩陣,求解線性方程組Ax=b。b的表達式如下:
備注:幻方矩陣的定義:如果一個數(shù)組具有相同行列且每行,每列和對角線上的和都一樣,則成這些數(shù)組則成為魔方矩陣,又叫幻方矩陣。
解:
MATLAB代碼如下:
clc;clear; A=magic(4); %生成四階的幻方矩陣 b=[34;34;34;34]; x=A\b; x'
運行結(jié)果:
ans =
- 0.980392156862745
- 0.941176470588235
- 1.058823529411765
- 1.019607843137255
分析:
4階的幻方矩陣是奇異矩陣,奇異矩陣也能求解,但是MATLAB會生成警告,警告信息如下:
警告: 矩陣接近奇異值,或者縮放錯誤。結(jié)果可能不準確。RCOND = 4.625929e-18。;
例題3
對線性方程組Ax=b進行求解,并分析是否有誤差。A,b矩陣如下:
,
解:
MATLAB代碼如下:
clc;clear; A=[1 2 0;0 4 3]; b=[8;18]; x=A\b E=norm(b-A*x) %求范數(shù)誤差
運行結(jié)果:
x =
- 0
- 3.999999999999997
- 0.666666666666670
E =
7.160723346098895e-15
分析:此方程屬于欠定方程,MATLAB采用最小二乘法進行求解,所以會出現(xiàn)誤差
另外最后舉一個利用稀疏矩陣對簡單線性方程組Ax=B進行求解。
MATLAB代碼:
clc;clear; A=sparse([0 2 0 1 0;4 -1 -1 0 0;0 0 0 3 -6;-2 0 0 0 2;0 0 4 2 0]);%稀疏矩陣 B=sparse([8;-1;-18;8;20]); %稀疏矩陣 x=A\B E1=norm(B-A*x) %范數(shù)誤差
運行結(jié)果:
x =
- (1,1) 1.000000000000000
- (2,1) 2.000000000000000
- (3,1) 3.000000000000000
- (4,1) 4.000000000000001
- (5,1) 5.000000000000001
E1 =
4.189529226675416e-15
二. 直接求解:判斷求解
給定矩陣A,B如下:
形成解的判定矩陣C如下:
線性方程組有解的判斷定理將分成三種情況:
2.1 m=n且rank(A)=rank(C)=n
此時,方程組有唯一的解
例題4
求解如下方程組:
解:
MATLAB代碼如下:
clc;clear; A=[1 2 3 4;4 3 2 1;1 3 2 4;4 1 3 2]; B=[5 1;4 2;3 3;2 4]; C=[A B]; %判定前提條件 rank_A=rank(A) rank_C=rank(C) %求解 x1=inv(A)*B %計算范數(shù)誤差 E1=norm(A*x1-B) %計算精確解 x2=inv(sym(A))*B %計算范數(shù)誤差 E2=norm(A*x2-B)
運行結(jié)果:
- rank_A =4
- rank_C = 4
x1 =
- -1.800000000000000 2.399999999999999
- 1.866666666666666 -1.266666666666667
- 3.866666666666667 -3.266666666666667
- -2.133333333333333 2.733333333333333
E1 =7.291088482824584e-15
x2 =
- [ -9/5, 12/5]
- [28/15, -19/15]
- [ 58/15, -49/15]
- [ -32/15, 41/15]
E2 =0
2.2 rank(A)=rank(C)=r<n
此時方程組有無窮多個解。將原方程組可以轉(zhuǎn)化為對應(yīng)的齊次方程組的解,形式如下:
此時需要求取A矩陣的化零矩陣,MATLAB格式如下:
Z=null(A)
或者求取A矩陣的化零矩陣的規(guī)范形式,MATLAB格式如下:
Z=null(A,'r')
例題5
求解以下方程組:
解:
MATLAB代碼如下:
clc;clear; A=[1 2 3 4;2 2 1 1;2 4 6 8;4 4 2 2]; B=[1;3;2;6]; C=[A B]; %判斷可解性 rank=[rank(A),rank(C)] %求解出規(guī)范化的化零空間 Z=null(sym(A)) %求特解 x0=sym(pinv(A)*B) %驗證得出的解 a1=randn(1); a2=rand(1); %a1和a2是不同分布的隨機數(shù) x1=a1*Z(:,1)+a2*Z(:,2)+x0; E=norm(double(A*x1-B)) %將通解表示出來 syms a3 a4; x2=a3*Z(:,1)+a4*Z(:,2)+x0
運行結(jié)果:
rank =2
分析:說明有解,且解不止一個
Z =
- [ 2, 3]
- [ -5/2, -7/2]
- [ 1, 0]
- [ 0, 1]
分析:此結(jié)果可以解釋為,如下等式:
x0 =
- 125/131
- 96/131
- -10/131
- -39/131
E =0
分析:此方法求出的結(jié)果沒有誤差
x2 =
- 2*a3 + 3*a4 + 125/131
- 96/131 - (7*a4)/2 - (5*a3)/2
- a3 - 10/131
- a4 - 39/131
分析:此結(jié)果可以寫成通式如下:
2.3
此時只能采用摩爾-彭羅斯(Moore-Penrose)廣義逆求解出其最小二乘法,MATLAB格式如下:
x=pinv(A)*B
這種方法只能使誤差范數(shù)測度||Ax-B||取的最小值,并不能完全符合原始的代數(shù)方程。
三. 矩陣求逆解線性方程組
通過反轉(zhuǎn)系數(shù)矩陣A,也可以實現(xiàn)對線性方程組Ax=b的求解。不幸的是,與之前反斜杠計算的方法相比,此方法的運算速度會更慢,殘差也相對較大。但其實,它也有它的優(yōu)點,我們來看一道例題。
例題 6
利用八階的幻方矩陣,來比較MATLAB中x1=A\b和x2=pinv(A)*b兩種求解方法的精確度。
解:
MATLAB代碼如下:
clc;clear; A=magic(8); A1=A(:,1:6); %行數(shù)全取,列數(shù)保留1~6列,具體原因見"分析" rank_A1=rank(A1) b=260*ones(8,1); %260是原A矩陣的幻數(shù)和 %使用反斜杠法求解 x1=A1\b norm_x1=norm(x1) E1=norm(A1*x1-b) %使用pinv()函數(shù)求解 x2=pinv(A1)*b norm_x2=norm(x2) E2=norm(A1*x2-b)
運行結(jié)果:
rank_A1 =3
警告: 秩虧,秩 = 3,tol = 1.882938e-13。
(由于A1矩陣非方陣,所以運算秩時矩陣出現(xiàn)警告,與矩陣秩相關(guān)的介紹可以看以下文章)
基于MATLAB的矩陣性質(zhì):行列式,秩,跡,范數(shù),特征多項式與矩陣多項式_嘮嗑!的博客-CSDN博客
x1 =
- 2.999999999999998
- 4.000000000000000
- 0
- 0
- 1.000000000000002
- 0
- norm_x1 =5.099019513592784
- E1 =1.392373714442771e-13
x2 =
- 1.153846153846151
- 1.461538461538463
- 1.384615384615385
- 1.384615384615383
- 1.461538461538461
- 1.153846153846153
- norm_x2 =3.281650616569467
- E2 =4.019436694230464e-13
分析:
- (1)將A矩陣轉(zhuǎn)換為A1,因為當只有六列時,方程仍然是一致的,依舊存在解,而且解并非全由1組成。另外,由于矩陣低秩,所以會有無數(shù)個解。
- (2)根據(jù)范數(shù)誤差計算的結(jié)果,兩種求解方法精度一致
(3)解x1的特殊之處在于它只有三個非零元素;解x2的特殊之處在于它的norm(x2)非常小
相關(guān)文章
matlab怎么畫函數(shù)圖像? MATLAB繪制函數(shù)圖像的實例教程
有時候我們在使用matlab的時候,想畫函數(shù)圖像,怎么畫呢,?繪制方法很簡單,下面來分享一下2025-01-11MATLAB如何調(diào)用function? 一文看懂functionfunction函數(shù)的試用技巧
MATLAB的functions函數(shù),用于查詢和調(diào)試函數(shù)句柄信息,該怎么了解并使用這個函數(shù)呢?詳細請看下文介紹2025-01-10matlab提示內(nèi)存不足怎么辦? MATLAB內(nèi)存不足及MAT文件版本過低解決方案
Matlab作為一款科學計算軟件,其內(nèi)存管理變得尤為重要,當遇到“內(nèi)存不足”的錯誤提示時,許多用戶可能會感到困擾,下面我們就來看看詳細解決辦法2025-01-10- 很多小伙伴還不了解matlab怎么生成圖像,其實很簡單的我們只要準備好的圖像復(fù)制到matlab的工作目錄中,然后入返回指令[data=imread在主界面將圖像轉(zhuǎn)換為數(shù)據(jù)就可以了,詳2024-01-20
matlab怎么刪除某一行錯誤代碼 matlab把錯誤行刪掉的技巧
用戶在使用matlab時多打一行代碼或者出現(xiàn)一行錯誤的代碼要怎么刪除,其實只要選中要刪除的代碼直接注釋或者delete刪除即可,詳細請看下文介紹2024-01-20- 很多用戶在不需要用到matlab這個軟件之后,想把它給卸載掉,但是不知道怎么卸載,不會操作,要想把這個軟件卸載干凈,可以直接到控制面板中把它卸載,詳細請看下文介紹2024-01-20
matlab怎么將數(shù)據(jù)從大到小排序? matlab正序或倒敘數(shù)據(jù)排序的技巧
matlab怎么將數(shù)據(jù)從大到小排序?matlab中的數(shù)據(jù)想要排序,該怎么從大到小排序,或者從小到大排序呢?下面我們就來看看matlab正序或倒敘數(shù)據(jù)排序的技巧2023-10-26matlab積分函數(shù)怎么寫? matlab求積分的教程
matlab積分函數(shù)怎么寫?matlab中想要求積分,該怎么操作呢?下面我們就來看看matlab求積分的教程2023-10-26matlab中for循環(huán)怎么用? MATLAB里for函數(shù)依次讀取的用法
matlab中for循環(huán)怎么用?matlab中想要使用for循環(huán)函數(shù),該怎么操作呢?下面我們就來看看MATLAB里for函數(shù)依次讀取的用法2023-10-26- Matlab常用快捷鍵有哪些?我們今天來看看一些常用的Matlab快捷操作命令掌握這些操作可以極大地提高使用Matlab時的效率,詳細請看下文介紹2023-10-26