matlab語言的常量與變量
matlab語言的變量命名規則
- 由一個字母引導,后面可以為其他字符。
- 區分大小寫 如Abc ≠ ABc
matlab的保留常量
以下為系統保留常量,自己定義的變量不能與他們重名
eps %表示機器的精度,其值通常在10^-16左右
i %虛數單位,表示根號-1
j %同i
pi %表示常量π
NaN %不定式,如∞/∞或者0/0的結果
Inf %表示無窮大
lastwarn
lasterr
matlab語言的數據結構
數值型數據
matlab主要使用雙精度的數據結構,滿足IEEE標准,單個數值型數據大小8字節占64位。64位中有11個指數位,52個數據位,1個符號位。表示范圍大概是±1.7x10^{308}。我們可以使用double()命令將變量轉換為雙精度數據結構。
擴展:其他數據類型
- 單精度數據結構single()32位
- uint8(),常用於圖像表示和處理
- int8(),int16(),int32(),uint16(),uint32()
符號型數據結構
-
matlab同時還支持符號型的數據結構,我們可以使用sym(A)這個命令把數值型變量A轉換為符號型數據
-
使用syms聲明符號變量
-
顯示符號變量的任何精度(前n位數值) vpa(A) vpa(A,n) 比如顯示pi的前n位數值
-
顯示符號變量的一些屬性 assumptions()
-
設置符號變量類型 assume(),assumeAlso()
舉個栗子:定義一個大於等於-1且小於5的實數
syms x real;
assume(x>=-1);
assumeAlso(x<5);
再以1/3的存儲內容舉個例子說明符號型數值與雙精度數值的區別
再來個例子:使用符號型數據結構表示數值12345678901234567890(20位,雙精度無法表示,需要使用符號型)
錯誤的方法 實際上計算機還是先將參數轉變成雙精度的數據類型再轉變成符號型的數據,轉換有偏差
正確的方法 將參數用字符串表示,再轉換為字符型
矩陣與向量的輸入
matlab的基本語句結構
直接賦值語句:
variable=expression;
將表達式運算得到的結果賦值給變量,賦值語句的結尾加分號可以阻止運算結果的顯示。如果未指定變量,則表達式的值被賦予保留變量ans。
實數矩陣輸入方法
復數矩陣輸入方法
復數元素1+9i之間不能有空格,不然會引起歧義
函數調用語句
[returned_arguments] = function_name(input_arguments)
[U S V] = svd[X] %函數調用舉例
函數可以通過不同的方式被調用,比如:
- 內核函數,*.m函數
- 匿名函數,inline函數 (不建議使用)
- 重載函數,私有函數等
冒號表達式
matlab下有個很重要的表達式 : 它是定義行向量的有效方法。例如,定義一個從s1到s3,間隔為s2的向量。默認間隔為1。
v = s1:s2:s3
但用這種方法時如果選擇的步距不合適,那么生成的行向量可能就不會包含s3,像下面這種情況:
如果像同時包含s1和s3呢?可以使用下面這個命令linspace(0,pi,50),表示從0到pi生成一個行向量,中間有50個數據點。
linspace(0,pi,50)
另外,如果輸入的布局s2為負數,顯然是錯誤的,但matlab仍然會執行該語句,結果返回一個1x0的空向量
子矩陣提取
B = A(1:2:end,:) %提取矩陣A的奇數行並賦值給B
C = A([1 1 1 1],:) %將A的第一行提取4次並賦值給C
結果如下圖所示
矩陣的代數運算
矩陣的轉置
C = A.' %一般轉置
C = A' %Hermite轉置
矩陣的加減法
C = A + B
C = A - B
它可以直接將兩個矩陣A、B相加減。如果其中一個是標量,那就會把這個標量遍加/減到另一個矩陣上。如果矩陣維數不匹配,系統會報錯。
矩陣的乘法
C = A*B
系統會自動檢測維數是否匹配,不匹配會報錯。
矩陣的除法
矩陣的除法 | 求解線性方程組 | MATLAB解法 | 最小二乘解 |
---|---|---|---|
矩陣左除 | AX=B | X=A\B(注意除號的方向!) 如果原方程不可解則得到最小二乘解 |
$X = A^{-1}B $ |
矩陣右除 | XA=B | X=B/A |
Matlab提供了兩種除法運算:左除(\)和右除(/)。
一般情況下,x=a\b是方程ax =b的解,而x=b/a是方程xa=b的解。
例:a=[1 2 3; 4 2 6; 7 4 9]
b=[4; 1; 2];
x=a\b
則顯示:x=
-1.5000 2.0000 0.5000
如果a為非奇異矩陣,則a\b和b/a可通過a的逆矩陣與b陣得到:
a\b = inv(a)*b b/a = b*inv(a)
什么是奇異矩陣
首先,看這個矩陣是不是方陣(即行數和列數相等的矩陣。若行數和列數不相等,那就談不上奇異矩陣和非奇異矩陣)。 然后,再看此方陣的行列式|A|是否等於0,若等於0,稱矩陣A為奇異矩陣;若不等於0,稱矩陣A為非奇異矩陣。 同時,由|A|≠0可知矩陣A可逆,這樣可以得出另外一個重要結論:可逆矩陣就是非奇異矩陣,非奇異矩陣也是可逆矩陣。 如果A為奇異矩陣,則AX=0有無窮解,AX=b有無窮解或者無解。如果A為非奇異矩陣,則AX=0有且只有唯一零解,AX=b有唯一解。
矩陣的翻轉
矩陣的乘方
前提:矩陣為方陣
數學描述:
matlab描述
F = A^x
求矩陣A的全部三次方艮,並檢驗結果
matlab代碼:
A = [1 2 3;4 5 6;7 8 0];
C = A^(1/3); %結果可以得到一個復數型的三次方根
e = norm(A-C^3); %所求值e為根C的誤差,通常用一個范數表示
運行結果 從圖中可以看出,誤差是很小的
但實際上A還有另外兩個根,如何得到呢?另外兩個根可以通過這個已知根旋轉得到,我們先求出旋轉常數j1,然后把它乘到得到的這個個根C上就能得到另外的兩個根A1 A2
j1=exp(sqrt(-1)*2*pi/3)
A1 = C*j1
A2 = C*j1^2
norm(A-A1^3),norm(A-A2^3)
點運算
點運算是矩陣對應元素直接進行的運算
矩陣的點乘法
矩陣的點乘方
比如求矩陣A的點乘方
那么點運算的用處在哪呢?答案是可以用來繪制函數曲線圖.比方說我們想要做一條y=x^2的曲線,那么我們需要先生成一個x向量,然后對x向量的每一個值單獨做平方運算.這不正是點乘方運算嗎?由此我們可以用以下代碼實現
x=--5:5
y=x.^2
矩陣的其他運算
矩陣的邏輯運算
- 與運算
A&B
- 或運算
A|B
- 非運算
B = ~A
- 異或運算
xor(A,B)
比較運算
各種允許的比較關系有:>
,>=
,<
,<=
,==
,~=
,find()
,all()
,any()
解析結果的化簡變換
函數simplify()
用於數學公式的化簡
s1 = simplify(s);
其他常用化簡函數
numden() 提取表達式分子分母
collect() 合並同類項
expand() 多項式展開
factor() 因式分解
例:化簡多項式
syms s;
P = (s+3)^2*(s^2+3*s+2)*(s^3+12*s^2+48*s+64); %輸入多項式
P1 = simplify(P) %得到最簡形式
P3 = factor(P),P4 = prod(P3) %用factor得到所有的公因式,再用prod將所有的公因數相乘
運行結果:
變量替換及轉成Latex表示
將函數表達式f中的x1
全部替換為x1*
應用舉例:試用s=(z-1)/(z+1)
對P(s)進行雙線性變換
syms s z;
P = (s+3)^2*(s^2+3*s+2)*(s^3+12*s^2+48*s+64); %輸入多項式
P1 = simplify(subs(P,s,(z-1)/(z+1))) %得到變換后的最簡形式
P3 = latex(P1) %轉為latex形式
運行結果
latex顯示如下:
matlab基本數論運算
函數 | 調用格式 | 功能 |
---|---|---|
floor() |
n=floor(x) | 向下(負無窮大)取整 |
ceil() |
n=ceil(x) | 向上取整 |
round() |
n=round(x) | 四舍五入取整 |
fix() |
n=fix(x) | 向0方向取整 |
rat() |
[n,d]=rat(x) | 將有理數轉化為分式 |
rem() |
B=rem(A,C) | 求余數 |
gcd() |
k=gcd(n,m) | 求最大公約數 |
lcm() |
k=lcm(n,m) | 求最小公倍數 |
factor() |
factor(n) | 因數分解 |
isprime() |
v1=isprime(v) | 判斷是否是素數(質數) |
perms() |
perms(1:5)或 perms('abcde') |
求全排列 |
流程結構
利用matlab的流程結構,我們可以編寫出復雜程序,實現更高級的功能。目前matlab支持的流程結構有:
- 循環結構
- 轉移結構
- 開關結構
- 試探結構
matlab的循環結構
for循環結構
for i=v , 循環體 , end
其執行機制為:v為一個向量
,循環變量i每次從v向量中取一個數值,執行一次循環體的內容,如此下去直到執行完v向量中的所有分量。v的內容可以任意排列
while循環結構
whlie(判斷條件) %滿足條件則進入循環
循環體,
end
用循環求解最小的m,使下式成立
\[\sum_{i=1}^{m} {i>10000} \]s=0;m=0; while(s<=10000),m=m+1; s=s+m; end,s,m
運行結果
條件轉移結構
if(condition)
statement group
end
%或者
if(conditon1)
statement group 1
elseif(conditon2)
statement group 2
else
statement group 3
end
用法與c語言中的if
else
基本差不多
開關結構
switch switch expression
case expression 1 ,statement 1
case {expression 2,expression 3 ,··· expression m},statement 2
otherwise statement n
end
與c語言也有些相似,不同之處在於執行完statement之后程序會自動結束,而在c語言中要使用break才能實現這樣的功能。如果開關表達式witch expression
滿足{expression 2,expression 3 ,··· expression m}
其中之一,那么程序就會執行這個段落,執行完后跳出此開關結構。
試探結構
try,
statement group 1,
catch,
statement group 2,
end
程序會先嘗試執行語句段1,如果不出現錯誤,那么這個結構就執行結束了。如果出現錯誤,程序就轉到語句段2去執行,執行完后結束此結構
函數編寫
matlab編程有腳本程序
和M-函數
,但更推薦使用M-函數
的形式,這是因為函數更靈活,可以應用到各種場合而不需要反復修改源程序,只需要改變參數就好。
matlab語言函數的基本結構
函數可以看作是一個信息處理單元,它接受輸入變量,然后經過處理計算,再返回相應的變量到上一級
函數的程序結構
function [return vars] = funname(input vars)
comments led by %
input and output variables check
main body of the function
通常,函數名funname
要起一個有意義的名字,最后這個函數將存為一個.m
文件(通常以函數名命名
例編寫一個函數生成nxm的Hilbert矩陣
\[h_{i,j}=\frac {1} {i+j-1} \]要求
- 輸入變量n,m,輸出變量H
- 若只給出一個輸入參數,則自動生成方陣
- 在函數中給出合適的幫助信息
- 檢測輸入和返回變量的個數nargin,nargout
function H=myhilb(n,m) %定義函數myhilb
if nargin == 1 ,m=n; end %如果輸入變量只有一個,那么強制將m=n
for i=1:n,for j=1:m %使用循環生成Hilbert矩陣
H(i,j)=1/(i+j-1);
end,end
函數的遞歸調用
例題:利用函數的遞歸調用計算階乘
原理 n! = n(n-1)!
遞歸調用函數編寫如下
function k = my_fact(n)
if nargin ~= 1,
error('error:只允許輸入一個變量');
end
if(abs(n-floor(n))>eps | n<0),
error('error:n不能為非負整數');
end %前半部分用於檢測n是否為非負整數
if n>1 k=n*my_fact(n-1); %如果n>1,則使用遞歸調用求階乘
elseif any([0 1]==n),k=1;end %如果n為0或1,則設置一個出口,讓k=1
除了我們自己編寫函數外,matlab自帶也有計算階乘的函數:factorial(n)
,prod(1:n)
,gamma(1+n)
,gamma(1+sym(n))
,factotial(sym(n))
。
可變輸入輸出參數個數
如何再matlab中編寫可變輸入輸出參數個數的函數呢?首先我們需要先了解一下輸入變量是怎么傳遞到函數里去的。
由圖可知,輸入輸出變量是存儲在varargin
,varargout
兩個數組中的,我們可以使用花括號{}
來提取某一個變量
例:任意多輸入變量
conv()
可以計算兩個多項式的積,試使用varargin實現任意個多項式的積
matlab編程
%思路 每次從varargin中提取出一個做累積
function a=convs(varargin),a=1;
for i=1:nargin, a=conv(a,varargin{i});end
下面我們來具體試一下自己編寫的函數
P=[1 2 4 0 5];
Q=[1 2];
F=[1 2 3]; %P Q F為三個多項式的系數
%傳統求法
E=conv(conv(P,Q),F)
%自編寫函數求解
D=convs(P,Q,F)
%我還能求更多 ^_+
G=convs(P,Q,F,[1,1],[1,3],[1,1])
運行結果
inline函數與匿名函數
inline函數會造成功能重疊,目前不建議大家使用。下面介紹一下匿名函數的使用
匿名函數的形式為:f=@(list of variables)function_contents
。它的好處是不用單獨編寫一個.m文件,但只適合處理簡單的函數。
二維曲線的繪制
二維圖形繪制基本語句
假設有兩個序列t=t1,t2,···,tn y=y1,y2,···,yn。構成向量t=[t1,t2,···,tn],y=[y1,y2,···,yn]。想要利用這兩組數據繪圖,只需要輸入plot(t,y)
即可
例:繪制函數曲線
繪制函數\[y=sin(tan(x))-tan(sin(x)),x\in[-\pi,\pi] \]
matlab代碼
x=[-pi:0.005:pi];
y=sin(tan(x))-tan(sin(x));
plot(x,y)
運行結果
例:繪制分段函數
繪制飽和函數方程
\[y= \begin{cases} 1.1sign(x),|x|>1.1 \\ x, \qquad\qquad |x|\leq1.1 \end{cases} \]
matlab繪圖語句(互斥條件才能這么寫)
x=[-2:0.02:2];
%互斥條件下將兩段表達式分別與其條件做點乘運算再相加即可得函數圖像
y=1.1*sign(x).*(abs(x)>1.1)+x.*(abs(x)<=1.1);
plot(x,y)
其他調用格式
情況1中y有m行,最終以t為x軸生成m條曲線。情況2、3也會得出多條曲線
更一般的調用格式
也可以重新設置曲線樣式
多縱軸曲線繪制方法
特殊二維圖形
其他二維圖形繪制語句
不同的繪制函數
常用函數 | 意義 | 常用函數 | 意義 |
---|---|---|---|
bar(x,y) |
二維條形圖 | comet(x,y) |
彗星狀軌跡圖 |
compass(x,y) |
羅盤圖 | comet(x,y,ym,yM) |
誤差限圖形 |
feather(x,y) |
羽毛狀圖 | fill(x,y,c) |
二維填充函數 |
hist(y,n) |
直方圖 | loglog(x,y) |
對數圖 |
polar(x,y) |
極坐標圖 | quiver(x,y) |
磁力線圖 |
stairs(x,y) |
階梯圖形 | stem(x,y) |
火柴桿圖 |
semilogx(x,y) |
x-半對數圖 | semilogy(x,y) |
y-半對數圖 |
例 用不同的函數繪制正弦函數圖形
t=0:.2:2*pi; y=sin(t); subplot(2,2,1), stairs(t,y), subplot(2,2,2), stem(t,y), subplot(2,2,3), bar(t,y), subplot(2,2,4), semilogx(t,y)
其中subplot(x,y,n)
函數用於在窗口的不同位置繪制圖形,x表示行,y表示列,n表示第幾幅圖
運行結果
隱函數繪制及應用
隱函數即滿足f(x,y)=0
方程的x,y之間的關系式。由於很多隱函數無法求出x,y之間的顯式關系,所以無法定義x向量再求出y向量從而進行繪制。matlab提供的ezplot(隱函數表達式,[xm,xM,ym,yM])
可以直接繪制隱函數曲線
例:繪制出下列隱函數的曲線
\[f(x,y)=x^2sin(x+y^2)+y^2e^{x+y}+5cos(x^2+y)=0 \]實現代碼
h=ezplot('x^2 *sin(x+y^2) +y^2*exp(x+y)+5*cos(x^2+y)',[-10 10]); %設置定義域為-10到10 set(h,'Color','b')
運行結果截圖
數據文件的存儲和讀取
matlab提供save
和load
命令來對外部文件讀寫數據
save mydat A B C %以二進制形式存儲ABC的數據
save /ascii mydat.dat A B C %以可讀形式存儲
X = load(filename) %讀取文件中的數據到x
X = xlsread(filename,range) %從excel文件中讀取數據 range可以為 'B6:C67'
xlswrite(filename,range) %向excel文件中寫入數據
例:已知excel文件census.xls給出某省人口數,第5-67行存儲數據,B列存儲年份,C列存儲人口。試用matlab對兩列數據繪圖
> X=xlsread('census.xls','B5:C67'); %讀入數據后X為一個兩列的矩陣 t=X(:,1); %第一列所有行數據給t p=X(:,2); %第二列所有行數據給p plot(t,p) %繪圖
三維圖形表示
三維曲線繪制
二維曲線繪制函數plot()
可以擴展到三維曲線的繪制中。這時可以使用plot3()
函數來繪制三維曲線,其調用格式如下,其中選項和二維曲線繪制的完全一致。當然也有其它三維曲線繪制函數,大部分是由二維曲線繪制函數擴展而來。
plot3(x,y,z)
plot3(x1,y1,z1,選項1,x2,y2,z2,選項2,x3,y3,z3,選項3)
例:試繪制如下參數方程的三維曲線
\[\begin{cases} x(t)=t^3sin(3t)e^{-t} \\ y(t)=t^3cos(3t)e^{-t} \\ z=t^2 \end{cases} \]matlab代碼
t=0:0.01:2*pi; %先構造t向量,注意下面的點運算 x=t.^3.*sin(3*t).*exp(-t); y=t.^3.*cos(3*t).*exp(-t); z=t.^2; plot3(x,y,z), grid %繪制三維曲線
運行結果截圖
三維曲面繪制
如果已知二元函數z=f(x,y)
,則可以繪制出該函數的三維曲面圖。在繪制三維曲面圖之前,應該先調用meshgrid()
函數生成網格矩陣數據x和y,這樣就可以按函數公式用點運算的方式計算出z矩陣,之后就可以用mesh()
或surf()
函數進行三維圖形繪制。具體的函數調用格式為
[x,y]=meshgrid(v1,v2) %生成網格數據,v1,v2為x、y軸的定義域及步距
z = ··· ,如z=x.*y %計算二元函數的z矩陣
surf(x,y,z)或mesh(x,y,z)%surf繪制表面圖,mesh繪制網格圖
例繪制出下列函數的三位表面圖形,定義域自取 $$ z=f(x,y)=(x^2-2x)e^{-x^2-y^2-xy} $$[x,y]=meshgrid(-3:0.1:3,-2:0.1:2); z=(x.^2-2*x).*exp(-x.^2-y.^2-x.*y); surf(x,y,z)
運行結果截圖
參數方程的表面圖
使用matlab繪制的方法如下
ezsurf(fx,fy,fz,[um,uM,vm,vM]) %如果不指定,默認區間為[-2pi,2pi]
例:繪制傳說中的莫比烏斯帶 $$ \begin{cases} x=cosu + v cosu cos\frac{u}{2} \\ y=sinu + v sinu cos\frac{u}{2} \\ z=vsin\frac{u}{2} \end{cases} $$ matlab代碼syms u v; x=cos(u)+v*cos(u)*cos(u/2); y=sin(u)+v*sin(u)*cos(u/2); z=v*sin(u/2); ezsurf(x,y,z,[0,2*pi,-0.5,0.5])
運行結果截圖
球面的繪制
繪制球面的話首先使用matlab函數生成數據
[x,y,z] = sphere(n)
函數會生成3個(n+1)x(n+1)的矩陣,根據這三個矩陣,我們就能畫出它的圖形了
例:繪制兩個球,一個圓心在原點,半徑為1,另一個圓心(0.9,-0.8,0.6),半徑為0.3
matlab代碼
[x,y,z]=sphere(50); surf(x,y,z), hold on, %繪制默認球體,並不讓其消失 x1=0.3*x+0.9; y1=0.3*y-0.8; z1=0.3*z+0.6; surf(x1,y1,z1) %繪制第二個球體
運行結果:
柱面的繪制
曲線r沿z軸旋轉一周可以得到廣泛意義下的柱面,在matlab中,我們可以使用[x,y,z]=cylinder(r,n)
生成柱面
例:試繪制如下柱面
\[r(z)=e^{-0.5Z^2}sinz,\quad z \in (-1,3) \]matlab代碼:
z0=-1:0.1:3; r=exp(-z0.^2/2).*sin(z0); [x,y,z]=cylinder(r); z=-1+4*z; surf(x,y,z)
運行結果截圖
特殊三維圖形
等高線繪制
- 直接繪制
contour(x,y,z,n) %n為等高線條數
- 帶有標志的等高線
[C,h]=cotour(x,y,z,n) %返回等高線的句柄及要標注的值
clabel(C,h) %疊加到等高線圖上
- 繪制三維等高線
contour3(x,y,z,n)
contourf(x,y,z,n)
三維隱函數的繪制
matlab本身並沒有提供三維隱函數繪制函數,我們需要到math work網站上下載ezimplot3()
函數。
ezimplot3(fun,[xm,xM,ym,yM,zm,zM])
默認的范圍是[-2pi,2pi]
三維圖形的視角設置
除了使用工具欄中的三維旋轉
按鈕,matlab還提供了視角設置函數view(α,β)
,也可以使用[α,β]=view(3)
讀出當前三維視角數據。視角的定義如下圖所示:
例:試繪制函數 $ z=f(x,y)=(x2-2x)e{-x2-y2-xy} $ 的三視圖
matlab代碼
[x,y] = meshgrid(-3:0.1:3,-2:0.1:2); z=(x.^2-2*x).*exp(-x.^2-y.^2-x.*y); subplot(224), surf(x,y,z), subplot(221), surf(x,y,z), view(0,90); %俯視圖 subplot(222), surf(x,y,z), view(90,0); %側視圖,右視圖 subplot(223), surf(x,y,z), view(0,0); %正視圖
運行結果截圖
三維曲面的旋轉
三維曲面旋轉函數,其中h為三維曲面的句柄,v為旋轉基軸,α為旋轉角度。如果繪圖時h=surf(x,y,z
,那么h就是該三維圖像的句柄。v是一個三維向量
,如v=[1 1 1],那么轉軸就是向量(1,1,1)所在的軸線。
rotate(h,v,α)
rotate(1,0,0) %繞x軸正方向旋轉
例:把上例的圖形繞x軸正方向旋轉360度的動畫,每0.01s旋轉一度,使用循環結構旋轉
matlab代碼
fingure; h=surf(x,y,z); r_ax=[1 0 0]; axis tight, for i=0:360, rotate(h,r_ax,1); pause(0.01), end