Ch51基本应用领域.docx
- 文档编号:9956661
- 上传时间:2023-02-07
- 格式:DOCX
- 页数:29
- 大小:231.21KB
Ch51基本应用领域.docx
《Ch51基本应用领域.docx》由会员分享,可在线阅读,更多相关《Ch51基本应用领域.docx(29页珍藏版)》请在冰豆网上搜索。
Ch51基本应用领域
第五章MATLAB基本应用领域
MATLAB的应用领域非常广泛,从最基本的线性代数、泛函分析,到应用广泛的信号处理、控制系统、通信系统,直到最新技术领域神经网络、小波理论等。
本章主要介绍MATLAB在基本应用领域:
线性代数、多项式与内插、数据分析与统计、泛函分析、常微分方程求解中的应用。
§5.1线性代数
5.1.1矩阵运算—参Ch1-基本操作
Ex.1设A=[13;24]判断A中元素是否大于2。
>>A=[13;24];
>>A>2
ans=
01
11
Ex.2找出A中小于0的元素的位置。
>>A=[1-23;45-4;5-67]
A=
1-23
45-4
5-67
>>B=find(A<0)
B=
4
6
8
5.1.2矢量范数和矩阵范数
1.矢量x的p-范数定义为
特别地,有
向量的1-范数:
向量的2-范数:
向量的-范数:
一般用norm(x,p)函数实现,p缺省时为p=2.如
>>v=[20-1];
>>n=[norm(v,1),norm(v),norm(v,inf)]
n=3.00002.23612.0000
2.矩阵A的p-范数定义为
特别地,有
矩阵的1-范数:
向量的2-范数:
向量的-范数:
一般用norm(A,p)函数实现,p缺省时为p=2.如
>>A=fix(10*rand(3,2))
A=
20
66
83
>>N=[norm(A,1),norm(A),norm(A,inf)]
N=16.000011.889512.0000
5.1.3线性代数方程求解
1.当A为非奇异时,解为
例1设矩阵
>>A=magic(3)
A=
816
357
492
举4种情形求解:
(1)Ax=b,b=[1:
3]’,x为未知列向量
>>x=A\b%给出方程组的解
x=
0.0500
0.3000
0.0500
>>A*x%验证解的正确性
ans=
1.0000
2.0000
3.0000
(2)xA=b,b=[1:
3],x为未知行向量
>>b=[1:
3]%验证解的正确性
b=123
>>x=b/A%给出方程组的解
x=-0.03330.4667-0.0333
>>x*A
ans=1.00002.00003.0000
(3)AX=B
>>B=Pascal(3)
B=
111
123
136
>>X=A\B
X=
0.06670.05000.0972
0.06670.30000.6389
0.06670.0500-0.0694
>>A*X%验证解的正确性
ans=
1.00001.00001.0000
1.00002.00003.0000
1.00003.00006.0000
(4)XA=B%矩阵方程
>>X=B/A%给出方程组的解
X=
0.06670.06670.0667
-0.03330.4667-0.0333
-0.15281.0556-0.2361
>>X*A%验证解的正确性
ans=
1.00001.00001.0000
1.00002.00003.0000
1.00003.00006.0000
2.当A为mn(m>n)矩阵时,AX=B中方程个数多于未知量个数,应采用最小二乘法来求解。
例如,对一维测量数据,根据数据图形趋势(如下图中虚线上红圆点所示)选取延迟指数函数
来拟合这组数据(这将是6个方程两个未知量的线性方程组(也称超定方程组?
)):
>>t=[0.3.81.11.62.3]';
>>y=[.82.72.63.60.55.50]';
>>A=[ones(size(t)exp(-t));
>>c=A\y
c=
0.4760
0.3413
>>T=[0:
.1:
2.5]';
Y=[ones(size(T))exp(-T)]*c;
plot(t,y,'--',T,Y,'-',t,y,'o')
title('最小二乘法曲线拟合')
xlabel('\itt'),ylabel('\ity')
以上得到拟合函数
,其曲线与数据关系图为
3.当方程组有无穷多解情况(参Ch2,§2.9符号方程求解)
练习:
求解
.
解(按第2种情形,以最小二乘拟合求得:
)
>>A=[11;12;14;18;116];
>>b=[3.03,5.05,9.09,16.98,33.04]';
>>x=A\b
x=
1.0467
1.9986
5.1.4矩阵求逆
求解AX=B的方程组时,用X=A\B比用X=inv(A)*B所用内存、时间、误差小等特点更优。
pinv(A)用于计算非方阵的伪逆,例如:
>>c=[94;28;67]
c=
94
28
67
>>x=pinv(c)
x=
0.1159-0.07290.0171
-0.05340.11520.0418
则可使x*c为单位阵,即
>>q=x*c
q=
1.00000.0000
1.01.0000
但应注意c*x并非单位阵,即
>>p=c*x
p=
0.8293-0.19580.3213
-0.19580.77540.3685
0.32130.36850.3952
5.1.5LU、QR分解
若有A=LU,则可由[L,U]=lu(A)求得。
如:
>>A=[123;456;426]
A=
123
456
426
>>[L,U]=lu(A)
L=
0.2500-0.25001.0000
1.000000
1.00001.00000
U=
4.00005.00006.0000
0-3.00000
001.5000
通过LU分解,可很容易得到:
det(A)=det(L)*det(U);
inv(A)=inv(U)*inv(L).
求线性方程Ax=b时,可得到
x=U\(L\b)
这种方法运算速度更快。
正交矩阵或具有正交列的矩阵,其所有列的长度为1,且与其它列正交。
即如果Q为正交矩阵,则
Q’Q=I.
通过正交或QR分解,可将任意二维矩阵分解成一个正交阵和一个上三角阵的乘积,即A=QR。
如
>>A=[94;28;67]
A=
94
28
67
>>[Q,R]=qr(A)
Q=
-0.81820.3999-0.4131
-0.1818-0.8616-0.4739
-0.5455-0.31260.7777
R=
-11.0000-8.5455
0-7.4817
00
5.1.6矩阵求幂和矩阵指数
矩阵求幂如A2、B3可很容易求出:
A^2,B^3.元素对元素的求幂,可输入:
A.^2,B.^3.
函数dqrtm(A)可求出
。
函数expm(A)可计算出矩阵A的指数,即eA.这些函数对求解微分方程是很有用的。
例如,要求解微分方程:
其解为
因此可输入:
>>x0=[1;1;1];
>>x=[];
>>X=[];
>>fort=0:
.01:
1
X=[Xexpm(t*A)*x0];
end
plot3(X(1,:
),X(2,:
),X(3,:
),'-o')
gridon
5.1.7特征值
>>A=[0-6-1;62-16;-520-10];
>>lambda=eig(A)
lambda=
3.0710
-2.4645+17.6008i
-2.4645-17.6008i
>>[V,D]=eig(A)
V=
-0.83260.2003-0.1394i0.2003+0.1394i
-0.3553-0.2110-0.6447i-0.2110+0.6447i
-0.4248-0.6930-0.6930
D=
-3.071000
0-2.4645+17.6008i0
00-2.4645-17.6008i
5.1.8奇异值分解——P164
§5.2多项式与内插
多项式函数如roots,poly,polyval,polyvalm,residue,polyfit,polyder,conv,deconv等,它们使对多项式的操作变得便捷迅速。
5.2.1多项式表示
用行矢量表示,其元素按幂指数降序排列,例如
或表示成
p=[10–2–5];
5.2.2多项式的根
为求多项式的根,即p(x)=0的解,可利用roots函数:
>>r=roots(p)
r=
2.0946
-1.0473+1.1359i
-1.0473-1.1359i
利用poly函数可从多项式的根中恢复出多项式:
>>p2=poly(r)
p2=
1.00000-2.0000-5.0000
p=[10–2–5];
5.2.3特征多项式
poly函数还可以用于计算矩阵的特征多项式的系数:
>>A=[12-1;345;-192];
>>poly(A)
ans=1.0000-7.0000-38.000090.0000
5.2.4多项式计算
polyval函数可计算出多项式在指定点处的值,例如:
>>y1=polyval(p,4)
y1=51
5.2.5卷积和去卷积
多项式的乘和除就相当于卷积和去卷积操作,这可由函数conv和deconv实现。
例如:
设
>>a=[133];b=[456];
>>c=deconv(c,b)
c=
413282718
通过多项式除法,可以得到:
这可由deconv函数求出:
>>[q,r]=deconv(c,b)
q=
123
r=
00000
5.2.6多项式求导
polyder函数可用于求出单个多项式的导数,也可用于求两个多项式之积或之比的导数,这可由下列示例说明:
>>p=[10-2-5];
>>q=polyder(p)%求p的导数
q=30-2
>>a=[135];b=[246];
>>c=polyder(a,b)%求积a*b的导数
c=
8305638
>>[q,r]=polyder(a,b)%求商a/b的导数
q=
-2-8-2
r=
416404836
最后求得的q(s)/r(s)为a(s)/b(s)的导数。
5.2.7多项式曲线拟合
polyfit函数可在最小二乘意义下找出一多项式来拟合给定的一组数据。
其调用格式为
p=polyfit(x,y,n)
其中x,y为给定的数据,n为多项式的阶数。
例如,用三阶多项式来拟合下列数据:
>>x=[12345];y=[5.543.1128290.7498.4];
p=polyfit(x,y,3);
x2=1:
.1:
5;
y2=polyval(p,x2);
figure
(1)
subplot(2,1,1)
plot(x,y,'o',x2,y2)
gridon
title('多项式曲线拟合')
5.2.8部分分式展开
residue函数可将有理多项式进行部分分式展开:
其调用格式为
[r,p,k]=residue(b,a)
residue函数还可将部分分式形式的多项式还原成有理多项式,其调用格式为
[b,a]=residue(r,p,k).
Ex.
>>a=[100-2];b=[1-1];
>>[q,r]=deconv(a,b)
q=
111
r=
000-1
得到
>>[r,p,k]=residue(a,b)
r=
-1
p=
1
k=
111
得到
Ex.
>>a=[1];b=[1-32];
>>[r,p,k]=residue(a,b)
r=
1
-1
p=
2
1
k=
[]
5.2.9一维内插
MATLAB中提供了两类一维内插:
多项式内插和其于FFT的内插。
1.多项式内插
interp1函数可完成一维内插,其调用格式为
yi=interp1(x,y,xi,method)
其中x,y为给定的数据对,xi为要内插的点矢量,method用于指定内插方法,可取:
nearest(最邻近内插):
将内插点设置成最接近于已有数据点的值;
linear(线性内插):
连接已有数据点作线性逼近。
这是interp1函数的缺省设置;
spline(三次样条内插):
利用一系列样条函数获得内插数据点,从而确定已有数据点之间的函数;
cubic(三次曲线内插):
通过y拟合三次曲线函数,从而确定内插点的值。
所有这四种方法都要求x中的数据为单调,每种方法不要求x为均匀间隔,但如果x已经为均匀间隔,则在method之前加上*,可使执行速度加快。
这四种方法对执行速度、内存要求及得到的平滑度是不同的。
总的说来,按nearst、linear、cubic、spline顺序,内存要求从小到大,执行速度由快到慢,平漫度由差到好。
2.基于FFT的内插
interpft完成基于FFT的一维内插,其调用格式为
y=interpft(x,n)
其中x中包含等间隔取样的周期函数值,n为要求得到的点数。
5.2.10二维内插——P168
§5.3数据分析与统计
MATLAB提供的数据分析与统计函数都是面向列的,即矩阵中的每一列代表一个变量的多项观测值,其列数相应于变量数,行数相应于测量点数。
Max和min函数可求出数据的最大值和最小值,mean和std函数可求出数据的均值和标准差,sum和prod函数可求出数据元素和与数据元素积.例如,对MATLAB内含的城市24小时的车流量数据count.dat,可作分析:
>>loadcount.dat
>>mx=max(count)
mx=
114145257
>>mu=mean(count)
mu=
32.000046.541765.5833
>>sigma=std(count)
sigma=
25.370341.405768.0281
对于有些函数还可给出位置,例如,在求出最小值的同时,可得到最小值所在的位置(行号);
>>[mx,indx]=min(count)
mx=
797
indx=
22324
5.32.1协方差和协相关系数
cov函数可以求出单个变量的协方差,而corrcoef可求出两个变量之间的相关系数.例如
>>cv=cov(count)
cv=
1.0e+003*
0.64370.98021.6567
0.98021.71442.6908
1.65672.69084.6278
>>cr=corrcoef(count)
cr=
1.00000.93310.9599
0.93311.00000.9553
0.95990.95531.0000
5.3.2数据预处理
MATLAB中遇到超出范围的数据是均用NaN(非数值)表示,而且在任何运算中,只要包含NaN,就将它传递到结果中,因此在对数据进行分析前,应对数据中出现的NaN作剔除处理.例如:
>>a=[1,2,3;5NaN8;742];
>>sum(a)
ans=13NaN13
在矢量x中删除NaN元素,可有下列四种方法:
(1)i=find(~isnan(x));x=x(i);
(2)x=x(find(~isnan(x)));
(3)x=x(~isnan(x));
(4)x(isnan(x))=[];
在矩阵X中删除NaN所在的行,可输入:
X(any(isnan(x)’),:
)=[]
经过这种预处理后的数据,可进行各种分析和统计操作.
5.3.3回归和曲线拟合
对给定数据进行拟合,可采用多项式回归,也可采用其它信号形式的回归,其基本原理是最小二乘法,这在MATLAB中显得轻而易举.
例1设通过测量得到一组时间t与变量y的数据:
t=[0.3.81.11.62.3]’;
y=[.5.821.141.251.351.40]’;
分别采用多项式
进行回归,可得到两种不同的结果.MATLAB程序如下:
X1=[ones(size(t))tt.^2];
a=X1\y;%解出y1的系数
X2=[ones(size(t))exp(-t)t.*exp(-t)];
b=X2\y;%解出y2的系数
T=[0:
.1:
2.5]';%为作出拟合曲线,将区间均分,步长为0.1
Y1=[ones(size(T))TT.^2]*a;%计算每个分点的函数值,即曲线上纵坐标
Y2=[ones(size(T))exp(-T)T.*exp(-T)]*b;
figure
(1)
subplot(2,2,1)
plot(T,Y1,'-',t,y,'o'),gridon
title('多项式回归')
subplot(2,2,2)
plot(T,Y2,'-',t,y,'o'),gridon
title('指数函数回归')
例2已知变量y与x1、x2有关,测得一组数据为
>>x1=[.2.5.6.81.01.1]';
>>x2=[.1.3.4.91.11.4]';
>>y=[.17.26.28.23.27.24]';
采用
来拟合,则有
>>a=X\y
a=
0.1018
0.4844
-0.2847
因此数据的拟合模型为
5.3.4滤波——P172
5.3.5傅里叶分析FFT(快速傅里叶变换)——P173
§5.4泛函分析
MATLABA提供了可对函数进行操作的函数,它们称之为泛函,如找出函数在区间上的极小值、求函数零值点、计算函数积分等,这些都属于泛函。
5.4.1函数在MATLAB中的表示
在MATLAB中,数学函数可用M文件表示,例如函数
在MATLAB中可表示成humps.m:
functiony=humps(x)
y=1./(x-0.3).^2+0.01)+1./((x-0.9).^2+0.04)-6;
这样,humps可用作为某些函数的输入变量,从而构成泛函的计算。
5.4.2数学函数的绘图
fplot函数可绘制出指定函数的图形,它可以指定函数自变量和函数值的范围,还可在同一张图上绘制出多个图形。
例如,下面程序段执行后,可得到如下图线。
>>figure
(1)
subplot(2,2,1),fplot('humps',[-5,5]),gridon
subplot(2,2,2),fplot('humps',[-55-1025]),gridon
subplot(2,2,3),fplot('5*sin(x)',[-1,1]),gridon
subplot(2,2,4),fplot('[5*sin(x),humps(x)]',[-1,1]),gridon
5.4.3函数极小值点和零值点
fmin函数用于求出指定单变量函数在特定区间上的极小值点,fmins函数用于求出多变量函数的极小值点。
Fzero函数可求出指定函数在特定区间上零值点。
例如,要求出函数humps在[0.3,1]区间上的极小值点,可输入:
>>x=fmin('humps',0.3,1)
Warning:
FMINisobsoleteandhasbeenreplacedbyFMINBND.
FMINnowcallsFMINBNDwhichusesthefollowingsyntax:
[X,FVAL,EXITFLAG,OUTPUT]=FMINBND(FUN,x1,x2,OPTIONS,P1,P2,...)
UseOPTIMSETtodefineoptimizationoptions,ortype
'editfmin'toviewthecodeusedhere.FMINwillbe
removedinthefuture;pleaseuseFMINBNDwiththenewsyntax.
>InC:
\MATLAB6p1\toolbox\matlab\funfun\fmin.matline62
x=
1.6370
>>humps(x)%求出相应极小值
ans=
11.2528
从上图中可以看出,函数humps在[-1,0]上有一零值点,这时可用两种方式来得到:
>>a=fzero('humps',-0.2)
a=
-0.1316
或
>>a=fzero('humps',[-1,0])
a=
-0.1316
5.4.4数值积分
quad和quad8函数可以求出指定函数在指定区间上的积分。
例如,q1=quad(‘humps’,0,1)可求出humps函数在[0,1]上的积分值,当然采用quad8函数可得到相同的结果,只是采用的积分算法有所不同。
1.如quad,quad8:
>>q=quad('exp(-x)+x.^2',0,1)%注意x2应写成x.^2,即按元素计算函数对应值序列
q=
0.9655
或
>>q=quad8('exp(-x)+x.^2',0,1)
Warning:
QUAD8isobsolete.QUADLisitsrecommendedreplacement.
>InC:
\MATLAB6p1\toolbox\matlab\funfun\quad8.matline35
q=
0.9655
2.dlquad函数可以计算函数的二重积分,例如要计算
,可输入:
>>q2=dblquad('f(x,y)',x1,x2,y1,y2)
注:
书上的上面这种写法有问题。
dblquad函数的详细用法,可参看帮助文件:
>>helpdblquad。
从中找到应用的例子,再用来求解:
>>dblquad('sin(x)+y-1'),1,2,1,4)
ans=
7.3693
或先用
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- Ch51 基本 应用领域