第 5章 MATLAB数值计算
5.1 特殊矩阵
5.2 矩阵分析
5.3 矩阵分解与线性方程组求解
5.4 数据处理与多项式计算
5.5 傅立叶分析
5.6 数值微积分
5.7 常微分方程的数值求解
5.8 非线性方程的数值求解
5.9 稀疏矩阵
5.1 特殊矩阵
5.1.1对角阵与三角阵
1,矩阵的对角元素
(1)提取矩阵的对角线元素设 A为 m× n矩阵,diag(A)函数用于提取矩阵 A主对角线元素产生一个具有 min(m,n)个元素的列向量 。
diag(A)函数还有更进一步的形式 diag(A,k),其功能是提取第 k条对角线的元素 。
(2)构造对角矩阵设 V为具有 m个元素的向量,diag(V)将产生一个 m× m对角矩阵,其主对角线元素即为向量 V的元素 。
diag(V)函数也有更进一步的形式 diag(V,k),其功能是产生一个 n× n(n=m+)对角阵,其第 k条对角线的元素即为向量
V的元素 。
例 5.1 先建立 5× 5矩阵 A,然后将 A的第 1行元素乘以 1,第 2行乘以 2,…,第 5行乘以 5。
命令如下：
A=[17,0,1,0,15;23,5,7,14,16;4,0,13,0,22;10,12,19,21,3;11,18,25,2,19];
D=diag([1,2,3,4,5]);
D*A
2,矩阵的三角阵
(1)下三角矩阵求矩阵 A的下三角阵的 MATLAB函数是 tril(A)。
tril(A)函数也有更进一步的一种形式 tril(A,k)，
其功能是求矩阵 A的第 k条对角线以下的元素 。
(2)上三角矩阵在 MATLAB中,提取矩阵 A的上三角矩阵的函数是 triu(A)和 triu(A,k),其用法与提取下三角矩阵的函数 tril(A)和 tril(A,k)完全相同 。
5.1.2 特殊矩阵的生成
1,魔方矩阵函数 magic(n),其功能是生成一个 n阶魔方阵 。
例 5.2 将 101~125等 25个数填入一个 5行 5列的表格中,使其每行每列及对角线的和均为 565。
命令如下：
B=100+magic(5)
2,范得蒙矩阵函数 vander(V)生成以向量 V为基础向量的范得蒙矩阵 。
3,希尔伯特矩阵生成希尔伯特矩阵的函数是 hilb(n)。 MATLAB中,有一个专门求希尔伯特矩阵的逆的函数 invhilb(n),其功能是求 n阶的希尔伯特矩阵的逆矩阵 。
4,托普利兹矩阵生成托普利兹矩阵的函数是 toeplitz(x,y),它生成一个以 x为第 1列,y为第 1行的托普利兹矩阵 。 这里 x,y均为向量,二者不必等长 。
5,友矩阵生成友矩阵的函数是,compan(P),生成多项式 P的友矩阵 。 P是一个多项式的系数向量,高次幂系数排在前,
低次幂排在后 。
6,帕斯卡矩阵函数 pascal(n)生成一个 n阶的帕斯卡矩阵 。
例 5.3求 (x+y)5的展开式 。
在 MATLAB命令窗口,输入命令：
pascal(6)
ans =
1 1 1 1 1 1
1 2 3 4 5 6
1 3 6 10 15 21
1 4 10 20 35 56
1 5 15 35 70 126
1 6 21 56 126 252
其次对角线上的元素 1,5,10,10,5,1即为展开式的系数 。
5.2 矩阵分析
5.2.1 矩阵结构变换
1,矩阵的转置转置运算符是单撇号 (')。
2,矩阵的旋转矩阵的旋转利用函数 rot90(A,k),功能是将矩阵 A旋转 90o
的 k倍,当 k为 1时可省略 。
3,矩阵的左右翻转对矩阵 A实施左右翻转的函数是 fliplr(A)。
4,矩阵的上下翻转对矩阵 A实施上下翻转的函数是 flipud(A)。
5.2.2矩阵的逆与伪逆
1,矩阵的逆求一个矩阵的逆非常容易 。 求方阵 A的逆可调用函数 inv(A)。
例 5.4用求逆矩阵的方法解线性方程组 。
命令如下：
A=[1,2,3;1,4,9;1,8,27];b=[5,–2,6]';
x=inv(A)*b
一般情况下,用左除比求矩阵的逆的方法更有效,即 x=A\b。
2,矩阵的伪逆
MATLAB中,求一个矩阵伪逆的函数是 pinv(A)。
例 5.5 求 A的伪逆,并将结果送 B。
命令如下：
A=[3,1,1,1;1,3,1,1;1,1,3,1];
B=pinv(A)
例 5.6 求矩阵 A的伪逆 。
在 MATLAB命令窗口,输入命令：
A=[0,0,0;0,1,0;0,0,1];
pinv(A)
5.2.3 方阵的行列式求方阵 A所对应的行列式的值的函数是 det(A)。
例 5.7用克莱姆 (Cramer)方法求解线性方程组 。
程序如下：
D=[2,2,-1,1;4,3,-1,2;8,5,-3,4;3,3,-2,2]; %定义系数矩阵
b=[4;6;12;6]; %定义常数项向量
D1=[b,D(:,2:4)]; %用方程组的右端向量置换 D的第 1列
D2=[D(:,1:1),b,D(:,3:4)]; %用方程组的右端向量置换 D的第 2列
D3=[D(:,1:2),b,D(:,4:4)]; %用方程组的右端向量置换 D的第 3列
D4=[D(:,1:3),b]; %用方程组的右端向量置换 D的第 4列
DD=det(D);
x1=det(D1)/DD;
x2=det(D2)/DD;
x3=det(D3)/DD;
x4=det(D4)/DD;
[x1,x2,x3,x4]
5.2.4 矩阵的秩
MATLAB中,求矩阵秩的函数是 rank(A)。 例如,
求例 5.7中方程组系数矩阵 D的秩,命令是：
D=[2,2,-1,1;4,3,-1,2;8,5,-3,4;3,3,-2,2];
r=rank(D)
r =
4
说明 D是一个满秩矩阵 。
5.2.5 向量和矩阵的范数
1,计算向量 3种常用范数的函数
(1)norm(V)或 norm(V,2) 计算向量 V的 2—范数
(2)norm(V,1) 计算向量 V的 1—范数
(3)norm(V,inf) 计算向量 V的 ∞—范数例 5.8 已知 V,求 V的 3种范数 。
命令如下：
V=[-1,1/2,1];
v1=norm(V,1) %求 V的 1—范数
v2=norm(V) %求 V的 2—范数
vinf=norm(V,inf) %求 ∞—范数
2,矩阵的范数及其计算函数
MATLAB中提供了求 3种矩阵范数的函数,其函数调用格式与求向量的范数的函数完全相同例 5.9 求矩阵 A的三种范数 。
命令如下：
A=[17,0,1,0,15;23,5,7,14,16;4,0,13,0,22;10,12,19,21,3;11,18,25,2,19];
a1=norm(A,1) %求 A的 1—范数
a2=norm(A) %求 A的 2—范数
ainf=norm(A,inf) %求 A的 ∞—范数
5.2.6 矩阵的条件数和迹
1,的条件数
MATLAB中,计算矩阵 A的 3种条件数的函数是：
(1)cond(A,1) 计算 A的 1—范数下的条件数
(2)cond(A)或 cond(A,2) 计算 A的 2—范数数下的条件数
(3)cond(A,inf) 计算 A的 ∞—范数下的条件数例 5.10 求矩阵 X的三种条件数 。
命令如下：
A=[2,2,3;4,5,-6;7,8,9];
C1=cond(A,1)
C2=cond(A)
C3=cond(A,inf)
2,矩阵的迹
MATLAB中,求矩阵的迹的函数是 trace(A)。 例如,
X=[2 2 3;4 5 -6;7 8 9];
trace(X)
ans =
16
5.2.7 矩阵的特征值与特征向量
MATLAB中,计算矩阵 A的特征值和特征向量的函数是
eig(A),常用的调用格式有 3种：
(1)E=eig(A) 求矩阵 A的全部特征值,构成向量 E。
(2)[V,D]=eig(A) 求矩阵 A的全部特征值,构成对角阵 D，
并求 A的特征向量构成 V的列向量 。
(3)[V,D]=eig(A,'nobalance') 与第 2种格式类似,但第 2种格式中先对 A作相似变换后求矩阵 A的特征值和特征向量,
而格式 3直接求矩阵 A的特征值和特征向量 。
例 5.11 用 3种不同的格式求 A的特征值和特征向量 。
命令如下：
A=[1,2,2;1,-1,1;4,-12,1];
E=eig(A)
[V,D]=eig(A)
[V,D]=eig(A,'nobalance')
例 5.12用求特征值的方法解方程 。
命令如下：
p=[3,-7,0,5,2,-18];
A=compan(p); %A的友矩阵
x1=eig(A) %求 A的特征值
x2=roots(p) %直接多项式 p的零点两种方法求得的方程的根是完全一致的,实际上,
roots函数正是应用求友矩阵的特征值的方法来求方程的根 。
5.2.8 MATLAB在三维向量中的应用
1,向量共线或共面的判断例 5.13 设 X=(1,1,1),Y=(-1,2,1),Z=(2,2,2),判断这三个向量的共线共面问题 。
命令如下：
X=[1,1,1];Y=[-1,2,1];Z=[2,2,2];
XY=[X;Y];YZ=[Y;Z];ZX=[Z;X];XYZ=[X;Y;Z];
rank(XY)
rank(YZ)
rank(ZX)
rank(XYZ)
2,向量方向余弦的计算例 5.14设向量 V=(5,-3,2),求 V的方向余弦 。
建立一个函数文件 direct.m：
function f=f(v)
r=norm(v);
ifr==0
f=0
else
f=[v(1)/r,v(2)/r,v(3)/r];
end
return
在 MATLAB命令窗口,输入命令：
v=[5,-3,2];
f=direct(v)
3,向量的夹角例 5.15 设 U=(1,0,0),V=(0,1,0),求 U,V间的夹角 θ。
命令如下：
U=[1,0,0];V=[0,1,0];
r1=norm(U);r2=norm(V);
UV=U*V';cosd=UV/r1/r2;
D=acos(cosd)
4,两点间的距离例 5.16 设 U=(1,0,0),V=(0,1,0),求 U,V两点间的距离 。
命令如下：
U=[1,0,0];V=[0,1,0];UV=U-V;
D=norm(UV)
5,向量的向量积例 5.17设 U=(2,-3,1),V=(3,0,4),求 U× V。
命令如下：
U=[2,-3,1];V=[3,0,4];W=eye(3);
A1=[W(1,:);U;V];A2=[W(2,:);U;V];A3=[W(3,:);U;V];
UV=[det(A1),det(A2),det(A3)]
UV=
-12 -5 9
6,向量的混合积例 5.18 设 U=(0,0,2),V=(3,0,5),W=(1,1,0),求以这三个向量构成的六面体的体积 。
命令如下：
U=[0,0,2];V=[3,0,5];W=[1,1,0];
A=[U;V;W];
det(A)
ans=
6
7,点到平面的距离例 5.19求原点到平面 X+Y+Z=1的距离 。
命令如下：
u=[0,0,0];v=[1,1,1]; %
A=B=C=1,u1=u2=u3=0,D=-1
r=abs(u*v'-1)/norm(v,2)
r =
0.5774
5.3 矩阵分解与线性方程组求解
5.3.1矩阵分解
1,实对称矩阵的 QDQ分解例 5.20设对称矩阵 A,对 A进行 QDQ分解 。
命令如下：
A=[2,1,4,6;1,2,1,5;4,1,3,4;6,5,4,2];
[Q,D]=eig(A)
Q*D*Q'
ans =
2.0000 1.0000 4.0000 6.0000
1.0000 2.0000 1.0000 5.0000
4.0000 1.0000 3.0000 4.0000
6.0000 5.0000 4.0000 2.0000
结果与 A相等,说明确实将 A分解为了 QDQ'的乘积 。
例 5.21求下列二次型的标准形式及变换矩阵 。
命令如下：
A=[1,2,1;2,1,1;1,1,3;];
[Q,D]=eig(A)
进一步作线性变换即得关于 u,v,w的标准二次型：
2,矩阵的 LU分解
MATLAB中,完成 LU分解的函数是：
(1)[L,U]=lu(A) 将方阵 A分解为交换下三角矩阵 L
和上三角矩阵 U,使 A=LU。
(2)[L,U,P]=lu(A) 将方阵 A分解为下三角矩阵 L和上三角矩阵 U,使 PA=LU。
例 5.22用 LU分解求方程组的根 。
3,矩阵的 QR分解对矩阵 A进行 QR分解的函数是 [Q,R]=qr(A),根据方阵 A,求一个正交矩阵 Q和一个上三角矩阵 R，
使 A=Q*R。 例如,对矩阵 A进行 QR分解的命令是：
A=[2,1,-2;1,2,1;2,5,3];
[Q,R]=qr(A)
5.3.2 线性方程组求解
1,线性方程组解的一般讨论解线性方程组的一般函数文件如下：
function [x,y]=line_solution(A,b)
[m,n]=size(A);y=[];
if norm(b)>0 %非齐次方程组
if rank(A)==rank([a,b]) %方程组相容
if rank(A)==m %有唯一解
x=A\b;
else %方程组有无穷多个解,基础解系
disp('原方程组有有无穷个解,其齐次方程组的基础解系为 y,特解为 x');
y=null(A,'r');
x=A\b;
end
else %方程组不相容,给出最小二乘法解
disp('方程组的最小二乘法解是,');
x=A\b;
end
else %齐次方程组
if rank(A)>=n %列满秩
x=zero(m,1) %0解
else %非 0解
disp('方程组有无穷个解,基础解系为 x');
x=null(A,'r');
2,应用举例例 5.23求线性方程组的解 。
在 MATLAB命令窗口,输入命令：
A=[2,2,-1,1;4,3,-1,2;8,5,-3,4;3,3,-2,2];b=[4,6,12,6]';
[x,y]=line_solution(A,b) %调用自定义函数例 5.24求下列线性方程组的解 。
在 MATLAB命令窗口,输入命令：
A=[2,7,3,1;3,5,2,2;9,4,1,7];b=[6,4,2]';
[x,y]=line_solution(A,b)
5.4 数据处理与多项式计算
5.4.1 数据统计与分析
1,求矩阵最大和最小元素
(1)求向量的最大最小元素
① y=max(X) 返回向量 X的最大元素存入 y。
② [y,I]=max(X) 返回向量 X的最大元素存入 y,最大元素的序号存入 I。
(2)求矩阵的最大和最小元素
① max(A) 返回一个行向量,向量的第 i个元素是 A矩阵的第 i列上的最大元素 。
② [Y,U]=max(A) 返回两个行向量,Y向量记录 A的每列的最大元素,U向量记录每列最大元素的行号 。
③ max(A,[],dim) dim取 1或 2。 dim取 1时,该函数和 max(A)
完全相同 。 dim取 2时,该函数返回一个列向量,其第 i
个元素是 A矩阵的第 i行上的最大元素 。
(3)两个向量或矩阵对应元素的比较
① U=max(A,B) A,B是两个同型的向量或矩阵 。 结果 U是与 A,B同型的向量或矩阵,U的每个元素等于 A,B对应元素的较大者 。
② U=max(A,n) n是一个标量 。 结果 U是与 A同型的向量或矩阵,U的每个元素等于 A对应元素和 n中的较大者 。
min函数的用法和 max完全相同 。
例 5.25 求矩阵 A的每行及每列的最大和最小元素,并求整个矩阵的最大和最小元 。
命令如下：
A=[13,-56,78;25,63,-235;78,25,563;1,0,-1];
max(A,[],2) %求每行最大元素
min(A,[],2) %求每行最小元素
max(A) %求每列最大元素
min(A) %求每列最小元素
max(max(A)) %求整个矩阵的最大元素
min(min(A)) %求整个矩阵的最小元素
2,求矩阵的平均值和中值求矩阵和向量元素的平均值的函数是 mean,求中值的函数是 median。 它们的调用方法和 max函数完全相同 。
3,矩阵元素求和与求积矩阵和向量求和与求积的基本函数是 sum和 prod,其使用方法和 max类似 。
例 5.26求矩阵 A的每行元素的乘积和全部元素的乘积 。
命令如下：
A=[1,2,3,4;5,6,7,8;9,10,11,12];
S=prod(A,2)
prod(S) %求 A的全部元素的乘积
4,矩阵元素累加和与累乘积
MATLAB中,使用 cumsum和 cumprod函数能方便地求得向量和矩阵元素的累加和与累乘积向量,函数的用法和 sum及 prod相同例 5.27求向量 X=(1！,2！,3！,…,10！ )。
命令如下：
X=cumprod(1:10)
5,标准方差
MATLAB中,提供了计算数据序列的标准方差的函数 std。
对于向量 X,std(X)返回一个标准方差 。 对于矩阵 A，
std(A)返回一个行向量,它的各个元素便是矩阵 A各列或各行的标准方差 。 std函数的一般调用格式为：
std(A,FLAG,dim)
其中 dim取 1或 2。 当 dim=1时,求各列元素的标准方差；
当 dim=2时,则求各行元素的标准方差 。 FLAG取 0或 1。
6,元素排序
MATLAB中对向量 X是排序函数是 sort(X),函数返回一个对 X中的元素按升序排列的新向量 。
sort函数也可以对矩阵 A的各列 (或行 )重新排序,
其调用格式为：
[Y,I]=sort(A,dim)
其中 dim指明对 A的列还是行进行排序,若 dim=1，
则按列排,若 dim=2,则按行排 。 Y是排序后的矩阵,而 I记录 Y中的元素在 A中位置 。
例 5.28对矩阵做各种排序 。
命令如下：
A=[1,-8,5;4,12,6;13,7,-13];
sort(A) %对 A的每列按升序排序
-sort(-A,2) %对 A的每行按降序排序
[X,I]=sort(A) %对 A按列排序,并将每个元素所在行号送矩阵 I
5.4.2 数值插值
1,一维数值插值
interp1函数调用格式为：
Y1=interp1(X,Y,X1,'method')
函数根据 X,Y的值,计算函数在 X1处的值 。 X,Y是两个等长的已知向量,分别描述采样点和样本值,X1是一个向量或标量,描述欲插值的点,Y1是一个与 X1等长的插值结果 。 method是插值方法,允许的取值有
'linear'(线性插值 ),'nearest'(最近插值 ),'spline'(三次样条插值 ),'cubic'（ 三次多项式插值 ）,缺省值是 'linear'。
例 5.29用不同的插值方法计算 sin(x)在 π/2点的值 。
这是一个一维插值问题 。 在 MATLAB命令窗口,输入命令：
X=0:0.2:pi;Y=sin(X); %给出 X,Y
interp1(X,Y,pi/2) %用缺省方法 (即线性插值方法 )计算
sin(π/2)
interp1(X,Y,pi/2,'nearest') %用最近方法计算 sin(π/2)
interp1(X,Y,pi/2,'linear') %用线性方法计算 sin(π/2)
interp1(X,Y,pi/2,'spline') %用三次样条方法计算 sin(π/2)
interp1(X,Y,pi/2,'cubic') %用三次多项式方法计算 sin(π/2)
MATLAB 中 有 一 个 专 门 的 三 次 样 条 插 值 函 数
Y1=spline(X,Y,X1),其功能及使用方法与函数
Y1=interp1(X,Y,X1,'spline')完全相同 。
例 5.30 已知检测参数 f随时间 t的采样结果,用数值插值法计算
t=2,7,12,17,22,17,32,37,42,47,52,57时 f的值 。
这是一个一维数值插值问题,命令如下：
T=0:5:65;
X=2:5:57;
F=[3.2015,2.2560,879.5,1835.9,2968.8,4136.2,5237.9,6152.7,...
6725.3,6848.3,6403.5,6824.7,7328.5,7857.6];
F1=interp1(T,F,X) %用线性方法插值
F1=interp1(T,F,X,'nearest') %用最近方法插值
F1=interp1(T,F,X,'spline') %用三次样条方法插值
F1=interp1(T,F,X,'cubic') %用三次多项式方法插值
2,二维数值插值
MATLAB中,提供了解决二维插值问题的函数 。
其调用格式为：
Z1=interp2(X,Y,Z,X1,Y1,'method')
其中 X,Y是两个向量,分别描述两个参数的采样点,Z是与参数采样点对应的采样变量的样本值,X1,Y1是两个向量或标量,描述欲插值的点 。 method的取值与一维插值函数相同 。
例 5.31设 Z=x2+y2,对 Z函数在 (0,1)× (0,2)区域内进行插值 。
命令如下：
x=0:0.1:10;y=0:0.2:20;
[X,Y]=meshgrid(x,y);
Z=X.^2+Y.^2;
interp2(x,y,Z,0.5,0.5) %对函数在 (0.5,0.5)点进行插值
interp2(x,y,Z,[0.5 0.6],0.4) %对函数在 (0.5,0.4)点和 (0.6,0.4)点进行插值
interp2(x,y,Z,[0.5 0.6],[0.4 0.5]) %对函数在 (0.5,0.4)点和 (0.6,0.5)点进行插值
interp2(x,y,Z,[0.5 0.6]',[0.4 0.5])%对函数在 (0.5,0.4),(0.6,0.4),(0.5,0.5)和 (0.6,0.5)点进行插值
3,三维数值插值对三维函数插值的函数是 interp3,其使用方法和 interp2相同 。 其调用格式为：
W1=interp3(X,Y,Z,W,X1,Y1,Z1,'method')
函数返回三维插值结果 。 其中 X,Y,Z是三个向量,分别描述三个参数的采样点,W是与参数采样点对应的采样变量的样本值,X1,Y1,Z1是三个向量或标量,
描述欲插值的点 。 method是插值方法,可选,其缺省值是 ‘ line'。 method的取值与一,二维插值函数相同 。
5.4.3 曲线拟合
MATLAB中,提供了解决使用最小二乘法进行曲线拟合的函数 。 调用格式为：
[P,S]=polyfit(X,Y,m)
函数根据采样点 X和采样点函数值 Y,产生一个 m
次多项式 P及其在采样点的误差向量 S。
其中 X,Y是两个等长的向量,P是一个长度为
m+1的向量 。
例 5.32 用一个 5次多项式在区间 [0,2π]内逼近函数 sin(x)。
命令如下：
X=linspace(0,2*pi,50);Y=sin(X);
[P,S]=polyfit(X,Y,5) %得到 5次多项式的系数和误差
plot(X,Y,'k*',X,polyval(P,X),'k-')
5.4.4 多项式计算
1,多项式的建立已知一个多项式的全部根 X求多项式系数的函数是 poly(X)，
该函数返回以 X为全部根的一个多项式 P,当 X是一个长度为 m的向量时,P是一个长度为 m+1的向量 。
2,多项式求根求多项式 p(x)的根的函数是 roots(P),这里,P是 p(x)的系数向量,该函数返回方程 p(x)=0的全部根 (含重根,复根 )。
3,多项式求值求多项式 p(x)在某点或某些点的函数值的函数是 polyval(P,x)。
若 x为一数值,则求多项式在该点的值；若 x为向量或矩阵,
则对向量或矩阵中的每个元素求其多项式的值 。
例 5.33 已知一个多项式,计算：
(1)计算 f(x)=0 的全部根 。
(2)由方程 f(x)=0的根构造一个多项式 g(x),并与 f(x)进行对比 。
(3)计算 f(5),f(7.8),f(9.6),f(12.3)的值 。
命令如下：
P=[3,0,4,-5,-7.2,5];
X=roots(P) %求方程 f(x)=0的根
G=poly(X) %求多项式 g(x)
X0=[5,7.8,9.6,12.3];
f=polyval(P,X0) %求多项式 f(x)在给定点的值多项式求值还有一个函数是 polyvalm,其调用格式与 polyval
相同,但含义不同 。 polyvalm函数要求 x为方阵,它以方阵为自变量求多项式的值 。
4,多项式的四则运算
(1)多项式的加减法
(2)多项式的乘法函数 conv(P1,P2)用于求多项式 P1和 P2的乘积 。
(3)多项式的除法函数 [Q,r]=deconv(P1,P2)用于对多项式 P1和 P2作除法运算 。 其中 Q返回多项式 P1除以 P2的商式,r返回 P1除以 P2的余式 。 这里,Q和 r仍是多项式系数向量 。
deconv是 conv的逆函数,即有 P1=conv(P2,Q)+r。
例 5.34设有两个多项式,计算：
(1)求 f(x)+g(x),f(x)-g(x)。
(2)求 f(x)·g(x),f(x)/g(x)。
在 MATLAB命令窗口,输入命令：
f=[3,-5,2,-7,5,6];g=[3,5,-3];g1=[0,0,0,g];
f+g1 %求 f(x)+g(x)
f-g1 %求 f(x)-g(x)
conv(f,g) %求 f(x)*g(x)
[Q,r]=deconv(f,g) %求 f(x)/g(x),商式送 Q,余式送 r。
5,多项式的导函数对多项式求导数的函数是：
p=polyder(P) 求多项式 P的导函数
p=polyder(P,Q) 求 P*Q的导函数
[p,q]=polyder(P,Q) 求 P/Q的导函数,导函数的分子存入 p，
分母存入 q。
例 5.35求有理分式的导数 。
命令如下：
P=[3,5,0,-8,1,-5];
Q=[10,5,0,0,6,0,0,7,-1,0,-100];
[p,q]=polyder(P,Q)
5.4.5 函数的最大值与最小值
MATLAB中用于求最小值的函数是：
fmin(f,a,b) 求单变量函数 f(x)在区间 (a,b)上的最小值点 。
fmins(F,X0) 求多变量函数 F(x)在估计值 X0附近的最小值点 。
MATLAB没有专门提供求函数最大值点的函数,但只要注意到 -f(x)在区间 (a,b)上的最小值点就是 f(x)在 (a,b)的最大值点,所以 fmin(-f,a,b)返回函数 f(x)在区间 (a,b)上的最大值 。
例 5.36 求函数 f(x)在区间 (-10,1)和 (1,10)上的最小值点 。
首先建立函数文件 fx.m：
function f=f(x)
f=x-1/x+5;
return
再在 MATLAB命令窗口,输入命令：
fmin('fx',-10,-1) %求函数在区间 (-10,-1)内的最小值点
fmin(f,1,10) %求函数在区间 (1,10)内的最小值点 。 注意函数名 f不用加 '
例 5.37 设有函数 f(x,y,z),求函数 f在 (0.5,0.5,0.5)附近的最小值 。
建立函数文件 fxyz.m：
function f=f(u)
x=u(1);y=u(2);z=u(3);
f=x+y.^2./x/4+z.^2./y+2./z;
return
在 MALAB命令窗口,输入命令：
U=fmins('fxyz',[0.5,0.5,0.5]) %求函数的最小值点
fxyz(U) %求函数的最小值
5.5 傅立叶分析
MATLAB中,提供了对向量 (或直接对矩阵的行或列 )进行离散傅立叶变换的函数,其调用格式是：
Y=fft(X,n,dim)
(1)当 X是一个向量时,返回对 X的离散傅立叶变换 。
(2)当 X是一个矩阵时,返回一个矩阵并送 Y,其列 (行 )是对 X的列 (行 )的离散傅立叶变换 。
例 5.38 求 X=(1,0,-3,5,2)的离散傅立叶逆变换 。
在 MATLAB命令窗口,输入命令：
X=[1,0,-3,5,2];
Y=fft(X) %对 X进行变换
3,离散傅立叶变换的逆变换
MATLAB中,对向量 (或直接对矩阵的行或列 )进行离散傅立叶逆变换的函数的调用方法是：
Y=ifft(X,n,dim)
函数对 X进行离散傅立叶逆变换 。 其中 X,n,dim的意义及用法和离散傅立叶变换函数 fft完全相同 。
例 5.39 对矩阵 A的列向量,行向量分别进行离散傅立叶变换,并对变换结果进行逆变换 。
命令如下：
A=[3,2,1,1;-5,1,0,1;3,2,1,5];
fftA=fft(A) %求 A的列向量的傅立叶变换
fftA2=fft(A,4,2) %求 A的行向量的傅立叶变换
ifft(fftA) %对矩阵 fftA的列向量进行傅立叶逆变换,
结果应等于 A
ifft(fftA2,4,2) %对矩阵 fftA2的行向量进行傅立叶逆变换,其结果应等于 A
5.6 数值微积分
5.6.1 数值微分
MATLAB中,没有直接提供求数值导数的函数,只有计算向前差分的函数 。
DX=diff(X) 计算向量 X的向前差分,DX(i)=X(i+1)-X(i)，
0<i<n。
DX=diff(X,n) 计算 X 的 n 阶 向 前 差 分,
diff(X,2)=diff(diff(X))。
DX=diff(A,n,dim) 计算矩阵 A的 n阶差分,dim=1时 (缺省状态 ),按列计算差分,dim=2,按行计算差分 。
例 5.40 求向量 sin(X)的 1~3阶差分 。 设 X由 [0,2π]间均匀分布的 10个点组成 。
命令如下：
X=linspace(0,2*pi,10);
Y=sin(X);
DY=diff(Y); %计算 Y的一阶差分
D2Y=diff(Y,2); %计算 Y的二阶差分,也可用命令
diff(DY)计算
D3Y=diff(Y,3); %计算 Y的三阶差分,也可用 diff(D2Y)
或 diff(DY,2)
例 5.41 用不同的方法求函数 f(x)的数值导数,并在同一个坐标系中做出 f'(x)的图象 。
程序如下：
f=inline('sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2');
g=inline('(3*x.^2+4*x-1)./sqrt(x.^3+2*x.^2-x+12)/2+1/6./(x+5).^(5/6)+5');
x=-3:0.01:3;
p=polyfit(x,f(x),5); %用 5次多项式 p拟合 f(x)
dp=polyder(p); %对拟合多项式 p求导数 dp
dpx=polyval(dp,x); %求 dp在假设点的函数值
dx=diff(f([x,3.01]))/0.01; %直接对 f(x)求数值导数
gx=g(x); %求函数 f的导函数 g在假设点的导数
plot(x,dpx,x,dx,'g.',x,gx,'r-'); %作图
5.6.2数值积分
(1)被积函数是一个解析式函数 quad(f,a,b,tol,trace)用于求被积函数 f(x)在 [a,b]上的定积分,tol是计算精度,缺省值是 0.001。 trace非 0时,画出积分图形 。 注意,调用 quad函数时,先要建立一个描述被积函数 f(x)的函数文件或语句函数 。 当被积函数
f含有一个以上的变量时,quad函数的调用格式为：
quad(f,a,b,tol,trace,g1,g2)
其中 f,a,b,tol,trace等参数的含义同前 。
数值积分函数还有一种形式 quad8,其用法与 quad完全相同 。
例 5.42 用两种不同的方法求积分 。
先建立一个函数文件 ex.m：
function ex=ex(x)
ex=exp(-x.^2); %注意应用点运算
return
然后,在 MATLAB命令窗口,输入命令：
quad('ex',0,1,1e-6) %注意函数名应加字符引号
quad8('ex',0,1,1e-6) %用另一函数求积分例 5.43用 trapz函数计算积分 。
在 MATLAB命令窗口,输入命令：
X=0:0.01:1;Y=exp(-X.^2);
trapz(X,Y)
(2)被积函数由一个表格定义
MATLAB中,对由表格形式定义的函数关系的求定积分问题用 trapz(X,Y)函数 。 其中向量 X,Y
定义函数关系 Y=f(X)。
(3)二重积分例 5.44计算二重积分 。
建立一个函数文件 fixy.m：
functionf=f(x,y)
f=exp(-x.^2-y.^2);
return
建立一个命令文件 ftxy1.m：
for i=1:20
int2(i)=quad('fixy',0,1,[],[],x(i)); %在二维函数 fixy中以 x=x(i)代入并对 y积分 。
end
在 MATLAB命令窗口,输入命令：
x=linspace(0,1,20);
ftxy1
trapz(x,int2)
实际上,MATLAB提供了计算二重积分的函数：
dblquad(f,a,b,c,d,tol,trace)
该函数求 f(x,y)在 [a,b]× [c,d]区域上的二重积分 。
参数 tol,trace的用法与函数 quad完全相同 。
如果直接使用这里介绍的二重积分函数 dblquad来求解本例就非常简单,命令如下：
g=inline('exp(-x.^2-y.^2)');
dblquad(g,0,1,0,1) %直接调用二重积分函数求解
5.7 常微分方程的数值求解基于龙格－库塔法,MATLAB提供了求常微分方程数值解的函数,一般调用格式为：
[X,Y]=ode23(f,[x0,xn],y0)
[X,Y]=ode45(f,[x0,xn],y0)
其中 X,Y是两个向量,X对应自变量 x在求解区间 [x1，
xn]的一组采样点,其采样密度是自适应的,无需指定；
Y是与 X对应的一组解,f是一个函数,[x0,xn]代表自变量的求解区间,y0=y(x0),由方程的初值给定 。 函数在求解区间 [x0,xn]内,自动设立采样点向量 X,并求出解函数 y在采样点 X处的样本值 。
例 5.45 求微分方程初值问题在 [1,3]区间内的数值解,并将结果与解析解进行比较 。
先建立一个该函数的 m文件 fxy1.m：
function f=f(x,y)
f=-2.*y./x+4*x %注意使用点运算符
return
再输入命令：
[X,Y]=ode45('fxy1',[1,3],2);
X' %显示自变量的一组采样点
Y' %显示求解函数与采样点对应的一组数值解
(X.^2+1./X.^2)' %显示求解函数与采样点对应的一组解析解例 5.46 求解初值问题在区间 [0,2]中的解 。
建立一个函数文件 fxy2.m：
function f=f(x,y)
f(2)=-x.*y(2)+x.^2-5;
f(1)=y(2);
f=f';
return
在 MATLAB命令窗口,输入命令：
[X,Y]=ode45('fxy2',[0,2],[5,6]);
[X,Y]
5.8 非线性方程的数值求解
1,单变量非线性方程求解
MATLAB中,提供了求解单变量方程的函数 fzero(f,x0,tol)，
该函数采用迭代法计算函数 f(x)的一个零点,迭代初值为 x0,当两次迭代结果小于 tol时停止迭代过程 。 tol的缺省值是 eps。
注意,在调用函数 fzero 之前,要使用 m文件建立自己要计算的函数 f(x),只有定义了函数 f(x)的 m文件后,才能在 fzero函数的参数中使用自定义函数名 。
例 5.47 求 f(x)=x-+5 在 x0=-5和 x0=1作为迭代初值时的零点 。
先编制一个函数文件 fz.m：
function f=f(x)
f=x-1/x+5;
然后,在 MATLAB命令窗口,输入命令：
fzero('fz',-5) %以 -5作为迭代初值
Zero found in the interval,[-4.8,-5.2].
fzero('fz',1)
2,非线性方程组求解函数 fsolve调用格式为：
X=fsolve(F,X0)
例 5.48 求方程组在 (1,1,1)附近的解并对结果进行验证 。
首先建立方程的函数文件 fxyz1.m：
function F=F(X)
x=X(1);y=X(2);z=X(3);
F(1)=sin(x)+y+z^2*exp(x);
F(2)=x+y*z;
F(3)=x*y*z;
在 MATLAB命令窗口,输入命令：
X=fsolve('fxyz1',[1,1,1]) %求解 X的三个分量 x,y,z
Y=fxyz1(X) %检验所求结果 X是否满足原方程组
norm(Y) %求 Y向量的模例 5.49 求圆和直线的两个交点 。
建立方程组函数文件 fxyz2.m：
function F=F(X)
x=X(1);y=X(2);z=X(3);
F(1)=x^2+y^2+z^2-9;
F(2)=3*x+5*y+6*z;
F(3)=x-3*y-6*z-1;
在 MATLAB命令窗口,输入命令：
X1=fsolve('fxyz2',[-1,1,-1]) %求直线与球面的第一个交点
X2=fsolve('fxyz2',[1,-1,1]) %求直线与球面的第二个交点
5.9 稀疏矩阵
5.9.1 矩阵存储方式
1,矩阵的完全存储模式
2,稀疏矩阵的存储方式
5.9.2 稀疏存储方式的产生与转化
1.将一个完全存储方式的转化为稀疏存储方式函数 B=sparse(A)将矩阵 A转化为稀疏存储方式的矩阵 B。
sparse函数还有其他一些格式：
sparse(m,n) 生成一个 m× n的所有元素都是 0的稀疏矩阵 。
sparse(u,v,S) u,v,S是三个等长的向量 。
此外,还有一些和稀疏矩阵操作有关的函数 。 例如
[U,V,S]=find(A) 返回矩阵 A中非 0元素的下标和元素 。 这里产生的 U,V、
S可作为 sparse(u,v,s)的参数 。
full(A) 返回和稀疏存储矩阵 A对应的完全存储方式矩阵 。
2,产生一个稀疏矩阵把要建立的稀疏矩阵的非 0元素及其所在行和列的位置表示出来后由 MATLAB自己产生其稀疏存储方式,这需要使用 spconvert函数 。 调用格式为：
B=spconvert(A)
其中 A为一个 m× 3或 m× 4的矩阵,其每行表示一个非 0元素,m是非 0元素的个数 。
3,单位稀疏矩阵的产生单位矩阵只有对角线元素为 1,其他元素都为 0,是一种具有稀疏特征的矩阵 。 我们知道,函数 eye产生一个完全存储方式的单位矩阵 。 MATLAB还有一个产生稀疏存储方式的单位矩阵的函数,这就是 speye。 函数
speye(m,n)返回一个 m× n的稀疏存储单位矩阵 。
5.1 特殊矩阵
5.2 矩阵分析
5.3 矩阵分解与线性方程组求解
5.4 数据处理与多项式计算
5.5 傅立叶分析
5.6 数值微积分
5.7 常微分方程的数值求解
5.8 非线性方程的数值求解
5.9 稀疏矩阵
5.1 特殊矩阵
5.1.1对角阵与三角阵
1,矩阵的对角元素
(1)提取矩阵的对角线元素设 A为 m× n矩阵,diag(A)函数用于提取矩阵 A主对角线元素产生一个具有 min(m,n)个元素的列向量 。
diag(A)函数还有更进一步的形式 diag(A,k),其功能是提取第 k条对角线的元素 。
(2)构造对角矩阵设 V为具有 m个元素的向量,diag(V)将产生一个 m× m对角矩阵,其主对角线元素即为向量 V的元素 。
diag(V)函数也有更进一步的形式 diag(V,k),其功能是产生一个 n× n(n=m+)对角阵,其第 k条对角线的元素即为向量
V的元素 。
例 5.1 先建立 5× 5矩阵 A,然后将 A的第 1行元素乘以 1,第 2行乘以 2,…,第 5行乘以 5。
命令如下：
A=[17,0,1,0,15;23,5,7,14,16;4,0,13,0,22;10,12,19,21,3;11,18,25,2,19];
D=diag([1,2,3,4,5]);
D*A
2,矩阵的三角阵
(1)下三角矩阵求矩阵 A的下三角阵的 MATLAB函数是 tril(A)。
tril(A)函数也有更进一步的一种形式 tril(A,k)，
其功能是求矩阵 A的第 k条对角线以下的元素 。
(2)上三角矩阵在 MATLAB中,提取矩阵 A的上三角矩阵的函数是 triu(A)和 triu(A,k),其用法与提取下三角矩阵的函数 tril(A)和 tril(A,k)完全相同 。
5.1.2 特殊矩阵的生成
1,魔方矩阵函数 magic(n),其功能是生成一个 n阶魔方阵 。
例 5.2 将 101~125等 25个数填入一个 5行 5列的表格中,使其每行每列及对角线的和均为 565。
命令如下：
B=100+magic(5)
2,范得蒙矩阵函数 vander(V)生成以向量 V为基础向量的范得蒙矩阵 。
3,希尔伯特矩阵生成希尔伯特矩阵的函数是 hilb(n)。 MATLAB中,有一个专门求希尔伯特矩阵的逆的函数 invhilb(n),其功能是求 n阶的希尔伯特矩阵的逆矩阵 。
4,托普利兹矩阵生成托普利兹矩阵的函数是 toeplitz(x,y),它生成一个以 x为第 1列,y为第 1行的托普利兹矩阵 。 这里 x,y均为向量,二者不必等长 。
5,友矩阵生成友矩阵的函数是,compan(P),生成多项式 P的友矩阵 。 P是一个多项式的系数向量,高次幂系数排在前,
低次幂排在后 。
6,帕斯卡矩阵函数 pascal(n)生成一个 n阶的帕斯卡矩阵 。
例 5.3求 (x+y)5的展开式 。
在 MATLAB命令窗口,输入命令：
pascal(6)
ans =
1 1 1 1 1 1
1 2 3 4 5 6
1 3 6 10 15 21
1 4 10 20 35 56
1 5 15 35 70 126
1 6 21 56 126 252
其次对角线上的元素 1,5,10,10,5,1即为展开式的系数 。
5.2 矩阵分析
5.2.1 矩阵结构变换
1,矩阵的转置转置运算符是单撇号 (')。
2,矩阵的旋转矩阵的旋转利用函数 rot90(A,k),功能是将矩阵 A旋转 90o
的 k倍,当 k为 1时可省略 。
3,矩阵的左右翻转对矩阵 A实施左右翻转的函数是 fliplr(A)。
4,矩阵的上下翻转对矩阵 A实施上下翻转的函数是 flipud(A)。
5.2.2矩阵的逆与伪逆
1,矩阵的逆求一个矩阵的逆非常容易 。 求方阵 A的逆可调用函数 inv(A)。
例 5.4用求逆矩阵的方法解线性方程组 。
命令如下：
A=[1,2,3;1,4,9;1,8,27];b=[5,–2,6]';
x=inv(A)*b
一般情况下,用左除比求矩阵的逆的方法更有效,即 x=A\b。
2,矩阵的伪逆
MATLAB中,求一个矩阵伪逆的函数是 pinv(A)。
例 5.5 求 A的伪逆,并将结果送 B。
命令如下：
A=[3,1,1,1;1,3,1,1;1,1,3,1];
B=pinv(A)
例 5.6 求矩阵 A的伪逆 。
在 MATLAB命令窗口,输入命令：
A=[0,0,0;0,1,0;0,0,1];
pinv(A)
5.2.3 方阵的行列式求方阵 A所对应的行列式的值的函数是 det(A)。
例 5.7用克莱姆 (Cramer)方法求解线性方程组 。
程序如下：
D=[2,2,-1,1;4,3,-1,2;8,5,-3,4;3,3,-2,2]; %定义系数矩阵
b=[4;6;12;6]; %定义常数项向量
D1=[b,D(:,2:4)]; %用方程组的右端向量置换 D的第 1列
D2=[D(:,1:1),b,D(:,3:4)]; %用方程组的右端向量置换 D的第 2列
D3=[D(:,1:2),b,D(:,4:4)]; %用方程组的右端向量置换 D的第 3列
D4=[D(:,1:3),b]; %用方程组的右端向量置换 D的第 4列
DD=det(D);
x1=det(D1)/DD;
x2=det(D2)/DD;
x3=det(D3)/DD;
x4=det(D4)/DD;
[x1,x2,x3,x4]
5.2.4 矩阵的秩
MATLAB中,求矩阵秩的函数是 rank(A)。 例如,
求例 5.7中方程组系数矩阵 D的秩,命令是：
D=[2,2,-1,1;4,3,-1,2;8,5,-3,4;3,3,-2,2];
r=rank(D)
r =
4
说明 D是一个满秩矩阵 。
5.2.5 向量和矩阵的范数
1,计算向量 3种常用范数的函数
(1)norm(V)或 norm(V,2) 计算向量 V的 2—范数
(2)norm(V,1) 计算向量 V的 1—范数
(3)norm(V,inf) 计算向量 V的 ∞—范数例 5.8 已知 V,求 V的 3种范数 。
命令如下：
V=[-1,1/2,1];
v1=norm(V,1) %求 V的 1—范数
v2=norm(V) %求 V的 2—范数
vinf=norm(V,inf) %求 ∞—范数
2,矩阵的范数及其计算函数
MATLAB中提供了求 3种矩阵范数的函数,其函数调用格式与求向量的范数的函数完全相同例 5.9 求矩阵 A的三种范数 。
命令如下：
A=[17,0,1,0,15;23,5,7,14,16;4,0,13,0,22;10,12,19,21,3;11,18,25,2,19];
a1=norm(A,1) %求 A的 1—范数
a2=norm(A) %求 A的 2—范数
ainf=norm(A,inf) %求 A的 ∞—范数
5.2.6 矩阵的条件数和迹
1,的条件数
MATLAB中,计算矩阵 A的 3种条件数的函数是：
(1)cond(A,1) 计算 A的 1—范数下的条件数
(2)cond(A)或 cond(A,2) 计算 A的 2—范数数下的条件数
(3)cond(A,inf) 计算 A的 ∞—范数下的条件数例 5.10 求矩阵 X的三种条件数 。
命令如下：
A=[2,2,3;4,5,-6;7,8,9];
C1=cond(A,1)
C2=cond(A)
C3=cond(A,inf)
2,矩阵的迹
MATLAB中,求矩阵的迹的函数是 trace(A)。 例如,
X=[2 2 3;4 5 -6;7 8 9];
trace(X)
ans =
16
5.2.7 矩阵的特征值与特征向量
MATLAB中,计算矩阵 A的特征值和特征向量的函数是
eig(A),常用的调用格式有 3种：
(1)E=eig(A) 求矩阵 A的全部特征值,构成向量 E。
(2)[V,D]=eig(A) 求矩阵 A的全部特征值,构成对角阵 D，
并求 A的特征向量构成 V的列向量 。
(3)[V,D]=eig(A,'nobalance') 与第 2种格式类似,但第 2种格式中先对 A作相似变换后求矩阵 A的特征值和特征向量,
而格式 3直接求矩阵 A的特征值和特征向量 。
例 5.11 用 3种不同的格式求 A的特征值和特征向量 。
命令如下：
A=[1,2,2;1,-1,1;4,-12,1];
E=eig(A)
[V,D]=eig(A)
[V,D]=eig(A,'nobalance')
例 5.12用求特征值的方法解方程 。
命令如下：
p=[3,-7,0,5,2,-18];
A=compan(p); %A的友矩阵
x1=eig(A) %求 A的特征值
x2=roots(p) %直接多项式 p的零点两种方法求得的方程的根是完全一致的,实际上,
roots函数正是应用求友矩阵的特征值的方法来求方程的根 。
5.2.8 MATLAB在三维向量中的应用
1,向量共线或共面的判断例 5.13 设 X=(1,1,1),Y=(-1,2,1),Z=(2,2,2),判断这三个向量的共线共面问题 。
命令如下：
X=[1,1,1];Y=[-1,2,1];Z=[2,2,2];
XY=[X;Y];YZ=[Y;Z];ZX=[Z;X];XYZ=[X;Y;Z];
rank(XY)
rank(YZ)
rank(ZX)
rank(XYZ)
2,向量方向余弦的计算例 5.14设向量 V=(5,-3,2),求 V的方向余弦 。
建立一个函数文件 direct.m：
function f=f(v)
r=norm(v);
ifr==0
f=0
else
f=[v(1)/r,v(2)/r,v(3)/r];
end
return
在 MATLAB命令窗口,输入命令：
v=[5,-3,2];
f=direct(v)
3,向量的夹角例 5.15 设 U=(1,0,0),V=(0,1,0),求 U,V间的夹角 θ。
命令如下：
U=[1,0,0];V=[0,1,0];
r1=norm(U);r2=norm(V);
UV=U*V';cosd=UV/r1/r2;
D=acos(cosd)
4,两点间的距离例 5.16 设 U=(1,0,0),V=(0,1,0),求 U,V两点间的距离 。
命令如下：
U=[1,0,0];V=[0,1,0];UV=U-V;
D=norm(UV)
5,向量的向量积例 5.17设 U=(2,-3,1),V=(3,0,4),求 U× V。
命令如下：
U=[2,-3,1];V=[3,0,4];W=eye(3);
A1=[W(1,:);U;V];A2=[W(2,:);U;V];A3=[W(3,:);U;V];
UV=[det(A1),det(A2),det(A3)]
UV=
-12 -5 9
6,向量的混合积例 5.18 设 U=(0,0,2),V=(3,0,5),W=(1,1,0),求以这三个向量构成的六面体的体积 。
命令如下：
U=[0,0,2];V=[3,0,5];W=[1,1,0];
A=[U;V;W];
det(A)
ans=
6
7,点到平面的距离例 5.19求原点到平面 X+Y+Z=1的距离 。
命令如下：
u=[0,0,0];v=[1,1,1]; %
A=B=C=1,u1=u2=u3=0,D=-1
r=abs(u*v'-1)/norm(v,2)
r =
0.5774
5.3 矩阵分解与线性方程组求解
5.3.1矩阵分解
1,实对称矩阵的 QDQ分解例 5.20设对称矩阵 A,对 A进行 QDQ分解 。
命令如下：
A=[2,1,4,6;1,2,1,5;4,1,3,4;6,5,4,2];
[Q,D]=eig(A)
Q*D*Q'
ans =
2.0000 1.0000 4.0000 6.0000
1.0000 2.0000 1.0000 5.0000
4.0000 1.0000 3.0000 4.0000
6.0000 5.0000 4.0000 2.0000
结果与 A相等,说明确实将 A分解为了 QDQ'的乘积 。
例 5.21求下列二次型的标准形式及变换矩阵 。
命令如下：
A=[1,2,1;2,1,1;1,1,3;];
[Q,D]=eig(A)
进一步作线性变换即得关于 u,v,w的标准二次型：
2,矩阵的 LU分解
MATLAB中,完成 LU分解的函数是：
(1)[L,U]=lu(A) 将方阵 A分解为交换下三角矩阵 L
和上三角矩阵 U,使 A=LU。
(2)[L,U,P]=lu(A) 将方阵 A分解为下三角矩阵 L和上三角矩阵 U,使 PA=LU。
例 5.22用 LU分解求方程组的根 。
3,矩阵的 QR分解对矩阵 A进行 QR分解的函数是 [Q,R]=qr(A),根据方阵 A,求一个正交矩阵 Q和一个上三角矩阵 R，
使 A=Q*R。 例如,对矩阵 A进行 QR分解的命令是：
A=[2,1,-2;1,2,1;2,5,3];
[Q,R]=qr(A)
5.3.2 线性方程组求解
1,线性方程组解的一般讨论解线性方程组的一般函数文件如下：
function [x,y]=line_solution(A,b)
[m,n]=size(A);y=[];
if norm(b)>0 %非齐次方程组
if rank(A)==rank([a,b]) %方程组相容
if rank(A)==m %有唯一解
x=A\b;
else %方程组有无穷多个解,基础解系
disp('原方程组有有无穷个解,其齐次方程组的基础解系为 y,特解为 x');
y=null(A,'r');
x=A\b;
end
else %方程组不相容,给出最小二乘法解
disp('方程组的最小二乘法解是,');
x=A\b;
end
else %齐次方程组
if rank(A)>=n %列满秩
x=zero(m,1) %0解
else %非 0解
disp('方程组有无穷个解,基础解系为 x');
x=null(A,'r');
2,应用举例例 5.23求线性方程组的解 。
在 MATLAB命令窗口,输入命令：
A=[2,2,-1,1;4,3,-1,2;8,5,-3,4;3,3,-2,2];b=[4,6,12,6]';
[x,y]=line_solution(A,b) %调用自定义函数例 5.24求下列线性方程组的解 。
在 MATLAB命令窗口,输入命令：
A=[2,7,3,1;3,5,2,2;9,4,1,7];b=[6,4,2]';
[x,y]=line_solution(A,b)
5.4 数据处理与多项式计算
5.4.1 数据统计与分析
1,求矩阵最大和最小元素
(1)求向量的最大最小元素
① y=max(X) 返回向量 X的最大元素存入 y。
② [y,I]=max(X) 返回向量 X的最大元素存入 y,最大元素的序号存入 I。
(2)求矩阵的最大和最小元素
① max(A) 返回一个行向量,向量的第 i个元素是 A矩阵的第 i列上的最大元素 。
② [Y,U]=max(A) 返回两个行向量,Y向量记录 A的每列的最大元素,U向量记录每列最大元素的行号 。
③ max(A,[],dim) dim取 1或 2。 dim取 1时,该函数和 max(A)
完全相同 。 dim取 2时,该函数返回一个列向量,其第 i
个元素是 A矩阵的第 i行上的最大元素 。
(3)两个向量或矩阵对应元素的比较
① U=max(A,B) A,B是两个同型的向量或矩阵 。 结果 U是与 A,B同型的向量或矩阵,U的每个元素等于 A,B对应元素的较大者 。
② U=max(A,n) n是一个标量 。 结果 U是与 A同型的向量或矩阵,U的每个元素等于 A对应元素和 n中的较大者 。
min函数的用法和 max完全相同 。
例 5.25 求矩阵 A的每行及每列的最大和最小元素,并求整个矩阵的最大和最小元 。
命令如下：
A=[13,-56,78;25,63,-235;78,25,563;1,0,-1];
max(A,[],2) %求每行最大元素
min(A,[],2) %求每行最小元素
max(A) %求每列最大元素
min(A) %求每列最小元素
max(max(A)) %求整个矩阵的最大元素
min(min(A)) %求整个矩阵的最小元素
2,求矩阵的平均值和中值求矩阵和向量元素的平均值的函数是 mean,求中值的函数是 median。 它们的调用方法和 max函数完全相同 。
3,矩阵元素求和与求积矩阵和向量求和与求积的基本函数是 sum和 prod,其使用方法和 max类似 。
例 5.26求矩阵 A的每行元素的乘积和全部元素的乘积 。
命令如下：
A=[1,2,3,4;5,6,7,8;9,10,11,12];
S=prod(A,2)
prod(S) %求 A的全部元素的乘积
4,矩阵元素累加和与累乘积
MATLAB中,使用 cumsum和 cumprod函数能方便地求得向量和矩阵元素的累加和与累乘积向量,函数的用法和 sum及 prod相同例 5.27求向量 X=(1！,2！,3！,…,10！ )。
命令如下：
X=cumprod(1:10)
5,标准方差
MATLAB中,提供了计算数据序列的标准方差的函数 std。
对于向量 X,std(X)返回一个标准方差 。 对于矩阵 A，
std(A)返回一个行向量,它的各个元素便是矩阵 A各列或各行的标准方差 。 std函数的一般调用格式为：
std(A,FLAG,dim)
其中 dim取 1或 2。 当 dim=1时,求各列元素的标准方差；
当 dim=2时,则求各行元素的标准方差 。 FLAG取 0或 1。
6,元素排序
MATLAB中对向量 X是排序函数是 sort(X),函数返回一个对 X中的元素按升序排列的新向量 。
sort函数也可以对矩阵 A的各列 (或行 )重新排序,
其调用格式为：
[Y,I]=sort(A,dim)
其中 dim指明对 A的列还是行进行排序,若 dim=1，
则按列排,若 dim=2,则按行排 。 Y是排序后的矩阵,而 I记录 Y中的元素在 A中位置 。
例 5.28对矩阵做各种排序 。
命令如下：
A=[1,-8,5;4,12,6;13,7,-13];
sort(A) %对 A的每列按升序排序
-sort(-A,2) %对 A的每行按降序排序
[X,I]=sort(A) %对 A按列排序,并将每个元素所在行号送矩阵 I
5.4.2 数值插值
1,一维数值插值
interp1函数调用格式为：
Y1=interp1(X,Y,X1,'method')
函数根据 X,Y的值,计算函数在 X1处的值 。 X,Y是两个等长的已知向量,分别描述采样点和样本值,X1是一个向量或标量,描述欲插值的点,Y1是一个与 X1等长的插值结果 。 method是插值方法,允许的取值有
'linear'(线性插值 ),'nearest'(最近插值 ),'spline'(三次样条插值 ),'cubic'（ 三次多项式插值 ）,缺省值是 'linear'。
例 5.29用不同的插值方法计算 sin(x)在 π/2点的值 。
这是一个一维插值问题 。 在 MATLAB命令窗口,输入命令：
X=0:0.2:pi;Y=sin(X); %给出 X,Y
interp1(X,Y,pi/2) %用缺省方法 (即线性插值方法 )计算
sin(π/2)
interp1(X,Y,pi/2,'nearest') %用最近方法计算 sin(π/2)
interp1(X,Y,pi/2,'linear') %用线性方法计算 sin(π/2)
interp1(X,Y,pi/2,'spline') %用三次样条方法计算 sin(π/2)
interp1(X,Y,pi/2,'cubic') %用三次多项式方法计算 sin(π/2)
MATLAB 中 有 一 个 专 门 的 三 次 样 条 插 值 函 数
Y1=spline(X,Y,X1),其功能及使用方法与函数
Y1=interp1(X,Y,X1,'spline')完全相同 。
例 5.30 已知检测参数 f随时间 t的采样结果,用数值插值法计算
t=2,7,12,17,22,17,32,37,42,47,52,57时 f的值 。
这是一个一维数值插值问题,命令如下：
T=0:5:65;
X=2:5:57;
F=[3.2015,2.2560,879.5,1835.9,2968.8,4136.2,5237.9,6152.7,...
6725.3,6848.3,6403.5,6824.7,7328.5,7857.6];
F1=interp1(T,F,X) %用线性方法插值
F1=interp1(T,F,X,'nearest') %用最近方法插值
F1=interp1(T,F,X,'spline') %用三次样条方法插值
F1=interp1(T,F,X,'cubic') %用三次多项式方法插值
2,二维数值插值
MATLAB中,提供了解决二维插值问题的函数 。
其调用格式为：
Z1=interp2(X,Y,Z,X1,Y1,'method')
其中 X,Y是两个向量,分别描述两个参数的采样点,Z是与参数采样点对应的采样变量的样本值,X1,Y1是两个向量或标量,描述欲插值的点 。 method的取值与一维插值函数相同 。
例 5.31设 Z=x2+y2,对 Z函数在 (0,1)× (0,2)区域内进行插值 。
命令如下：
x=0:0.1:10;y=0:0.2:20;
[X,Y]=meshgrid(x,y);
Z=X.^2+Y.^2;
interp2(x,y,Z,0.5,0.5) %对函数在 (0.5,0.5)点进行插值
interp2(x,y,Z,[0.5 0.6],0.4) %对函数在 (0.5,0.4)点和 (0.6,0.4)点进行插值
interp2(x,y,Z,[0.5 0.6],[0.4 0.5]) %对函数在 (0.5,0.4)点和 (0.6,0.5)点进行插值
interp2(x,y,Z,[0.5 0.6]',[0.4 0.5])%对函数在 (0.5,0.4),(0.6,0.4),(0.5,0.5)和 (0.6,0.5)点进行插值
3,三维数值插值对三维函数插值的函数是 interp3,其使用方法和 interp2相同 。 其调用格式为：
W1=interp3(X,Y,Z,W,X1,Y1,Z1,'method')
函数返回三维插值结果 。 其中 X,Y,Z是三个向量,分别描述三个参数的采样点,W是与参数采样点对应的采样变量的样本值,X1,Y1,Z1是三个向量或标量,
描述欲插值的点 。 method是插值方法,可选,其缺省值是 ‘ line'。 method的取值与一,二维插值函数相同 。
5.4.3 曲线拟合
MATLAB中,提供了解决使用最小二乘法进行曲线拟合的函数 。 调用格式为：
[P,S]=polyfit(X,Y,m)
函数根据采样点 X和采样点函数值 Y,产生一个 m
次多项式 P及其在采样点的误差向量 S。
其中 X,Y是两个等长的向量,P是一个长度为
m+1的向量 。
例 5.32 用一个 5次多项式在区间 [0,2π]内逼近函数 sin(x)。
命令如下：
X=linspace(0,2*pi,50);Y=sin(X);
[P,S]=polyfit(X,Y,5) %得到 5次多项式的系数和误差
plot(X,Y,'k*',X,polyval(P,X),'k-')
5.4.4 多项式计算
1,多项式的建立已知一个多项式的全部根 X求多项式系数的函数是 poly(X)，
该函数返回以 X为全部根的一个多项式 P,当 X是一个长度为 m的向量时,P是一个长度为 m+1的向量 。
2,多项式求根求多项式 p(x)的根的函数是 roots(P),这里,P是 p(x)的系数向量,该函数返回方程 p(x)=0的全部根 (含重根,复根 )。
3,多项式求值求多项式 p(x)在某点或某些点的函数值的函数是 polyval(P,x)。
若 x为一数值,则求多项式在该点的值；若 x为向量或矩阵,
则对向量或矩阵中的每个元素求其多项式的值 。
例 5.33 已知一个多项式,计算：
(1)计算 f(x)=0 的全部根 。
(2)由方程 f(x)=0的根构造一个多项式 g(x),并与 f(x)进行对比 。
(3)计算 f(5),f(7.8),f(9.6),f(12.3)的值 。
命令如下：
P=[3,0,4,-5,-7.2,5];
X=roots(P) %求方程 f(x)=0的根
G=poly(X) %求多项式 g(x)
X0=[5,7.8,9.6,12.3];
f=polyval(P,X0) %求多项式 f(x)在给定点的值多项式求值还有一个函数是 polyvalm,其调用格式与 polyval
相同,但含义不同 。 polyvalm函数要求 x为方阵,它以方阵为自变量求多项式的值 。
4,多项式的四则运算
(1)多项式的加减法
(2)多项式的乘法函数 conv(P1,P2)用于求多项式 P1和 P2的乘积 。
(3)多项式的除法函数 [Q,r]=deconv(P1,P2)用于对多项式 P1和 P2作除法运算 。 其中 Q返回多项式 P1除以 P2的商式,r返回 P1除以 P2的余式 。 这里,Q和 r仍是多项式系数向量 。
deconv是 conv的逆函数,即有 P1=conv(P2,Q)+r。
例 5.34设有两个多项式,计算：
(1)求 f(x)+g(x),f(x)-g(x)。
(2)求 f(x)·g(x),f(x)/g(x)。
在 MATLAB命令窗口,输入命令：
f=[3,-5,2,-7,5,6];g=[3,5,-3];g1=[0,0,0,g];
f+g1 %求 f(x)+g(x)
f-g1 %求 f(x)-g(x)
conv(f,g) %求 f(x)*g(x)
[Q,r]=deconv(f,g) %求 f(x)/g(x),商式送 Q,余式送 r。
5,多项式的导函数对多项式求导数的函数是：
p=polyder(P) 求多项式 P的导函数
p=polyder(P,Q) 求 P*Q的导函数
[p,q]=polyder(P,Q) 求 P/Q的导函数,导函数的分子存入 p，
分母存入 q。
例 5.35求有理分式的导数 。
命令如下：
P=[3,5,0,-8,1,-5];
Q=[10,5,0,0,6,0,0,7,-1,0,-100];
[p,q]=polyder(P,Q)
5.4.5 函数的最大值与最小值
MATLAB中用于求最小值的函数是：
fmin(f,a,b) 求单变量函数 f(x)在区间 (a,b)上的最小值点 。
fmins(F,X0) 求多变量函数 F(x)在估计值 X0附近的最小值点 。
MATLAB没有专门提供求函数最大值点的函数,但只要注意到 -f(x)在区间 (a,b)上的最小值点就是 f(x)在 (a,b)的最大值点,所以 fmin(-f,a,b)返回函数 f(x)在区间 (a,b)上的最大值 。
例 5.36 求函数 f(x)在区间 (-10,1)和 (1,10)上的最小值点 。
首先建立函数文件 fx.m：
function f=f(x)
f=x-1/x+5;
return
再在 MATLAB命令窗口,输入命令：
fmin('fx',-10,-1) %求函数在区间 (-10,-1)内的最小值点
fmin(f,1,10) %求函数在区间 (1,10)内的最小值点 。 注意函数名 f不用加 '
例 5.37 设有函数 f(x,y,z),求函数 f在 (0.5,0.5,0.5)附近的最小值 。
建立函数文件 fxyz.m：
function f=f(u)
x=u(1);y=u(2);z=u(3);
f=x+y.^2./x/4+z.^2./y+2./z;
return
在 MALAB命令窗口,输入命令：
U=fmins('fxyz',[0.5,0.5,0.5]) %求函数的最小值点
fxyz(U) %求函数的最小值
5.5 傅立叶分析
MATLAB中,提供了对向量 (或直接对矩阵的行或列 )进行离散傅立叶变换的函数,其调用格式是：
Y=fft(X,n,dim)
(1)当 X是一个向量时,返回对 X的离散傅立叶变换 。
(2)当 X是一个矩阵时,返回一个矩阵并送 Y,其列 (行 )是对 X的列 (行 )的离散傅立叶变换 。
例 5.38 求 X=(1,0,-3,5,2)的离散傅立叶逆变换 。
在 MATLAB命令窗口,输入命令：
X=[1,0,-3,5,2];
Y=fft(X) %对 X进行变换
3,离散傅立叶变换的逆变换
MATLAB中,对向量 (或直接对矩阵的行或列 )进行离散傅立叶逆变换的函数的调用方法是：
Y=ifft(X,n,dim)
函数对 X进行离散傅立叶逆变换 。 其中 X,n,dim的意义及用法和离散傅立叶变换函数 fft完全相同 。
例 5.39 对矩阵 A的列向量,行向量分别进行离散傅立叶变换,并对变换结果进行逆变换 。
命令如下：
A=[3,2,1,1;-5,1,0,1;3,2,1,5];
fftA=fft(A) %求 A的列向量的傅立叶变换
fftA2=fft(A,4,2) %求 A的行向量的傅立叶变换
ifft(fftA) %对矩阵 fftA的列向量进行傅立叶逆变换,
结果应等于 A
ifft(fftA2,4,2) %对矩阵 fftA2的行向量进行傅立叶逆变换,其结果应等于 A
5.6 数值微积分
5.6.1 数值微分
MATLAB中,没有直接提供求数值导数的函数,只有计算向前差分的函数 。
DX=diff(X) 计算向量 X的向前差分,DX(i)=X(i+1)-X(i)，
0<i<n。
DX=diff(X,n) 计算 X 的 n 阶 向 前 差 分,
diff(X,2)=diff(diff(X))。
DX=diff(A,n,dim) 计算矩阵 A的 n阶差分,dim=1时 (缺省状态 ),按列计算差分,dim=2,按行计算差分 。
例 5.40 求向量 sin(X)的 1~3阶差分 。 设 X由 [0,2π]间均匀分布的 10个点组成 。
命令如下：
X=linspace(0,2*pi,10);
Y=sin(X);
DY=diff(Y); %计算 Y的一阶差分
D2Y=diff(Y,2); %计算 Y的二阶差分,也可用命令
diff(DY)计算
D3Y=diff(Y,3); %计算 Y的三阶差分,也可用 diff(D2Y)
或 diff(DY,2)
例 5.41 用不同的方法求函数 f(x)的数值导数,并在同一个坐标系中做出 f'(x)的图象 。
程序如下：
f=inline('sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2');
g=inline('(3*x.^2+4*x-1)./sqrt(x.^3+2*x.^2-x+12)/2+1/6./(x+5).^(5/6)+5');
x=-3:0.01:3;
p=polyfit(x,f(x),5); %用 5次多项式 p拟合 f(x)
dp=polyder(p); %对拟合多项式 p求导数 dp
dpx=polyval(dp,x); %求 dp在假设点的函数值
dx=diff(f([x,3.01]))/0.01; %直接对 f(x)求数值导数
gx=g(x); %求函数 f的导函数 g在假设点的导数
plot(x,dpx,x,dx,'g.',x,gx,'r-'); %作图
5.6.2数值积分
(1)被积函数是一个解析式函数 quad(f,a,b,tol,trace)用于求被积函数 f(x)在 [a,b]上的定积分,tol是计算精度,缺省值是 0.001。 trace非 0时,画出积分图形 。 注意,调用 quad函数时,先要建立一个描述被积函数 f(x)的函数文件或语句函数 。 当被积函数
f含有一个以上的变量时,quad函数的调用格式为：
quad(f,a,b,tol,trace,g1,g2)
其中 f,a,b,tol,trace等参数的含义同前 。
数值积分函数还有一种形式 quad8,其用法与 quad完全相同 。
例 5.42 用两种不同的方法求积分 。
先建立一个函数文件 ex.m：
function ex=ex(x)
ex=exp(-x.^2); %注意应用点运算
return
然后,在 MATLAB命令窗口,输入命令：
quad('ex',0,1,1e-6) %注意函数名应加字符引号
quad8('ex',0,1,1e-6) %用另一函数求积分例 5.43用 trapz函数计算积分 。
在 MATLAB命令窗口,输入命令：
X=0:0.01:1;Y=exp(-X.^2);
trapz(X,Y)
(2)被积函数由一个表格定义
MATLAB中,对由表格形式定义的函数关系的求定积分问题用 trapz(X,Y)函数 。 其中向量 X,Y
定义函数关系 Y=f(X)。
(3)二重积分例 5.44计算二重积分 。
建立一个函数文件 fixy.m：
functionf=f(x,y)
f=exp(-x.^2-y.^2);
return
建立一个命令文件 ftxy1.m：
for i=1:20
int2(i)=quad('fixy',0,1,[],[],x(i)); %在二维函数 fixy中以 x=x(i)代入并对 y积分 。
end
在 MATLAB命令窗口,输入命令：
x=linspace(0,1,20);
ftxy1
trapz(x,int2)
实际上,MATLAB提供了计算二重积分的函数：
dblquad(f,a,b,c,d,tol,trace)
该函数求 f(x,y)在 [a,b]× [c,d]区域上的二重积分 。
参数 tol,trace的用法与函数 quad完全相同 。
如果直接使用这里介绍的二重积分函数 dblquad来求解本例就非常简单,命令如下：
g=inline('exp(-x.^2-y.^2)');
dblquad(g,0,1,0,1) %直接调用二重积分函数求解
5.7 常微分方程的数值求解基于龙格－库塔法,MATLAB提供了求常微分方程数值解的函数,一般调用格式为：
[X,Y]=ode23(f,[x0,xn],y0)
[X,Y]=ode45(f,[x0,xn],y0)
其中 X,Y是两个向量,X对应自变量 x在求解区间 [x1，
xn]的一组采样点,其采样密度是自适应的,无需指定；
Y是与 X对应的一组解,f是一个函数,[x0,xn]代表自变量的求解区间,y0=y(x0),由方程的初值给定 。 函数在求解区间 [x0,xn]内,自动设立采样点向量 X,并求出解函数 y在采样点 X处的样本值 。
例 5.45 求微分方程初值问题在 [1,3]区间内的数值解,并将结果与解析解进行比较 。
先建立一个该函数的 m文件 fxy1.m：
function f=f(x,y)
f=-2.*y./x+4*x %注意使用点运算符
return
再输入命令：
[X,Y]=ode45('fxy1',[1,3],2);
X' %显示自变量的一组采样点
Y' %显示求解函数与采样点对应的一组数值解
(X.^2+1./X.^2)' %显示求解函数与采样点对应的一组解析解例 5.46 求解初值问题在区间 [0,2]中的解 。
建立一个函数文件 fxy2.m：
function f=f(x,y)
f(2)=-x.*y(2)+x.^2-5;
f(1)=y(2);
f=f';
return
在 MATLAB命令窗口,输入命令：
[X,Y]=ode45('fxy2',[0,2],[5,6]);
[X,Y]
5.8 非线性方程的数值求解
1,单变量非线性方程求解
MATLAB中,提供了求解单变量方程的函数 fzero(f,x0,tol)，
该函数采用迭代法计算函数 f(x)的一个零点,迭代初值为 x0,当两次迭代结果小于 tol时停止迭代过程 。 tol的缺省值是 eps。
注意,在调用函数 fzero 之前,要使用 m文件建立自己要计算的函数 f(x),只有定义了函数 f(x)的 m文件后,才能在 fzero函数的参数中使用自定义函数名 。
例 5.47 求 f(x)=x-+5 在 x0=-5和 x0=1作为迭代初值时的零点 。
先编制一个函数文件 fz.m：
function f=f(x)
f=x-1/x+5;
然后,在 MATLAB命令窗口,输入命令：
fzero('fz',-5) %以 -5作为迭代初值
Zero found in the interval,[-4.8,-5.2].
fzero('fz',1)
2,非线性方程组求解函数 fsolve调用格式为：
X=fsolve(F,X0)
例 5.48 求方程组在 (1,1,1)附近的解并对结果进行验证 。
首先建立方程的函数文件 fxyz1.m：
function F=F(X)
x=X(1);y=X(2);z=X(3);
F(1)=sin(x)+y+z^2*exp(x);
F(2)=x+y*z;
F(3)=x*y*z;
在 MATLAB命令窗口,输入命令：
X=fsolve('fxyz1',[1,1,1]) %求解 X的三个分量 x,y,z
Y=fxyz1(X) %检验所求结果 X是否满足原方程组
norm(Y) %求 Y向量的模例 5.49 求圆和直线的两个交点 。
建立方程组函数文件 fxyz2.m：
function F=F(X)
x=X(1);y=X(2);z=X(3);
F(1)=x^2+y^2+z^2-9;
F(2)=3*x+5*y+6*z;
F(3)=x-3*y-6*z-1;
在 MATLAB命令窗口,输入命令：
X1=fsolve('fxyz2',[-1,1,-1]) %求直线与球面的第一个交点
X2=fsolve('fxyz2',[1,-1,1]) %求直线与球面的第二个交点
5.9 稀疏矩阵
5.9.1 矩阵存储方式
1,矩阵的完全存储模式
2,稀疏矩阵的存储方式
5.9.2 稀疏存储方式的产生与转化
1.将一个完全存储方式的转化为稀疏存储方式函数 B=sparse(A)将矩阵 A转化为稀疏存储方式的矩阵 B。
sparse函数还有其他一些格式：
sparse(m,n) 生成一个 m× n的所有元素都是 0的稀疏矩阵 。
sparse(u,v,S) u,v,S是三个等长的向量 。
此外,还有一些和稀疏矩阵操作有关的函数 。 例如
[U,V,S]=find(A) 返回矩阵 A中非 0元素的下标和元素 。 这里产生的 U,V、
S可作为 sparse(u,v,s)的参数 。
full(A) 返回和稀疏存储矩阵 A对应的完全存储方式矩阵 。
2,产生一个稀疏矩阵把要建立的稀疏矩阵的非 0元素及其所在行和列的位置表示出来后由 MATLAB自己产生其稀疏存储方式,这需要使用 spconvert函数 。 调用格式为：
B=spconvert(A)
其中 A为一个 m× 3或 m× 4的矩阵,其每行表示一个非 0元素,m是非 0元素的个数 。
3,单位稀疏矩阵的产生单位矩阵只有对角线元素为 1,其他元素都为 0,是一种具有稀疏特征的矩阵 。 我们知道,函数 eye产生一个完全存储方式的单位矩阵 。 MATLAB还有一个产生稀疏存储方式的单位矩阵的函数,这就是 speye。 函数
speye(m,n)返回一个 m× n的稀疏存储单位矩阵 。
