程式流程控制 方煒 台大生機系.

Presentation on theme: "程式流程控制 方煒 台大生機系."— Presentation transcript:

1 程式流程控制 方煒 台大生機系

2 迴圈指令 MATLAB 提供兩種迴圈指令： for 迴圈 while 迴圈

3 for 迴圈 for 迴圈的使用語法如下： 運算式 end 另一種 for 迴圈的使用語法如下： for 變數 = 向量,

4 程式流程控制- for loop 下列 for 迴圈會產生一個長度為 6 的調和數列（Harmonic Sequence）：
範例 1 : forLoop01.m x = zeros(1,6); % 變數 x 是一個 1×6 大小的零矩陣 for i = 1:6 x(i) = 1/i; end x % 顯示 x x = 上例中，矩陣 x 最初是一個 1×6 大小的零矩陣，在 for 迴圈中，變數 i 的值依次是 1 到 6，因此矩陣 x 的第 i 個元素的值依次被設為 1/i。 可用分數形式來顯示此數列： >> format rat >> disp(x) / / / / /6

5 For k=m:s:n Statements End

6 練習 Summing a series with a for loop
% script file: summing.m s=0; % set a variable to 0 so that 1/n^2 can be repeatedly added to it N=10000; % set the upper limit of the sum for n=1:N % add 1/n^2 to s each time, then put the answer back into s s=s+1/n^2; end fprintf(’ Sum = %g \n’,s) % print the answer % calculate the sum of the squares of the reciprocals of the % integers from 1 to 10,000 n=1:10000; sum(1./n.^2)

7 練習 Products with a for loop
% script file: product.m P=1; % set the first term in the product N=20; % set the upper limit of the product for n=2:N P=P*n; end fprintf(’ N! = %g \n’,P) % print the answer factorial(20) gamma(21)

8 練習 Recursion relations with for loops
%script file: recursion.m a(1)=1; % put the first element into the array N=19; % the first one is loaded, so let’s load 19 more for n=1:N a(n+1)=(2*n-1)/(2*n+1)*a(n); End format short disp(a) % display the resulting array of values format rat disp(a)

9 程式流程控制- for loop for 迴圈可以是多層或巢狀式（Nested）的，在下例中即產生一個 6×6 的Hilbert 矩陣 h，其中位於第 i 列、第 j 行的元素為 ： 範例 2 : forLoop02.m h = zeros(6); % 變數 x 是一個 6×6 大小的零矩陣 for i = 1:6 for j = 1:6 h(i,j) = 1/(i+j-1); end format rat h h = / / / / /6 1/ / / / / /7 1/ / / / / /8 1/ / / / / /9 1/ / / / / /10 1/ / / / / /11

10 程式流程控制- for loop 在下例中，for 迴圈列出先前產生的 Hilbert 矩陣的每一直行的平方和：
範例 3 : forLoop03.m format short % 回到預設形式來顯式所有數值 for i = h disp(norm(i)^2); % 印出每一行的平方和 End 1.4914 0.5118 0.2774 0.1787 0.1262 0.0944 在上例中，由於 h 是一個矩陣，因此每一次 i 的值就是矩陣 h 的一直行的內容。 for i=1:6,disp(sum(h(i,:).^2)), end

11 程式流程控制- for loop 若要跳出 for 迴圈，可用 break 指令。例如，若要找出最小的 n 值，滿足 n! > 10100，可輸入如下： 範例 4 : break01.m for i = 1:1000 if prod(1:i) > 1e100 fprintf('%g! = %e > 1e100\n', i, prod(1:i)); break; % 跳出 for 迴圈 end 70! = e+100 > 1e100

12 程式流程控制- for loop 在一個迴圈內若要直接跳至到此迴圈下一回合的執行，可使用 continue 指令。
範例 5 : continue01.m x = [ ]; posTotal = 0; for i = 1:length(x) if x(i)<0, continue; end % 若 x(i) 小於零，跳到此迴圈的下一回合 posTotal=posTotal+x(i); end posTotal % 顯示 posTotal 的值 posTotal = 9

13 while 迴圈 while 迴圈使用語法如下： while 條件式 運算式； end

14 程式流程控制– while loop 先前產生調和數列的例子，亦可用 while 迴圈改寫如下： 範例 6 : while01.m
x = zeros(1,6); i = 1; while i<=6 x(i) = 1/i; i = i+1; end x % 顯示 x x =

15 程式流程控制– while loop 若要用 while 指令找出最小的 n 值，使得 n! > 10100 ，可輸入如下：
範例 7 : while02.m n = 1; while prod(1:n) < 1e100 n = n+1; end fprintf('%g! = %e > 1e100\n', n, prod(1:n)); 70! = e+100 > 1e100 與前述的 for 迴圈相同，在任何時刻若要跳出 while 迴圈，亦可使用 break 指令；若要跳到下一回合的 while 迴圈，也可以使用 continue 指令。 無論是 for 或 while 迴圈，均會降低 MATLAB 的執行速度，因此盡量使用向量化的運算（Vectorized Operations）而盡量少用迴圈。 break 指令若用在多重迴圈中，每次只跳出包含break指令的最內部迴圈。

16 練習 x = 5;k = 0; while x < 25 k = k + 1; y(k) = 3*x; x = 2*x-1; end

17 練習 total = 0;k = 0; while total < 1e+4 k = k + 1;
total = 5*k^2 - 2*k + total; end disp(’The number of terms is:’) disp(k) disp(’The sum is:’) disp(total)

18 條件/分支指令 二種條件指令（Branching Command）: if-then-else 條件指令
switch - case - otherwise 條件指令

19 if-then-else條件指令 if - then - else，其使用語法： if 條件式 運算式一; else 運算式二; end
當條件式成立時，執行運算式一，否則，就執行運算式二。 若不使用運算式二，則可省略 else 和運算式二。

20 程式流程控制- if-then-else 在數值運算的過程中，若變數值為 NaN（即 Not A Number ）時，我們要立刻印出警告訊息，可輸入如下例： 範例8 : if01.m x = 0/0; if isnan(x) disp('Warning: NaN detected!'); end Warning: Divide by zero. … Warning: NaN detected! 在上例中，第一個警告訊息是 MATLAB 自動產生的，第二個警告訊息則是我們的程式碼產生的，其中 isnan(x) 可用於判斷 x 是否為 NaN，若是，則傳回 1（真），否則即傳回 0（偽）。

21 程式流程控制- if-then-else 可根據向量 y 的元素值為奇數或偶數，來顯示不同的訊息： 範例9 : if02.m
for i = 1:length(y) if rem(y(i), 2)==0 fprintf('y(%g) = %g is even.\n', i, y(i)); else fprintf('y(%g) = %g is odd.\n', i, y(i)); end y(1) = 0 is even. y(2) = 3 is odd. y(3) = 4 is even. y(4) = 1 is odd. y(5) = 6 is even. 上述的 if - then - else 為雙向條件，亦即程式只會執行「運算式一」或「運算式 二」，不會有第三種可能。

22 練習 %script file: logic.m clear; a=1;b=3;
% If the number a is positive set c to 1; if a is 0 or negative set c to 0 if a>0 c=1 else c=0 end % if either a or b is non-negative, add them to obtain c; % otherwise multiply a and b to obtain c if a>=0 | b>=0 % either non-negative c=a+b c=a*b % otherwise multiply them to obtain c

23 練習 x = [4,-9,25]; if x < 0 disp(’Some elements of x are negative.’)
else y = sqrt(x) end 結果會是甚麼？

24 練習 x = [4,-9,25]; if x >= 0 y = sqrt(x) else
disp(’Some elements of x are negative.’) end 結果會是甚麼？

25 elseif if logical expression 1 statement group 1
elseif logical expression 2 statement group 2 else statement group 3 end

26 練習 if x > 10 y = log(x) elseif x >= 0 y = sqrt(x) else
y = exp(x) - 1 end

27 練習

28 程式流程控制- if-then-else 若要進行更多向的條件，只需一再重覆 elseif 即可。例如，欲判斷 y 向量之元素是屬於 3n、3n＋1、 或 3n+2，可輸入如下： 範例10 : if03.m y = [ ]; for i = 1:length(y) if rem(y(i),3)==0 fprintf('y(%g)=%g is 3n.\n', i, y(i)); elseif rem(y(i), 3)==1 fprintf('y(%g)=%g is 3n+1.\n', i , y(i)); else fprintf('y(%g)=%g is 3n+2.\n', i , y(i)); end y(1)=3 is 3n. y(2)=4 is 3n+1. y(3)=5 is 3n+2. y(4)=9 is 3n. y(5)=2 is 3n+2.

29 練習 陣列 A中不小於0的元素取平方根，其他元素則加上50，求新的陣列 A =[0 -1 4;9 -14 25;-34 49 64];
for m=1:size(A,1) for n = 1:size(A,2) if A(m,n)>=0 B(m,n)=sqrt(A(m,n)); else B(m,n)=A(m,n)+50; end A =[0 -1 4; ; ]; C=(A>=0); A(C)=sqrt(A(C); A(~C)=A(~C)+50;

30 Implied loop 內涵迴圈 for k=1:21 x=(k-1)*5 x=[0:5:100]; y(k)=cos(x); 等效於
end x=[0:5:100]; y=cos(X) 等效於 x=[-10,-5,0,5,10]; m=0; for k=1:length(x) If x(k)>0 m=m+1; y(m)=k; end x=[-10,-5,0,5,10]; y=find(x>0) 等效於

31 switch-case-otherwise 指令
MATLAB 在第五版開始支援 switch-case-otherwise 的多向條件指令，其使用語法如下： switch expression case value(1) statement(1) case value(2) statement(2) case value(n-1) statement(n-1) otherwise statement(n) end 在上述語法中，expression 為一數值或字串，當其值和 value(k) 相等時，MATLAB 即執行 statement(k) 並跳出 switch 指令。若 expression 不等於 value(k)，k=1, 2, …, n-1，則 MATLAB 會執行 statement(n) 並跳出 switch 指令。

32 程式流程控制 – switch case 欲根據月份來判斷其季別，可輸入如下： 範例11 : switch01.m
for month = 1:12 switch month case {3,4,5} season = 'Spring'; case {6,7,8} season = 'Summer'; case {9,10,11} season = 'Autumn'; case {12,1,2} season = 'Winter'; end fprintf('Month %d ===> %s.\n', month, season); Month 1 ===> Winter. Month 12 ===> Winter.

33 程式流程控制 - switch case 如果 expression 是字串，那麼若要在 case 之後比對多個字串，就必需使用字串的異值陣列（Cell Array of Strings）： 範例12 : switch02.m month = {'Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep'}; for i = 1:length(month) switch month{i} case {'Mar','Apr','May'} season = 'Spring'; case {'Jun','Jul','Aug'} season = 'Summer'; case {'Sep','Oct','Nov'} season = 'Autumn'; case {'Dec','Jan','Feb'} season = 'Winter'; end fprintf('%s is %s.\n', month{i}, season);

34 程式流程控制- switch case 上述範例output如下: Jan is Winter. Feb is Winter.
Mar is Spring. Apr is Spring. May is Spring. Jun is Summer. Jul is Summer. Aug is Summer. Sep is Autumn. MATLAB 的 switch 指令和 C 語言的 switch 指令略有差別：在 C 語言的 switch 敘述內，每個 case 敘述需加上 break 以跳出該 switch 敘述，而在 MATLAB 則不必多此一舉。 一般而言，switch–case–otherwise 的執行效率優於 if–then–else 。

35 補充 - Using fzero function f=fz(x) % evaluate the function fz(x) whose
% roots are being sought f=exp(-x)-x; %**************************************** % Here is the matlab code that uses fz.m to find % a zero of f(x)=0 near the guess x=.7 % Note that sign is used to tell Matlab that % the name of an M-file is being passed into fzero %***************************************

36 補充：inline, Example 1 This example creates a simple inline function to square a number. >>g = inline('t^2') g = Inline function: g(t) = t^2 >>g(3) ans = 9 You can convert the result to a string using the char function. >>char(g) t^2

37 補充：inline, Example 2 f = inline('3*sin(2*x.^2)') argnames(f)
Inline function: f(x) = 3*sin(2*x.^2) argnames(f) ans = 'x' formula(f) 3*sin(2*x.^2)

38 補充：inline, Example 3 >>f = inline('sin(alpha*x)')
Inline function: f(alpha,x) = sin(alpha*x) If inline does not return the desired function variables or if the function variables are in the wrong order, you can specify the desired variables explicitly with the inline argument list. >>g = inline('sin(alpha*x)','x','alpha') g = g(x,alpha) = sin(alpha*x)

39 補充 – inline , Example 4 func=inline(’exp(-x)-x’,’x’); x=0:.01:1;
f=func(x); plot(x,f,’r-’)

40 切線法求根

41 切線法 (Secant method) 求根 % set chk, the error, to 1 so it won’t trigger
% script file: secant01.m clear;close all; %************************************ % Define the function as an in line function func=inline(’exp(-x)-x’,’x’); % First plot the function x=0:.01:2; f=func(x); plot(x,f,’r-’,x,0*x,’b-’) % From the plot the solution is near x=.6 % Secant method to solve the exp(-x)-x = 0 % Use an initial guess of x1=0.6 x1=0.6; % find f(x1) f1=func(x1); % find a nearby second guess x2=0.99*x1; % set chk, the error, to 1 so it won’t trigger % the while before the loop starts chk=1; % start the loop while chk>1e-8 % find f(x2) f2=func(x2); % find the new x from the straight line approximation and print it xnew = x2 - f2*(x2-x1)/(f2-f1) % find chk the error by seeing how closely f(x)=0 is approximated chk=abs(f2); % load the old x2 and f2 into x1 and f1; then put the new x into x2 x1=x2;f1=f2;x2=xnew; end

42 切線法 (Secant method) 求根 script file: secant02.m
func=inline('exp(x)-3*x','x'); x=0.6159 x=1.5121

43 function function showNearestPoint1(action)
global h0 h1 h2 h3 h4 x1 y1 if nargin==0 action='b'; end x1=linspace(0,4*pi,30); y1=sin(x1); plot(x1,y1,'r-o') axis equal %

44 switch(action) case 'b' set(gcf,'windowbuttonmotionfcn','showNearestPoint1 a') case 'a' currPt=get(gca,'CurrentPoint'); x2=currPt(1,1); y2=currPt(1,2); Minx=sqrt((x2-x1).^2+(y2-y1).^2); n=find(Minx==min(Minx)); line(x1(n),y1(n),'marker','o') location=[num2str(x1(n)) ',' num2str(y1(n))]; %text(x1(n),y1(n),loction);%顯示在點上 text(x2+.2,y2,location);%顯示在游標旁 end