欧美bbbwbbbw肥妇,免费乱码人妻系列日韩,一级黄片

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

  發(fā)布時間:2025-01-11 16:17:14   作者:佚名   我要評論
線性方程組是數(shù)學中一類重要的問題,廣泛應(yīng)用于科學、工程和經(jīng)濟等領(lǐng)域,在Matlab中,我們可以利用內(nèi)置的函數(shù)和工具箱來解決線性方程組,本文將介紹如何使用Matlab求解線性方程組,并提供相應(yīng)的源代碼示例

線性方程組的直接解法方法很多,包括Gauss消去法、選主元消去法、平方根法和追趕法等等。但是在MATLAB中,可以直接利用“\”或者“/”來解決問題。這兩種方法的內(nèi)部包含非常多的自適應(yīng)算法,比如對超定方程使用最小二乘法;對欠定方程給出誤差范數(shù)最小的一個解;對三對角陣方程組使用追趕法。

一. 直接求解:矩陣除法

對線性方程A*X=B求解,MATLAB調(diào)用格式如下:

x=A\B

調(diào)用此函數(shù)時,矩陣A、B必須具有相同的行數(shù)。如果矩陣A沒有正確縮放或者接近奇異矩陣,運行代碼時,MATLAB就會顯示警告信息。

對矩陣A可以分成如下三種情況:

  • 如果A是標量,那么A\B就等同于A.\B
  • 如果A是n\times n方陣,B是n行矩陣,那么A\B就是方程A*X=B的解
  • 如果A是m\times n矩陣,而且m\neq n,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的表達式如下:

b=\begin{bmatrix}34\\34\\34\\34 \end{bmatrix}

備注:幻方矩陣的定義:如果一個數(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矩陣如下:

A=\begin{bmatrix}1&2&0\\0&4&3 \end{bmatrix},b=\begin{bmatrix}8\\18 \end{bmatrix}

解:

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如下:

C=\begin{bmatrix}a_{11}&a_{12}&\cdots&a_{1n}&b_{11}&b_{12}&\cdots&b_{1p}\\ a_{21}&a_{22}&\cdots&a_{2n}&b_{21}&b_{22}&\cdots&b_{2p}\\ \vdots&\vdots&\ddots&\vdots&\vdots&\vdots&\ddots&\vdots\\ a_{m1}&a_{m2}&\cdots&a_{mn}&b_{m1}&b_{m2}&\cdots&b_{mp} \end{bmatrix}

線性方程組有解的判斷定理將分成三種情況:

2.1 m=n且rank(A)=rank(C)=n

此時,方程組有唯一的解 x=A^{-1}B

例題4

求解如下方程組:

\begin{bmatrix}1&2&3&4\\ 4&3&2&1\\ 1&3&2&4\\ 4&1&3&2 \end{bmatrix}X=\begin{bmatrix}5&1\\4&2\\3&3\\2&4 \end{bmatrix}

解:

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)的齊次方程組的解,形式如下:

\hat x=\alpha_1x_1+\alpha_2x_2+\ldots+\alpha_{n-r}x_{n-r},\quad \alpha_i,i=1,2,\cdots,n-r

此時需要求取A矩陣的化零矩陣,MATLAB格式如下:

Z=null(A)

 或者求取A矩陣的化零矩陣的規(guī)范形式,MATLAB格式如下:

Z=null(A,'r')

例題5

求解以下方程組:

\begin{bmatrix}1&2&3&4\\ 2&2&1&1\\ 2&4&6&8\\ 4&4&2&2 \end{bmatrix}X=\begin{bmatrix}1\\3\\2\\6 \end{bmatrix}

解:

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é)果可以解釋為,如下等式:

\begin{bmatrix}x_1\\x_2\\x_3\\x_4 \end{bmatrix}=\begin{bmatrix}2&3\\-2.5&-3.5\\1&0\\0&1 \end{bmatrix}\begin{bmatrix}x_3\\x_4 \end{bmatrix}

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é)果可以寫成通式如下:

\begin{bmatrix}x_1\\x_2\\x_3\\x_4 \end{bmatrix}=a_3\begin{bmatrix}2\\-2.5\\1\\0 \end{bmatrix}+a_4\begin{bmatrix}3\\-3.5\\0\\1 \end{bmatrix}+\begin{bmatrix}0.9542\\0.7328\\-0.0763\\-0.2977 \end{bmatrix}

2.3 rank(A)\leq rank(C)

此時只能采用摩爾-彭羅斯(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)文章

最新評論