f(x)=(x1+10x2)2+5(x3−x4)2+(x2−2x3)4+10(x1−x4)4
我們將這個函式實現為MATLAB程式碼，並使用FR共軛梯度法對其進行最佳化。首先需要定義目標函式及其梯度。然後，使用前面介紹的FR共軛梯度法進行最佳化。
目標函式和梯度的定義
我們需要先定義目標函式 f(x)f(x)f(x) 及其梯度 ∇f(x)\nabla f(x)∇f(x)。
1 function [x, fval, iter] = fr_cg(f, grad_f, x0, tol, max_iter) 2 % f: 待最佳化的目標函式控制代碼 3 % grad_f: 目標函式的梯度控制代碼 4 % x0: 初始點 5 % tol: 終止準則 6 % max_iter: 最大迭代次數 7 8 % 初始化 9 x = x0; 10 g = grad_f(x); 11 d = -g; 12 iter = 0; 13 14 while norm(g) > tol && iter < max_iter 15 % 計算步長alpha 16 alpha = line_search(f, grad_f, x, d); 17 18 % 更新x 19 x_new = x + alpha * d; 20 21 % 計算新的梯度 22 g_new = grad_f(x_new); 23 24 % 計算beta 25 beta = (g_new' * g_new) / (g' * g); 26 27 % 更新搜尋方向 28 d = -g_new + beta * d; 29 30 % 更新x和g 31 x = x_new; 32 g = g_new; 33 34 % 增加迭代次數 35 iter = iter + 1; 36 end 37 38 fval = f(x); 39 end 40 41 function alpha = line_search(f, grad_f, x, d) 42 % 使用簡單的線搜尋方法來確定步長alpha 43 alpha = 1; 44 rho = 0.8; 45 c = 1e-4; 46 47 while f(x + alpha * d) > f(x) + c * alpha * grad_f(x)' * d 48 alpha = rho * alpha; 49 end 50 end
1 % 定義目標函式 2 f = @(x) (x(1) + 10 * x(2))^2 + 5 * (x(3) - x(4))^2 + (x(2) - 2 * x(3))^4 + 10 * (x(1) - x(4))^4; 3 4 % 定義目標函式的梯度 5 grad_f = @(x) [2 * (x(1) + 10 * x(2)) + 40 * (x(1) - x(4))^3; 6 20 * (x(1) + 10 * x(2)) + 4 * (x(2) - 2 * x(3))^3; 7 10 * (x(3) - x(4)) - 8 * (x(2) - 2 * x(3))^3; 8 -10 * (x(3) - x(4)) - 40 * (x(1) - x(4))^3]; 9 10 % 選取初始點 11 initial_points = [0, 0, 0, 0; 1, 1, 1, 1; -1, -1, -1, -1]; 12 13 % 設定終止準則和最大迭代次數 14 tol = 1e-6; 15 max_iter = 1000; 16 17 % 執行實驗並輸出結果 18 for i = 1:size(initial_points, 1) 19 x0 = initial_points(i, :)'; 20 [x, fval, iter] = fr_cg(f, grad_f, x0, tol, max_iter); 21 fprintf('初始點: (%.1f, %.1f, %.1f, %.1f)\n', x0(1), x0(2), x0(3), x0(4)); 22 fprintf('最優解: (%.6f, %.6f, %.6f, %.6f)\n', x(1), x(2), x(3), x(4)); 23 fprintf('最優值: %.6f\n', fval); 24 fprintf('迭代次數: %d\n', iter); 25 fprintf('-----------------------------\n'); 26 end