《之四程序设计.ppt》由会员分享，可在线阅读，更多相关《之四程序设计.ppt（42页珍藏版）》请在三一文库上搜索。
1、MATLAB 软件使用简介,之四,程序设计,MATLAB提供了一个完善的程序设计语言环境，使我们能方便地编制复杂的程序，完成各种运算。本节先介绍关系、逻辑运算和条件循环语句等编程的重要手段，然后再介绍MATLAB所特有的M文件。M文件的语法与C语言十分相似，对广大C语言爱好这来说，编写相当容易,一、关系和逻辑运算,1关系运算,MATLAB中的关系运算符为：,=等于； 大于； =大于等于； =不等于。,关系运算是比较二个数之间的关系，若指出的关系成立，则得到逻辑值（结果）为1（表示真），否则为0（表示假）,如： A=1:5,B=5:-1:1,A = 1 2 3 4 5,B = 5 4 3 2 1
2、,C=A=4,D=A=B,C = 0 0 0 1 1,D = 0 0 1 0 0,注：关系运算与数值运算可以结合起来,x=(-3:3)/3;x=x+(x=0)*eps; y=sin(x)./x,y = 0.8415 0.9276 0.9816 1.0000 0.9816 0.9276 0.8415,2逻辑运算,MATLAB中的逻辑运算符为：,&与运算；|或运算；非运算。,前二个运算有二个运算对象，形如A&B或A|B第三个逻辑运算只有一个运算对象。运算规则为：,逻辑运算将任何非零元素视为1（真值），逻辑运算也可以作用于矩阵或数组。,a=1:9,b=9-a,a = 1 2 3 4 5 6 7 8
3、9,b = 8 7 6 5 4 3 2 1 0,c=(a4),d=(a=3)&(b6),c = 1 1 1 1 0 0 0 0 0,d = 0 0 0 1 1 1 1 1 1,MATLAB提供了一些关系和逻辑函数，常见的有all和any，用法：,Y=all(x)：若x为向量，当所有元素非零时Y=1，否则Y=0；若x 为矩阵，all作用与列元素，Y为列向量。,Y=any(x)：若x为向量，当有一个元素非零时Y=1，否则Y=0；若x 为矩阵，any作用与列元素，Y为列向量。,a=1,0,-5,0;-3,0,8,2,b=all(a) , c=any(a),d=all(b),e=any(c),a =
4、1 0 -5 0 -3 0 8 2,b = 1 0 1 0,c = 1 0 1 1,d = 0,e = 1,二、条件和循环语句,条件和循环语句属于控制语句，MATLAB控制语句有个：if,while,for,witch,它们都用end结束。, if语句,if end,if语句的最简单用法是：,如果关系表达式的值为，则执行语句；否则，执行end的后续命令。,if语句的另一种用法是：,if else end,多个选择时还可以用下列结构：,if elseif elseif else end,如果关系表达式的值为，则执行语句；否则(关系表达式值为)，执行语句，然后再执行end的后续命令。,如果关系表达
5、式j(j=1,2,n)的值为，则执行语句j,然后执行end的后续命令；否则跳过，不予执行，直至到执行语句n+1, 然后执行end的后续命令。,例如，可以用下列程序得到右图所示的分段函数。,x=-2:0.2:2;n=length(x); for i=1:n if x(i)=0 y(i)=0; elseif x(i)=1 y(i)=x(i); else y(i)=1; end end plot(x,y),for语句,循环语句for的一般形式为：,for： end,步长为时可以省略。对每个参数，语句都重复执行。,注：,当作多次循环时，for 语句可以嵌套使用，如下列程序可以生成Hilbert矩阵。,
6、for I=1:3 for j=1:4 a(I,j)=1/(I+j-1); end end format rat,a = 1 1/2 1/3 1/4 1/2 1/3 1/4 1/5 1/3 1/4 1/5 1/6,for语句的循环参数可以是任意数组或矩阵，循环参数依次取数组元素的值，或按矩阵的列依次取值。,while语句,循环语句for的循环判别条件通常是对循环次数的判别 。如果不知道循环参数，可用while循环来完成。While语句的判别控制可以是一个逻辑判别语句，因此它的适用范围会更广些。其表达形式为：,while end,关系表达式的值为（真时）时，语句被循环执行，直到关系表达式值为（假
7、时）时终止。如：,n=0;Eps=1; while (1+Eps)1 Eps=Eps/2; n=n+1; end Eps=2*Eps; n,format short e,Eps,n = 53 Eps = 2.2204e-016,这个例子给出了计算MATLAB中特殊常数的eps过程（我们用大写Eps以便与eps相区别）。Eps不断地被2除直到（1+ Eps）1为假时终止。这里要注意的是，MATLAB用16位数来表示数据，因此当Eps接近时它会认为（1+ Eps）1不成立。,MATLAB还提供了跳出循环系统的break语句，遇到此语句立即执行循环end后续语句。如上面计算eps的程序可以改为：,E
8、ps=1; for n=1:100 Eps=Eps/2; if (1+Eps)=1 Eps=2*Eps; break end end n,format short e,Eps,矩阵指数的计算,Expm(A)=I+A+A2/2!+A3/3!+,e=zeros(size(a); f=eye(size(a); k=1; while norm(e+f-e)0 e=e+f; f=a*f/k; k=k+1; end e, k,a=2,3;3,4 e=zeros(size(a); f=eye(size(a); k=1; while norm(e+f-e)0 e=e+f; f=a*f/k; k=k+1; en
9、d e, k,e = 162.7871 224.6754 224.6754 312.5707 k = 38,e=expm(a) e = 162.7871 224.6754 224.6754 312.5707,4 Switch语句,Switch语句根据表达式的值来执行相应语句。用法如下：,switch case value1 case value2 otherwise end,A=1,2,3,4;2,4,6,8; n,m=size(A); switch(m-n); case 0 disp(A是方阵); otherwise disp(A不是方阵); end,A=1,2,3,4;2,4,6,8 1,
10、3,5,7;1,1,1,1; n,m=size(A); switch(m-n); case 0 disp(A是方阵); otherwise disp(A不是方阵); end,A不是方阵,A是方阵,用循环语句生成一个n维向量，要求它的元素满足Fibonacci规则，即向量第k+2个元素满足 ak+2=ak+ak+1，其中k=1,2, ,且a1=1,a2=1,while(1) n=input(输入维数n(n2); if (n2) break; end end,A=1,1; for i=3:n A(i)=A(i-2)+A(i-1); end disp(A=); disp(A);,A=1 1 2 3
11、5 8 13 21 34 55 89 144 233 377 610,三、M文件,为了处理复杂的问题和大量数据，MATLAB提供了M文件以解决此问题。M文件有两种类型：文本M文件和函数M文件。使用文本文件将自动执行一系列命令（相当于DOS中的批处理文件），函数文件提供了MATLAB的外部函数，M文件可用各类文字编辑器或文字处理器建立。,M文件包含了一系列的MATLAB命令，存储在磁盘的某个目录下，M文件可调用其它M文件（包括自己本身）。,M文件也可在WINOWS中File菜单中选择New,再选择Mfile，这时MATLAB将打开一个文本编辑窗口，在此窗口下输入命令和数据。,1文本M文件,一个比
12、较复杂的程序常常反复调用，这时不妨建立一个文本文件把它储存起来，可以随时调用进行计算。M文件必须以m为扩展名。,其一般形式为：,M,如hilblm，ppm等。,值得注意的是，文本M文件中的变量都是全局变量，在执行过程中，文本M文件中的命令可以使用目前工作区中的变量，它所产生的变量也将成为工作区的一部分。,比如我们把生成Hilbert矩阵的程序写成如下的文本M文件：,a=ones(m,n); for I=1:m for j=1:n a(I,j)=1/(I+j-1); end end a=rats(a),命名为Hilbertm，那么当需要一个23 Hilbert矩阵时，可以在MATLAB工作区中进
13、行：,m=2,n=3,Hilbert,a,m = 2 n = 3 a = 1 1/2 1/3 1/2 1/3 1/4,2函数M文件,函数M文件是另一类M文件，我们可根据需要建立自己的函数文件，它能够像库函数一样方便地调用，从而极大地扩展MATLAB的功能。如果对一类特殊问题，建立起许多函数M文件，就能最终形成独立的工具箱。,函数M文件的第一行有特殊要求，其形式必须为：,function =（）,其他各行为从自变量计算因变量的语句，并最终将结果赋值于因变量，而这个M文件的文件名必须是m。下面给出一个函数文件的简单例子。,如果我们经常要调用这样的随机矩阵，其每个元素等概率地取从0到9的整数值，就不
14、妨建立如下的函数M文件：,function a=randint(m,n) % FANDINT Randomly generated integral matrix. % randint(m,n) retuns an m-by-n such matrix with entries between 0 and 9 . a=floor(10*rand(m,n);,x=randint(2,3),x = 9 6 8 2 4 7,函数M文件有许多因变量时，要用 将它们括起来。如：建立一个数的平均值和均方差函数，,function mean,stdev=stat(x) % STAT mean and sta
15、ndard deviation % For a vector x , stat(x) returns the mean and standard deviation of x % For a matrix x , stat(x) retuns two row vectors containing , respectively. m,n=size(x); if m=1 m=n; % handle case of a row vector. end mean =sum(x)/m; stdev=sqrt(sum (x.2)/m-mean.2);,x=2,4,-7,0,5,-1 ; xm,xd=sta
16、t(x),xm = 0.5000 xd = 3.9476,function A=Fi(n) A=1,1; for j=3:n A(j)=A(j-2)+A(j-1); end,建立生成Fiobnacci数列函数,A=Fi(10),A = 1 1 2 3 5 8 13 21 34 55,建立生成第n个Fibonacci数,function A=Fib(n) A=1,1; for j=3:n A(j)=A(j-2)+A(j-1); end A=A(n),f=Fib(20),f = 6765,应用举例：,逢山开路（1994年数学建模竞赛试题A）,要在一山区修建公路，首先测得一些地点的高程，数据见表1（
17、平面区域 ,表中数据为坐标点的高程，单位：米）。数据显示：在y=3200处有一东西走向的山峰；从坐标（2400，2400）到（4800，0）有一西北东南走向的山谷；在（2000，2800）附近有一山口湖，其最高水位略高于1350米，雨季在山谷形成一溪流，雨量大是溪流水面宽度W与（溪流最深处的）x坐标的关系可近似表示为：,公路从山脚（0，800）处开始，经居民点（4000，2000）至矿区（4000，4000）。已知路段工程成本及对路段坡度（上升高程与水平距离之比）的限制如表2。,1）试给出一种线路设计方案，包括原理、方法及比较精确的线路位置（含桥梁、隧道），并估算该方案的总成本。 2）如果居民
18、点改为的居民区，公路只须经过居民区即可，那么你的方案有什么变化。,表一 4800 1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150 4400 1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210 4000 1380 1410 1430 1450 1470 1320 1280 1200 1080 940 780 620 460 370 350 3600 1420 1430 1450 1480 1500 1550 1510 1430 13
19、00 1200 980 850 750 550 500 3200 1430 1450 1460 1500 1550 1600 1550 1600 1600 1600 1550 1500 1500 1550 1550 2800 950 1190 1370 1500 1200 1100 1550 1600 1550 1380 1070 900 1050 1150 1200 2400 910 1090 1270 1500 1200 1100 1350 1450 1200 1150 1010 880 1000 1050 1100,2000 880 1060 1230 1390 1500 1500 14
20、00 900 1100 1060 950 870 900 930 950 1600 830 980 1180 1320 1450 1420 1400 1300 700 900 850 840 380 780 750 1200 740 880 1080 1130 1250 1280 1230 1040 900 500 700 780 750 650 550 800 650 760 880 970 1020 1050 1020 830 800 700 300 500 550 480 350 400 510 620 730 800 850 870 850 780 720 650 500 200 30
21、0 350 320 0 370 470 550 600 670 690 670 620 580 450 400 300 100 150 250 y/x 0 400 800 1200 1600 2000 2400 2800 3200 3600 4000 4400 4800 5200 5600,表二,现在利用MATLAB的M文件画出地形图及等高线图,js94.m A=370 470 550 600 670 690 670 620 580 450 400 300 100 150 250 510 620 730 800 850 870 850 780 720 650 500 200 300 350 3
22、20 650 760 880 970 1020 1050 1020 830 800 700 300 500 550 480 350 740 880 1080 1130 1250 1280 1230 1040 900 500 700 780 750 650 550 830 980 1180 1320 1450 1420 1400 1300 700 900 850 840 380 780 750 880 1060 1230 1390 1500 1500 1400 900 1100 1060 950 870 900 930 950 910 1090 1270 1500 1200 1100 1350
23、1450 1200 1150 1010 880 1000 1050 1100 950 1190 1370 1500 1200 1100 1550 1600 1550 1380 1070 900 1050 1150 1200 1430 1450 1460 1500 1550 1600 1550 1600 1600 1600 1550 1500 1500 1550 1550 1420 1430 1450 1480 1500 1550 1510 1430 1300 1200 980 850 750 550 500 1380 1410 1430 1450 1470 1320 1280 1200 108
24、0 940 780 620 460 370 350 1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210 1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150; mesh(A);,A=370 470 550 600 670 690 670 620 580 450 400 300 100 150 250 510 620 730 800 850 870 850 780 720 650 500 200 300 350 320 650 760 8
25、80 970 1020 1050 1020 830 800 700 300 500 550 480 350 740 880 1080 1130 1250 1280 1230 1040 900 500 700 780 750 650 550 830 980 1180 1320 1450 1420 1400 1300 700 900 850 840 380 780 750 880 1060 1230 1390 1500 1500 1400 900 1100 1060 950 870 900 930 950 910 1090 1270 1500 1200 1100 1350 1450 1200 11
26、50 1010 880 1000 1050 1100 950 1190 1370 1500 1200 1100 1550 1600 1550 1380 1070 900 1050 1150 1200 1430 1450 1460 1500 1550 1600 1550 1600 1600 1600 1550 1500 1500 1550 1550 1420 1430 1450 1480 1500 1550 1510 1430 1300 1200 980 850 750 550 500 1380 1410 1430 1450 1470 1320 1280 1200 1080 940 780 62
27、0 460 370 350 1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210 1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150; contour(A),A=370 470 550 600 670 690 670 620 580 450 400 300 100 150 250 510 620 730 800 850 870 850 780 720 650 500 200 300 350 320 650 760 880 970 102
28、0 1050 1020 830 800 700 300 500 550 480 350 740 880 1080 1130 1250 1280 1230 1040 900 500 700 780 750 650 550 830 980 1180 1320 1450 1420 1400 1300 700 900 850 840 380 780 750 880 1060 1230 1390 1500 1500 1400 900 1100 1060 950 870 900 930 950 910 1090 1270 1500 1200 1100 1350 1450 1200 1150 1010 88
29、0 1000 1050 1100 950 1190 1370 1500 1200 1100 1550 1600 1550 1380 1070 900 1050 1150 1200 1430 1450 1460 1500 1550 1600 1550 1600 1600 1600 1550 1500 1500 1550 1550 1420 1430 1450 1480 1500 1550 1510 1430 1300 1200 980 850 750 550 500 1380 1410 1430 1450 1470 1320 1280 1200 1080 940 780 620 460 370
30、350 1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210 1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150; B=rot90(A);C=rot90(B);D=rot90(C); subplot(221),mesh(A);subplot(222),mesh(B);subplot(223),mesh(C); subplot(224),mesh(D);,应用举例2 Taylor级数,我们可以建立函数M文件对不同的n值画出Taylor多项
31、式的图形，并有此看出随着n的增大，Taylor多项式将趋于sinx.,function y=taylor(n) x=-2*pi:pi/100:2*pi; k=length(x); w=zeros(1,k); y=x; p=1;q=x; for j=1:n p=p*(2*j)*(2*j+1); q=q*(-1).*x.2; y=y+q/p; end; z=sin(x); plot(x,w,k,x,z,r,x,y,b);,taylor(6),subplot(221),taylor1(3);subplot(222), taylor1(5);subplot(223), taylor1(6);subpl
32、ot(224), taylor1(8);,傅立叶级数,傅立叶级数为,z(i)=1; elseif x(i) 0 z(i)=-1; elseif x(i)pi z(i)=1; elseif x(i)2*pi z(i)=-1; else z(i)=1; end end plot(x,w,k,x,z,r,x,q,b);,function y=fly(n) x=-3*pi:pi/100:3*pi; k=length(x); w=zeros(1,k); q=w; for j=1:n q=q+4*sin(2*j-1)*x)/(2*j-1)*pi); end; for i=1:k if x(i)=-2*pi
33、 z(i)=-1; elseif x(i)-pi,fly(3),fly(33),subplot(221),fly(3);subplot(222), fly(8);subplot(223), fly(16);subplot(224), fly(30);,function y=djf(n); x=0:0.0001:1; w=x.2; plot(x,w); t=0:1/n:1; hold on y=t.2; bar(t,y) hold off,djf(6);,dif(30);,djf(100);,subplot(221),djf(6);subplot(222), djf(10);subplot(223), djf(20);subplot(224), djf(100);,