Bài giảng Tối ưu hoá

ear, clc   f = inline(ʹx(1)^4 ‐ 16*x(1)^2 ‐ 5*x(1) + x(2)^4 ‐ 16*x(2)^2 ‐ 5*x(2)ʹ,ʹxʹ);  l = [‐5  ‐5];   u = [5  5]; %bien duoi/tren  %x0 = [‐0.5   ‐1];  x0 = [0   0];  kmax = 500;   q = 1;   tolfun = 1e‐10;  [xmin, fmin] = simannealing(f, x0, l, u, kmax, q, tolfun)  [xmin1, fmin1] = neldermead(f, x0, 1e‐5, 1e‐8, kmax)  [xmin2, fmin2] = fminsearch(f, x0)  Trong chương  trình  trên,  ta dùng  thêm các  thuật  toán khác để so sánh. Kết  quả là thuật toán SA có thể tìm được cực tiểu toàn cục. Tuy nhiên không phải  khi nào  thuật  toán cũng  thành công. Sự  thành công của  thuật  toán này phụ  thuộc giá trị đầu và may mắn, trong khi các thuật toán khác sự thành công chỉ  phụ thuộc giá trị đầu.  §9. THUẬT TOÁN DI TRUYỀN    Thuật toán di truyền (genetic algorithm ‐ GA) là một kỹ thuật tìm ngẫu  nhiên có định hướng, mô phỏng sự chọn lọc tự nhiên để có các cá thể sống sót  thích nghi nhất. Cũng như  thuật  toán SA, GA  cho phép  tìm  được  cực  tiểu  toàn cục ngay cả khi hàm đối tượng có nhiều cực trị gồm các điểm uốn, các  cực tiểu địa phương, các cực đại địa phương. Thuật toán di truyền lai gồm các  bước: khởi gán, chọn lọc, vượt qua, đột biến. Thuật toán gồm các bước sau:  • Cho giá  trị đầu  [xo] =  [xo1, xo2,...,xoN]  (N  là kích  thước của biến), biên  dưới [l] = [l1,...,lN], biên trên [u] = [u1,...,uN], kích thước của quần thể Np,  vec  tơ Nb =  [Nb1,..., NbN] gồm số bit các gán cho mỗi  thể hiện của mỗi  biến xi, xác suất sống sót Pc, xác suất đột biến Pm, tỉ lệ học η(0≤ η ≤ 1, thể   394 hiện học nhanh hay chậm), số lần lặp cực đại kmax. Chú ý là kích thước  của [xo], [u], [l] đều là N.  • Khởi tạo ngẫu nhiên số cá thể ban đầu của quần thể.   Cho [xo] = [xo], fo = f([xo]) và xây dựng theo cách ngẫu nhiên mảng các cá  thể ban đầu [X1] chứa Np trạng thái(trong vùng cho phép báo bởi biên  trên [u] và biên dưới [l]) bao gồm cả trạng thái ban đầu [xo] bằng ccáh  đặt:    [X1(1)] = [xo], [X1(k)] = [l] + rand.×([u] ‐ [l]) k = 2,..,Np    (1)  Sau  đó mã hoá mỗi  số  của mảng  quần  thể này  thnàh một  chuỗi nhị  phân bằng:  − = = ⎛ ⎞+⎜ ⎟⎝ ⎠∑ ∑ m 1 m 1 bi bi i 1 i 1 P n,1 N : N  = biểu diễn nhị phân của X1(n ,m) với Nbm bit                                     ( ) −= − −bmN 1X (n,m) l(m)2 1 u(m) l(m)       (2)  với n = 1,...,Np và m = 1,...,N  sao cho  toàn  thể mảng quần  thể  trở  thành mảng mà mỗi hàng  là một  nhiễn sắc thể được biểu diễn bằng chuỗi nhị phân  = ∑N bi i 1 N bit.  • Lặp từ k = 1 đến kmax:    ∗ giải mã mỗi số của mảng thành số thập phân bằng:  =kX (n,m) biểu diễn  thập phân  của  − = = ⎛ ⎞+⎜ ⎟⎝ ⎠∑ ∑ m 1 m 1 bi bi i 1 i 1 P n,1 N : N   với Nbm bit                  ( ) −= +−bmk N u(m) l(m)P (n,.) l(m) 2 1       (3)  n = 1,...,N; m = 1,...,N  và tính giá trị f(n) của hàm đối với mỗi hàng Xk(n, :) = [x(n)]  tương ứng với mỗi nhiễm sắc thể và tìm cực tiểu fmin = f(nb)  tương ứng với Xk(n, :) = [x(nb)]  ∗ nếu fmin = f(nb) < fo thì đặt fo = f(nb) và [xo] = [x(nb)]  ∗ biến đổi giá trị của hàm thành giá trị thích hợp bằng:    { }= −pN1 n=1f (n) Max f(n) f(n)           (4)  ∗ nếu  { } ≈pNn=1Max f(n) 0 , kết thúc quá trình và [xo]là giá trị tốt nhất.  Nếu không, để  tạo nhiều nhiễn sắc thể hơn quanh điểm tốt nhất  [x(nb)] cho thế hệ sau, ta dùng quy tắc chọn lọc:  395   [ ] [ ] [ ] [ ]{ }−= + η −1 b 1 b 1 b f (n ) f (n)x(n) x(n) x(n ) x(n) f (n )     (5)  để có được quần thể mới [Xk+1] có Xk+1(n,  :) = [x(n)] và mã hoá nó  để xây dựng mảng Pk+1 mới theo (2)  ∗ xáo trộn chỉ số hàng của mảng Pk+1  ∗ với xác suất tồn tại Pc, thay đổi phần cuối bắt đầu từ vài bit ngẫu  nhiên  của  các  số  trong 2  cặp nhiễm  sắc  thể ngẫu nhiên(hàng  cả  Pk+1)với các nhiễm sắc thể khác để có ma trận  +′k 1P   ∗ với xác suất đột biến Pm, nghịch đảo một bít ngẫu nhiên của mỗi  hàng biểu diễn bởi nhiễm sắc thể (hàng của  +′k 1P ) để tạo ra mảng  Pk+1  Lưu đồ thuật toán như sau:  Ta xây dựng hàm genetic() thực hiên thuật toán trên:  function [xo, fo] = genetic(f, x0, l, u)  Đột biến Khởi gán Đánh giá Nếu giá trị  hàm của các nhiễm sắc  thể bằng nhau  Kết thúc Chọn lọc Vượt qua 396 % Thuat toan Genetic Algorithm tim cuc tieu cua ham f(x) tg doan l <= x <= u  N = length(x0);  kmax = 100; %so lan lap(the he)  eta = 1;%ti le hoc(0 < eta < 1)  Pm = 0.01; %xac suat dot bien  Pc = 0.5; end %xac suat ton tai  Nb = 8*ones(1, N);   Np = 10;  %so nhiem sac the  %khoi gan  NNb = sum(Nb);  xo = x0(:)ʹ;   l = l(:)ʹ;   u = u(:)ʹ;  fo = feval(f, xo);  X(1, :) = xo;  for n = 2:Np      X(n, :) = l + rand(size(x0)).*(u ‐ l); %Pt.(1)  end   P = genencode(X, Nb, l, u); %Pt.(2)  for k = 1:kmax      X = gendecode(P, Nb, l, u); %Pt.(3)      for n = 1:Np          fX(n) = feval(f, X(n,:));       end      [fxb, nb] = min(fX);       if fxb < fo          fo = fxb;           xo = X(nb, :);       end      fX1 = max(fX) ‐ fX; %Pt.(4)      fXm = fX1(nb);      if fXm < eps          return;       end %ket thuc neu tat ca cac nhiem sac the nhu nhau     %Chon loc the h tiep theo      for n = 1:Np  397         X(n, :) = X(n, :) + eta*(fXm ‐ fX1(n))/fXm*(X(nb, :) ‐ X(n, :)); %Pt.(5)      end      P = genencode(X,Nb,l,u);      is = shuffle([1:Np]);      for n = 1:2:Np ‐ 1          if rand < Pc              P(is(n:n + 1), :) = crossover(P(is(n:n + 1), :), Nb);           end      end      %Dot bien      P = mutation(P, Nb, Pm);  end  function X = gendecode(P, Nb, l, u)  % giai ma   Np = size(P, 1);   N = length(Nb);   for n = 1:Np      b2 = 0;      for m = 1:N          b1 = b2 + 1;           b2 = b1 + Nb(m) ‐ 1; %Pt.(3)          X(n, m) = bin2dec(P(n,b1:b2))*(u(m) ‐ l(m))/(2^Nb(m) ‐ 1) + l(m);      end  end  function P = genencode(X, Nb, l, u)  Np = size(X,1);   N = length(Nb);   for n = 1:Np      b2 = 0;      for m = 1:N          b1 = b2 + 1;           b2 = b2 + Nb(m);          Xnm = (2^Nb(m)‐ 1)*(X(n, m) ‐ l(m))/(u(m) ‐ l(m)); %Pt.(2)          P(n, b1:b2) = dec2bin(Xnm, Nb(m));   398     end  end  function chrms = crossover(chrms2, Nb)  Nbb = length(Nb);  b2 = 0;  for m = 1:Nbb      b1 = b2 + 1;       bi = b1 + mod(floor(rand*Nb(m)), Nb(m));       b2 = b2 + Nb(m);      tmp = chrms2(1, bi:b2);      chrms(1, bi:b2) = chrms(2, bi:b2);      chrms(2, bi:b2) = tmp;  end  function P = mutation(P, Nb, Pm)   Nbb = length(Nb);  for n = 1:size(P,1)      b2 = 0;      for m = 1:Nbb          if rand < Pm              b1 = b2 + 1;               bi = b1 + mod(floor(rand*Nb(m)),Nb(m));               b2 = b2 + Nb(m);              P(n,bi) = ~P(n,bi);          end      end  end  function is = shuffle(is)   N = length(is);  for n = N:‐1:2      in = ceil(rand*(n ‐ 1));       tmp = is(in);      is(in) = is(n);       is(n) = tmp;   end  399 Để tìm cực tiểu của hàm ta dùng chương trình ctgenetic.m:  clear all, clc  f = inline(ʹx(1).^2 + 2*x(2).^2ʹ);  l = [‐5  ‐5 ];   u = [5  5]; %bien duoi/tren  x0 = [0   0];  [xmin, fmin] = genetic(f, x0, l, u)  §10. THUẬT TOÁN FIBONACCI  Trong thuật toán tỉ lệ vàng, hai lần tính giá trị của hàm được thực hiện  tại lần lặp đầu tiên và sau đó chỉ tính giá trị hàm một lần trong các lần lặp tiếp  theo. Giá trị của r là hằng số trong mỗi đoạn con và việc tìm điểm cực tiểu kết  thúc tại đoạn con thứ k có  − < δk ka b  hay  − < εk kf(b f(a ) . Phương pháp tìm  theo thuật toán Fibonacci khác phương pháp tỉ lệ vàng ở chỗ r không phải là  hằng số trên mỗi đoạn con. Ngoài ra số đoạn con (số bước lặp) được xác định  trước. Thuật toán tìm Fibonacci dựa trên dãy số Fibonacci được xác định bằng  phương trình:    fo = 0    f1 = 1    fn = fn‐1 + fn‐2 với n = 2,3,...  Như vậy dãy số Fibonacci là: 0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55,...     Giả sử ta có hàm f(x) có cực tiểu trên đoạn [a, b]. Như ở trong phương  pháp tỉ lệ vàng, 0.5 < ro < 1 được chọn sao cho cả hai điểm bên trong co và do  sẽ được dùng  trong  đoạn  con  tiếp  theo và như vậy  chỉ  cần  tính giá  trị  của  hàm một lần. Nếu f(co) < f(do) thì cực tiểu nằm trong đoạn [ao, do] và ta thay  =1 oa a  và  =1 ob d và tiếp tục tìm trong khoảng mới  [ ] [ ]=1 1 o oa ,b a ,d . Nếu f(co)  > f(do) thì cực tiểu nằm trong đoạn [co, bo] và ta thay a1 = co và b1 = bo và tiếp  tục tìm trong khoảng mới [ ] [ ]=1 1 o oa ,b c ,b  như hình vẽ.     ao  eo  co  do  bo ao eo co do  bo  400 Nếu f(co) < f(do) và chỉ muốn tính giá trị của hàm một lần trong đoạn [ao, bo] ta  sẽ chọn r1 (0.5 < r1 < 1) trong đoạn con [ ] [ ]=1 1 o oa ,b a ,b . Ta đã kí hiệu b1 = do và  do co ∈ [ao, do] nên ta có:    do ‐ co = b1 ‐ d1                  (1)  Tỉ số ro được chọn sao cho do ‐ ao = ro(bo ‐ ao) và co ‐ ao = (1 ‐ro)(bo ‐ ao) và thay  thế:    do ‐ co = (do ‐ ao) ‐ (co ‐ ao)    do ‐ co = ro(bo ‐ ao) ‐ (1 ‐ ro)(bo ‐ ao)    do ‐ co = (2ro ‐ 1)(bo ‐ ao)               (2)  và r1 được chọn sao cho:    b1 ‐ d1 = (1 ‐ r1)(b1 ‐ a1)                (3)  Thay (2) và (3) vào (1) ta có:    (2ro ‐ 1)(bo ‐ ao) = (1 ‐ r1)(b1 ‐ a1)            (4)  Như vậy đoạn [a, b] bị co ngắn bằng hệ số ro và (b1 ‐ a1) = ro(bo ‐ ao) và:     (2ro ‐ 1)(bo ‐ ao) = (1 ‐ r1)ro(bo ‐ ao)             (5)  Rút gọn ta có:    (2ro ‐ 1) = (1 ‐ r1)ro                  (6)  Từ (6) ta tính được r1:    −= o1 o 1 rr r                     (7)  Trong (7), thay  −= n 1o n fr f  ta có:  − − − − − − − −= = = n 1 n n 1 n 2n 1 n 1 n 1 n 1 n f1 f f ffr f f f f Ta rút ra rằng thuật toán tìm Fibonacci có thể bắt đầu bằng:    −= n 1o n fr f   − − = n 21 n 1 fr f và:  − − − = n 1 kk n k fr f , k = 1, 2,..., n ‐ 3  Bước cuối cùng là:    − = =2n 3 3 f 1r f 2 401 Thuật  toán  tìm Fibonacci gồm  (n  ‐ 2)  lần  tính. Đoạn  con  thứ  (k+1)  có  được  bằng cách giảm độ dài của đoạn  thứ k bằng hệ số  − − − = n 1 kk n k fr f . Sau  (n  ‐ 2)  lần  tính, độ dài của bước cuối cùng là:    − − − − − − = − = −Ln 1 n 2 n 3 3 3 2o o o o o o n n 1 n 2 4 3 n n f f f f f f 1(b a ) (b a ) (b a ) f f f f f f f Nếu sai số cho trước là ε, nghĩa là  − < εo o n 1 (b a ) f  và cần dùng n lần lặp với n là  số nguyên nhỏ nhất sao cho:    −> ε o o n b af                    (8)  Các điểm bên trong ck và dk được xác định bằng:    − − − ⎛ ⎞= + + −⎜ ⎟⎝ ⎠ n 1 k k k k k n k fc a 1 (b a ) f               (9)    − − − = + + −n 1 kk k k k n k fd a 1 (b a ) f               (10)  Ta xây dựng hàm fibonacci() để thực hiện thuật toán trên:  function [x, y] = fibonacci(f, a, b, n)  % Phuong phap Fibonacci de tim cuc tieu cua   % ham f trong (a, b) voi n buoc tinh  fn2 = 1;   fn1 = 1;   fn = fn1 + fn2;  for i = 3:n      fn2 = fn1;       fn1 = fn;       fn = fn1 + fn2;  end  l = (b ‐ a)*fn2/fn;  x1 = a + l;   x2 = b ‐ l;  f1 = feval(f, x1);   f2 = feval(f,x2);  fibn = fn;   ll1 = b ‐ a;  402 for j = 3:n      llj = ll1*fn2/fibn;      fn = fn1;       fn1 = fn2;       fn2 = fn ‐ fn1;      if f2 > f1         b = x2;          l = (b ‐ a)*fn2/fn;         x2 = x1;          x1 = a + l;          f2 = f1;          f1 = feval(f, x1);      else         a = x1;          l = (b ‐ a)*fn2/fn;         x1 = x2;          x2 = b ‐ l;          f1 = f2;         f2 = feval(f, x2);      end  end  x = x1; y = f1;  Để tìm cực tiểu của hàm trong đoạn (a, b) ta dùng chương trình ctfibonacci.m:  clear all, clc  f = inline(ʹ1.6*x^2 ‐ 3*x + 2ʹ);  a = ‐0.;   b = 1;  n = 50;  [x, y] = fibonacci(f, a, b, n)