matlab杂例.docx - 冰豆网

资源描述

matlab杂例.docx

《matlab杂例.docx》由会员分享，可在线阅读，更多相关《matlab杂例.docx（23页珍藏版）》请在冰豆网上搜索。

matlab杂例.docx

matlab杂例

公元2011一剑飞雪倾情奉献，悉心整理……………………

matlab中contour的用法

悬赏分：

0|解决时间：

2008-11-1915:

58|提问者：

xiaoxin2305

matlab中contour的用法和参数的说明，中文的，越详细越好。

我要画它等值线图。

给你一个实例吧。

x=1:

14;

y=1:

14;

[xx,yy]=meshgrid（y,x）;

z=[1234567891011121314;

1234567891011121314;

1234567891011121314;];

contour（xx,yy,z,15）;holdon

contour（yy,xx,z,15）结果如图。

。

用matlab求圆周率

如何用Matlab画矢量图

[xx,yy]=meshgrid（-2:

0.2:

2）;

zx=xx-2*yy;

zy=6*xx+yy;

quiver（xx,yy,zx,zy）结果如下

gradient的用法－－matlab的一个函数

怎么样用gradient绘制z=1/（x^2-2*x+4）+1/（y^3－2*y+4）的梯度场呢？

[x,y]=meshgrid（[-5:

0.5:

5]）

z=1./（x.^2-2*x+4）+1./（y.^3-2*y+4）

g=gradient（z）

subplot（1,2,1）,surf（x,y,z）

subplot（1,2,2）,surf（x,y,z,g）

MATLAB中subs（subs（））两个连起来怎么使用？

symsab

y=a^2+sin（b）;

subs（subs（y,a,2）,b,pi）%就是将a=2,b=pi代入y的表达式，分两次代。

matlab中加点什么用？

悬赏分：

0|解决时间：

2008-12-2413:

05|提问者：

139bbc

比如y=1./（1+x.^4）

为什么1后面也要加点X后面要加点？

在什么时候就不用加点了

点乘表示和矩阵相乘的时候用的。

。

6.4矩阵的数组运算：

矩阵乘：

A*B。

A的列数与B的行数要相等。

矩阵的数组乘：

A.*B。

表示为A、B矩阵的对应元素一一相乘，即Aij*Bij。

A与B的维度要相同。

矩阵的数组除：

A./B或者A.\B

举例：

x=-10:

y=sin（x）%correct

y=sin（x^2）%incorrect-Matrixxmustbesquare

y=sin（x.^2）%correct

clear

x=0:

250;|

y=x';

X=ones（size（y））*x;|

Y=y*ones（size（x））;

Z=40*X-X.^2+30*Y-Y.^2;

surf（X,Y,Z）

你忘记在X，Y后面加.了

x=[1,2,3;4,5,6;7,8,9]

surf（x）

inline是用来定义内联函数的

比如说：

y=inline（'sin（x）','x'）%第一个参数是表达式，第二个参数是函数变量

y（0）%计算sin（0）的值

y（pi）%计算sin（pi）的值

q=quad（y,0,1）;%计算sin（x）在0到1上的积分

试求双重积分

[说明]该函数形式比较简单，可首先定义函数f，接着直接调用dblquad（）求得该积分的数值解。

>>formatlong

>>fun=inline（'exp（-x.^2）.*cos（x+y.^2）','x','y'）;

>>Y_num=dblquad（fun,1,2,0,10）

Y_num=

.0519********

matlab控制运算精度用的是digits和vpa这两个函数

digits用于规定运算精度，比如：

digits（20）;

这个语句就规定了运算精度是20位有效数字。

但并不是规定了就可以使用，因为实际编程中，我们可能有些运算需要控制精度，而有些不需要控制。

vpa就用于解决这个问题，凡是用需要控制精度的，我们都对运算表达式使用vpa函数。

例如：

digits（5）;

a=vpa（sqrt

（2））;

这样a的值就是1.4142，而不是准确的1.4142135623730950488016887242097

又如：

digits（5）;

a=vpa（sqrt

（2））;

b=sqrt

（2）;

这样a的值是1.4142，b没有用vpa函数，所以b是1.4142135623730950488016887242097......

vpa函数对其中每一个运算都控制精度，并非只控制结果。

digits（11）;

a=vpa（2/3+4/7+5/9）;

b=2/3+4/7+5/9;

a的结果为1.7936507936，b的结果为1.793650793650794......也就是说，计算a的值的时候，先对2/3，4/7，5/9这三个运算都控制了精度，又对三个数相加的运算控制了精度。

而b的值是真实值，对它取11位有效数字的话，结果为1.7936507937，与a不同，就是说vpa并不是先把表达式的值用matlab本身的精度求出来，再取有效数字，而是每运算一次，都控制精度。

我举的例子不太好，因为加法不太会增加数字位数。

希望你能理解我的意思....

方阵可逆时，用inv就可以了，不可逆时，用pinv得到的是广义逆

42-1

3-12

1130

则A（:

1）即A的第一列

ans=

rref求增广矩阵的行最简形式，可得到最简同解方程组，

A=[6,-2,2,1;1,-1,0,1;2,0,1,3]

6-221

1-101

2013

>>B=[3,1,2]'

>>C=[AB]

6-2213

1-1011

20132

rref（C）

ans=

1.000000.500000.3571

01.00000.50000-0.2143

0001.00000.4286

求解方程组x1+2x2+3x3+4x4=1的解

2x1+2x2+x3+x4=3

2x1+4x2+6x3+8x4=2

4x1+4x2+2x3+2x4=6

>>A=[1,2,3,4;2,2,1,1;2,4,6,8;4,4,2,2]

1234

2211

2468

4422

>>B=[1;2;3;4]

>>C=[AB]

12341

22112

24683

44224

>>[rank（A）,rank（C）]

ans=

因为rank（A）不等于rank（C），故无解，只能利用函数pinv（）求取Moore-Penrose广义逆，从而得到原方程的最小二乘解。

x=pinv（A）*B

0.5466

0.4550

0.0443

-0.0473

>>norm（A*x-B）验证所得的解的精确性

ans=

0.4472

一个稀疏矩阵中有许多元素等于零,这便于矩阵的计算和保存.如果MATLAB把一个矩

　　阵当作稀疏矩阵,那么只需在m×3的矩阵中存储m个非零项.第1列是行下标,第2列是列下

　　标,第3列是非零元素值,不必保存零元素.如果存储一个浮点数要8个字

高斯

节,存储每个下标

　　要4个字节,那么整个矩阵在内存中存储需要16×m个字节.

　　A=eye（1000）;

　　得到一个1000×1000的单位矩阵,存储它需要8Mb空间.如果使用命令:

　　B=speye（1000）;

　　用一个1000×3的矩阵来代表,每行包含有一个行下标,列下标和元素本身.现在只需

　　16Kb的空间就可以存储1000×1000的单位矩阵,它只需要满单位矩阵的0.2%存储空间.对于

　　许多的广义矩阵也可这样来作.

在MATLAB中,用命令sparse来创建一个稀疏矩阵.

　　命令集87创建稀疏矩阵

　　sparse（A）由非零元素和下标建立稀疏矩阵A.如果A已是一个稀疏

　　矩阵,则返回A本身.

　　sparse（m,n）生成一个m×n的所有元素都是0的稀疏矩阵.

　　sparse（u,v,a）生成一个由长度相同的向量u,v和a定义的稀疏矩阵.其

　　中u和v是整数向量,a是一个实数或者复数向量.（ui,vi）

　　对应值ai,如果a中有零元素,则将这个元素排除在外.

　　稀疏矩阵的大小为max（u）×max（v）.

　　sparse（u,v,a,m,n）生成一个m×n的稀疏矩阵,（ui,vi）对应值ai.向量u,v和a

　　必须长度相同.

　　sparse（u,v,a,m,n,生成一个m×n的含有nzmax个非零元素的稀疏矩阵.（ui,

　　nzmax）vi）对应值ai.nzmax的值必须大于或者等于向量u和v的长度.

　　find（x）返回向量x中非零元素的下标.如果x=X是一

稀疏矩阵

个矩阵,那

　　么X的向量就作为一个长向量来考虑.

　　[u,v]=find（A）返回矩阵A中非零元素的下标.

　　[u,v,s]=find（A）返回矩阵A中非零元素的下标.用向量s中元素的值及u和v中

　　相应的下标,实际上就是向量u,v和s作为命令sparse的参数.

　　spconvert（D）将一个有三列的矩阵转换成一个稀疏矩阵.D中的第1列作

　　为行的下标,第2列作为列的下标,最后一列作为元素值.

　　而且可以使用命令full将稀疏矩阵转换成一个满矩阵.

　　命令集88转换成满矩阵

　　full（S）将稀疏矩阵S转换成一个满矩阵.

　　a）创建一个5×5的单位矩阵:

　　A=eye（5）

　　将矩阵A转换成稀疏矩阵B:

非线性方程的数值解法有：

二分法，迭代法，

非线性方程的数符号法有：

solve（）

求ax2+bx+c=0和e-x=sin

的解

x=solve（'a*x^2+b*x+c'）

[1/2/a*（-b+（b^2-4*a*c）^（1/2））]

[1/2/a*（-b-（b^2-4*a*c）^（1/2））]

x_start=solve（'exp（-x）=sin（pi*x/2）','x'）

x_start=

1.3935273560689078034994453392807/pi

fzero（）函数可以求出一维变量的零点

求非线性方程x3-x2-1=0的根。

调用格式x=fzero（fun,x0）

fzero（fun,x0,option）

[x,fval]=fzero（…）

x是方程的零点，fval是计算终止时的函数值，fun是方程的函数，x0为初始点；option是选择项，它包括Display和Tolx

Display为显示迭代情况，它有如下参数：

off表示不显示；iter表示显示；final表示只显示最终的结果。

Tolx表示x的终止精度

options=optimset（'Display','iter'）;

x=fzero（'x^3-x^2-1',1.5,options）

Func-countxf（x）Procedure

11.50.125initial

21.45757-0.0278754search

Lookingforazerointheinterval[1.4576,1.5]

31.46531-0.000918587interpolation

41.465573.90703e-007interpolation

51.46557-9.8848e-011interpolation

61.46557-4.44089e-016interpolation

71.465572.22045e-015interpolation

Zerofoundintheinterval:

[1.4576,1.5].

1.4656

fsolve（）主要是求非线性方程组的解。

调用格式x=fsolve（fun,x0）

fsolve（fun,x0,option）

[x,fval]=fsolve（…）

x是方程的零点，fval是计算终止时的函数值，fun是方程的函数，x0为初始点；option是选择项，它包括Display和Tolx、jacobian、maxfunevals、maxiter;

Display为显示迭代情况，它有如下参数：

off表示不显示；iter表示显示；final表示只显示最终的结果。

Tolx表示x的终止精度maxfunevals是调用函数的最大次数，maxiter是指最大迭代次数

例子x

-5=0

（x

+1）x

-（3x

+1）=0取初始点x

=（1,1）

编写m文件

functiondx=ch52ff1（x）

dx=[x

（1）^2+x

（2）^2-5;（x

（1）+1）*x

（2）-（3*x

（1）+1）];

>>x0=[1;1]；

>>options=optimset（'jacobian','off','display','iter'）;

>>x=fsolve（@ch52ff1,x0,options）

NormofFirst-order

IterationFunc-countf（x）stepoptimalityCG-iterations

14131100

272.738281.274758.021

3100.01560420.3344890.5651

4135.87406e-0070.0276840.003061

5162.51492e-0150.0002160212.23e-0071

Optimizationterminatedsuccessfully:

RelativefunctionvaluechangingbylessthanOPTIONS.TolFun

1.0000

2.0000

求矩阵x使其满足方程x*x*x=

并设初始解向量为x=[1,1;1,1].

编写m文件

functionF=ch52ff3（x）

F=x*x*x-[1,2;3,4];

x0=ones（2,2）

x0=

>>opitions=optimset（'display','off'）;

>>[x,fval]=fsolve（@ch52ff3,x0,options）

NormofFirst-order

IterationFunc-countf（x）stepoptimalityCG-iterations

16141240

2110.889160.6123724.162

3160.1396680.2034870.4132

4210.1328920.01910470.01222

5260.13270.0116260.04542

6310.1324920.005730250.01982

7360.1322670.01056430.05212

8410.1320040.007126310.02692

9460.1316920.01102060.05812

10510.13130.009806410.03572

11560.1308150.0127040.06412

12610.1301610.0141530.04842

13660.1293280.01583040.07352

14710.1280880.02145330.06832

15760.1264760.02132180.08882

16810.1236890.03514270.1032

17860.1199930.03169210.1152

18910.1118930.06515880.172

19960.101140.05384540.1632

201010.07078370.1365650.2832

211060.03973250.09172980.1972

221110.006965990.1399730.1522

231160.0002326340.03764330.02092

241211.04115e-0050.007835430.01582

251263.78107e-0070.001320610.0007452

261313.7973e-0080.0002858570.0002712

Optimizationterminatedsuccessfully:

RelativefunctionvaluechangingbylessthanOPTIONS.TolFun

-0.12910.8602

1.29031.1612

fval=

1.0e-003*

0.1541-0.1163

0.0109-0.0243

Ezplot符号绘图函数，在图形窗口绘出函数的图形。

对于符号函数，MATLAB提供了一个非常简单的作图指令：

ezplot（）函数。

通过这个命令，可以在图型窗口绘制出符号函数图形。

Ezplot（）函数有以下调用格式。

ezplot（f）对于符号函数f=f（x），按照x默认的范围：

－2*pi

的图形。

如：

f=’cox（x）’ezplot（f）;

ezplot（f,[a,b]）在图形窗口中绘制符号函数f=f（x）的图形，x的范围由[a,b]确定，即a

ezplot（f）对于符号函数f=f（x,y），ezplot（f）在图形窗口中绘制符号方程f（x,y）=0的图形。

x和y的取值范围为-2*pi

ezplot（f,[xmin,xmax,ymin,ymax]）对于符号函数f=f（x,y）,在没有指定y的范围时，y的范围和x的一致。

用axisequal；可以维持两坐标的一致。

ezplot（x,y）对于符号函数x=x（t），y=y（t），ezplot（x,y）在图形窗口中绘制符号方程x=x（t）

y=y（t）的图形。

t的取值范围为：

ezplot（x,y,[tmin,tmax]）设定t的取值范围：

tmin

fnplt绘制样条函数图形

matlab的样条插值工具箱中提供csapi（）函数来定义一个三次样条函数类，调用格式为S=csapi（x,y）,x和y分别为n阶向量

lsqcurvefit的使用

[xresnorm]=lsqcurvefit（fun,x0,xdata,ydata,...）

fun 是我们需要拟合的函数，这是重点

x0 是我们对函数中各参数的猜想值，这也是重点

xdata 则是横轴坐标的值

ydata 是纵轴的值

观察数据，发现该数据一个指数衰减的函数和一个三角正弦函数的叠加。

然后用psd观察频谱，它在0.4004处有最大峰，然后在0点和其他各处也有。

因为我们拟合要得到的结果是初始相位，只对0.404处的信号感兴趣。

以下是我们需要拟合出来的一个函数。

例

试用最小二乘法求一个多项式，使其与已知数据相拟合。

x=[1,3:

8]'

y=[10542112]'

>>plot（x,y）

>>axis（[0.590.511]）

可设拟合曲线为二次多项式y=

functiony=fitline_1（a,x）

y=a

（1）+a

（2）*x+a（3）*x.^2

>>a=lsqcurvefit（'fitline_1',[1;1;1],x,y）；

>>a'

ans=

13.2916-3.51600.2575

>>y1=fitline_1（a,x）;

10.0331

5.0613

3.3480

2.1498

1.4666

1.2984

1.6454

>>plot（x,y,x,y1,'o'）

>>holdon

>>plot（x,y1）拟合结果

求解下面的微分方程的通解2y''（t）+4y'（t）+y+6t=

>>y=dsolve（['2*D2y+y+4*Dy+6*t=exp（-6*t）*cos（4*t）']）;

>>simple（y）

simplify:

17/13378*exp（（-6+4*i）*t）+17/13378*exp（（-6-4*i）*t）+24-6*t+40/6689*i*exp（（-6+4*i）*t）-40/6

展开阅读全文