Bài giảng Tối ưu hóa trong thiết kế cơ khí - Chương 11: Ứng dụng Matlab giải các bài toán tối ưu hóa

Giới thiệu về Optimization Toolbox 10

STT Loại bài toán Hàm sử

dụng

1 Cực tiểu hóa hàm 1 biến số

(Scalar Minimization)

fminbnd

2 Cực tiểu hóa hàm nhiều biến số không có ràng buộc

(Unconstrained Minimization)

fminunc

fminsearch

3 Quy hoạch tuyến tính

(Linear Programming)

linprog

4 Quy hoạch bậc hai

(Quadratic Programming)

quadprog

5

Cực tiểu hóa hàm phi tuyến với các ràng buộc

tuyến tính và phi tuyến

(Constrained Minimization)

fmincon

6 Cực tiểu hóa nửa vô hạn

(Semi-Infinite Minimization) fseminf

pdf 51 trang kimcuc 10001
Bạn đang xem 20 trang mẫu của tài liệu "Bài giảng Tối ưu hóa trong thiết kế cơ khí - Chương 11: Ứng dụng Matlab giải các bài toán tối ưu hóa", để tải tài liệu gốc về máy hãy click vào nút Download ở trên

Tóm tắt nội dung tài liệu: Bài giảng Tối ưu hóa trong thiết kế cơ khí - Chương 11: Ứng dụng Matlab giải các bài toán tối ưu hóa

Bài giảng Tối ưu hóa trong thiết kế cơ khí - Chương 11: Ứng dụng Matlab giải các bài toán tối ưu hóa
Trường Đại học Công nghiệp thành phố Hồ Chí Minh 
Khoa Công nghệ Cơ khí 
CHƯƠNG 11: 
ỨNG DỤNG MATLAB GIẢI CÁC BÀI 
TOÁN TỐI ƯU HÓA 
Thời lượng: 3 tiết 
2 
Thống nhất phiên bản MATLAB 
Trường câu lệnh 
Các biến 
tính toán 
xong 
Trường lịch 
sử các câu 
lệnh đã 
dùng 
Chọn thư 
mục làm 
việc 
3 
Làm quen với MATLAB 
Dấu >> Ký tự nhắc mặc định trong MATLAB, đầu câu lệnh >> 
Dấu ; 
Dấu chấm phẩy ở cuối dòng sẽ tránh việc in kết quả ra trường 
câu lệnh (khi ta không cần in kết quả ra cho cửa sổ ngắn gọn) 
>> a=2; 
Dấu  
Dấu ba chấm ở cuối dòng cho phép tiếp tục code ở dòng tiếp 
theo 
>> a= 
2 
help 
tên_câu_lệnh 
Hiển thị các thông tin chi tiết về câu lệnh mà người dùng cần 
sử dụng 
>>help linprog 
Chữ cái viết thường và viết hoa được phân biệt khác nhau trong MATLAB 
>> a=2 
>>A=2 
>>A+a 
MATLAB coi tất cả các biến đều ở dạng mảng (arrays) 
Tên biến 
Trong MATLAB tên biến được bắt đầu từ chữ cái và có chiều 
dài tối đa 31 ký tự bao gồm chữ cái (in hoa và viết thường là 
khác nhau), số và dấu gạch dưới 
>>A_a_bc_9=12 
Các phép 
toán thông 
dụng 
+ 
- 
* 
/ 
^ 
>>2+3*6^2/4-7 
Tránh trùng 
với các tên 
biến tích hợp 
của hệ thống, 
như các hằng 
số, tên hàm 
pi 
sin 
cos 
v.v.. 
>>pi 
4 
Ma trận trong MATLAB 
Đối tượng Cách thực hiện Ví dụ 
Véctơ hàng Ngoặc vuông, Dấu cách A=[1 2 3 4] 
Véctơ cột 
Cách 1: Ngoặc vuông và xuống dòng 
B=[5 
2 
] 
Cách 2: Ngoặc vuông và dấu chấm phẩy B=[5;2;3] 
Cách 3: Đảo véc tơ hàng thành cột bằng dấu ‘ B=[5 2 3]' 
Ma trận [mxn] 
Cách 1: áp dụng véc tơ hàng và cột: dấu cách 
và xuống dòng 
A=[1 2 3 
4 5 6 
7 8 9] 
Cách 2: Dấu cách và dấu chấm phẩy A=[1 2 3;4 5 6; 7 
8 9] 
Ma trận 1 đơn vị ở 
đường chéo 
- eye(m,n) – ma trận mxn 
- eye(n) – ma trận vuông kích thước n 
eye(4) 
eye(3,4) 
Ma trận toàn 1 ones(m,n) hoặc ones(n) ones(3) 
Ma trận toàn 0 zeros(m,b) hoặc zeros(n) zeros(4) 
Tạo ra 1 dãy cấp số cộng M:icr:N 
M:N – khi icr=1 theo mặc định 
100:-7:50 
50:100 
5 
Ma trận trong MATLAB (tiếp) 
Đối tượng Cách thực hiện Ví dụ 
Phần tử trong ma trận Dấu ngoặc tròn: 
- A(i,j) chọn phần tử hàng i, cột j của ma 
trận A 
- A(m1:m2,n1:n2) Chọn các phần tử từ 
hàng m1 đến hàng m2, cột n1 đến cột n2 
A(2,3) 
A(2:3,1:3) 
6 
Kịch bản (Scripts) 
Là một dạng M-file đơn giản nhất, không có biến vào và biến ra. Nó 
chỉ gồm một chuỗi các trình tự câu lệnh. 
7 
Hàm số (Function) 
function [O_1,O_2,,O_m] = function_name (I_1,I_2,,I_n) 
%----------------------------------------------------------------------- 
% Ở dưới ghi function body 
end 
8 
function [ p ] = func1( x ) 
%func1 of this function goes here 
%Detailed explanation goes here 
p=x^3-2*x+cos(x); 
end 
9 
function [x1,x2] = PTB2(a,b,c) 
% Ham PTB2 dung de giai phuong trinh bac 2 co dang a*x^2+b*x+c=0 
% Tham bien dau vao la 3 he so a, b, c trong do a !=0 
Delta = b^2-4*a*c; 
x1 = (-b+sqrt(Delta))/(2*a); 
x2 = (-b-sqrt(Delta))/(2*a); 
[x1 x2] 
end 
[x1, x2] = feval('PTB2',1,2,-5) 
10 Giới thiệu về Optimization Toolbox 
STT Loại bài toán Hàm sử 
dụng 
1 
Cực tiểu hóa hàm 1 biến số 
(Scalar Minimization) 
fminbnd 
2 
Cực tiểu hóa hàm nhiều biến số không có ràng buộc 
(Unconstrained Minimization) 
fminunc 
fminsearch 
3 
Quy hoạch tuyến tính 
(Linear Programming) 
linprog 
4 
Quy hoạch bậc hai 
(Quadratic Programming) 
quadprog 
5 
Cực tiểu hóa hàm phi tuyến với các ràng buộc 
tuyến tính và phi tuyến 
(Constrained Minimization) 
fmincon 
6 
Cực tiểu hóa nửa vô hạn 
(Semi-Infinite Minimization) 
fseminf 
11 Các thuật toán của các công cụ 
12 Các thuật toán của các công cụ 
13 
CỰC TIỂU HÓA HÀM MỘT BIẾN SỐ 
 min;f x a x b 
[x,fval,exitflag,output] = fminbnd(@Objfun,a,b,options) 
x – xuất ra giá trị x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
• exitflag=1: có nghĩa là hàm đã hội tụ tại điểm lời giải x nếu Stopping Criteria = 'TolX‘ 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu hoặc số lượng vòng lặp đã đạt 
ngưỡng cho phép nếu Stopping Criteria = 'MaxIter' hoặc 'MaxFunEvals' 
• exitflag=-1: có nghĩa là thuật toán bị dừng vì hàm đầu ra 
• exitflag=-2: có nghĩa là khoảng giá trị bị sai (a>b) 
output – xuất ra các thông tin về số vòng lặp tính toán, số lần tính hàm số, các thuật toán 
tại các bước tính và thông báo cuối cùng 
Objfun – tên của M-file xác định hàm mục tiêu 
a,b – giá trị 2 biên của biến x 
options – các thuộc tính cần thiết cho việc giải bài toán bằng hàm fminbnd, được xác định 
trước dòng cú pháp 
14 
CỰC TIỂU HÓA HÀM MỘT BIẾN SỐ 
Tìm cực tiểu hàm số  2
0.75 1
0.65 0.65 arctan ; 0;0.5
1
f x x x
x x
1) Bước 1: Tạo 1 thư mục cho bài toán, ví dụ fminbnd1 
Vào trong 
thư mục 
15 2) Bước 2: Tạo M-file hàm mục tiêu, ví dụ Objfun.m 
function f = Objfun( x ) 
%OBJFUN Summary of this function goes here 
% Detailed explanation goes here 
f= 0.65 - (0.75/(1+x^2))- 0.65*x*atan(1/x); 
end 
Sửa lại code 
của hàm 
16 3) Bước 3: Tạo M-file Script lời giải, ví dụ Solve.m 
clear;clc;format long;warning('off'); 
% Nhap mien xac dinh cua bien x: a<=x<=b 
a=0; 
b=0.5; 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
Alg='interior-point'; 
%Alg='levenberg-marquardt'; 
%Alg='trust-region-dogleg'; 
%Alg='lm-line-search'; 
% Chon Stopping Criteria 
GTN=1e-4;GTL=1e4; 
% Neu Stopping Criteria la sai so cua tham bien 
StCr='TolX'; 
% Neu Stopping Criteria la so luong Iterations 
%StCr='MaxIter'; 
% Neu Stopping Criteria la so luong tinh cac ham so 
%StCr='MaxFunEvals'; 
if strcmp(StCr,'MaxIter')==1 || strcmp(StCr,'MaxFunEvals')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 
 GT=GTN; 
end 
options = 
optimset('Algorithm',Alg,'Display','iter',StCr,GT,'PlotFcns',@optimplotfval); 
[x,fval,exitflag,output]=fminbnd(@Objfun,a,b,options) 
17 
0 1 2 3 4 5 6 7 8 9
-0.32
-0.31
-0.3
-0.29
-0.28
-0.27
-0.26
-0.25
-0.24
-0.23
Iteration
F
u
n
c
ti
o
n
 v
a
lu
e
Current Function Value: -0.31002
18 
CỰC TIỂU HÓA HÀM NHIỀU BIẾN SỐ KHÔNG 
RÀNG BUỘC 
 1 2min; , , , nf x x x x x
[x,fval,exitflag,output, grad,hessian] = fminunc(@Objfun,x0,options) 
x – xuất ra véctơ tham biến x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
• exitflag>0: có nghĩa là hàm đã hội tụ tại điểm lời giải x 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu hoặc số lượng vòng lặp đã đạt 
ngưỡng cho phép nếu Stopping Criteria = 'MaxIter' hoặc 'MaxFunEvals' 
• exitflag<0: Hàm không hội tụ tại điểm lời giải 
output – xuất ra các thông tin về số vòng lặp tính toán, số lần tính hàm số, các thuật toán 
tại các bước tính, v.v 
grad – xuất ra véctơ Gradient của hàm số tại điểm lời giải 
hessian – xuất ra ma trận Hessian của hàm số tại điểm lời giải 
x0 – véc tơ tham biến khởi đầu 
Objfun - tên của M-file xác định hàm mục tiêu 
19 
Tìm cực tiểu hàm số: 1 2 1 2
1 2
50 20
, minf f x x x x
x x
 x
1) Bước 1: Tạo 1 thư mục cho bài toán, ví dụ MultiUncon1 
2) Bước 2: Tạo M-file hàm mục tiêu, ví dụ Objfun.m 
function f = Objfun( x ) 
%OBJFUN Summary of this function goes 
here 
% Detailed explanation goes here 
f = x(1)*x(2) + 50/x(1) + 20/x(2); 
end 
3) Bước 3: Tạo M-file Script lời giải, ví dụ Solvefminunc.m 
20 clear;clc;format long;warning('off'); 
% Nhap vecto tham bien khoi dau: 
x0 = [1,1]; 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
%Alg='interior-point'; 
%Alg='levenberg-marquardt'; 
Alg='trust-region-dogleg'; 
%Alg='lm-line-search'; 
% Chon Stopping Criteria 
GTN=1e-4;GTL=1e4; 
% Neu Stopping Criteria la sai so cua tham bien 
%StCr='TolX'; 
% Neu Stopping Criteria la sai so cua ham so 
StCr='TolFun'; 
% Neu Stopping Criteria la so luong Iterations 
%StCr='MaxIter'; 
% Neu Stopping Criteria la so luong tinh cac ham so 
%StCr='MaxFunEvals'; 
if strcmp(StCr,'MaxIter')==1 || strcmp(StCr,'MaxFunEvals')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 || strcmp(StCr,'TolFun')==1 
 GT=GTN; 
end 
options = 
optimset('Algorithm',Alg,'Display','iter',StCr,GT,'PlotFcns',@optimplotfval); 
[x,fval,exitflag,output, grad,hessian] = fminunc(@Objfun,x0,options) 
21 
22 
0 2 4 6 8 10 12 14
30
35
40
45
50
55
60
65
70
75
Iteration
F
u
n
c
ti
o
n
 v
a
lu
e
Current Function Value: 30
23 
CỰC TIỂU HÓA HÀM NHIỀU BIẾN SỐ KHÔNG 
RÀNG BUỘC 
 1 2min; , , , nf x x x x x
[x,fval,exitflag,output] = fminsearch(@Objfun,x0,options) 
x – xuất ra véctơ tham biến x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
• exitflag>0: có nghĩa là hàm đã hội tụ tại điểm lời giải x 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu hoặc số lượng vòng lặp đã đạt 
ngưỡng cho phép nếu Stopping Criteria = 'MaxIter' hoặc 'MaxFunEvals' 
• exitflag<0: Hàm không hội tụ tại điểm lời giải 
output – xuất ra các thông tin về số vòng lặp tính toán, số lần tính hàm số, các thuật toán 
tại các bước tính, v.v 
x0 – véc tơ tham biến khởi đầu 
Objfun - tên của M-file xác định hàm mục tiêu 
Cú pháp này sử dung pp Nelder-Mead và chỉ tìm được cực tiểu local 
24 4) Bước 4: Tạo M-file Script lời giải, ví dụ Solvefminsearch.m 
clear;clc;format long;warning('off'); 
% Nhap vecto tham bien khoi dau: 
x0 = [1,1]; 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
%Alg='interior-point'; 
%Alg='levenberg-marquardt'; 
Alg='trust-region-dogleg'; 
%Alg='lm-line-search'; 
% Chon Stopping Criteria 
GTN=1e-4;GTL=1e4; 
% Neu Stopping Criteria la sai so cua tham bien 
%StCr='TolX'; 
% Neu Stopping Criteria la sai so cua ham so 
StCr='TolFun'; 
% Neu Stopping Criteria la so luong Iterations 
%StCr='MaxIter'; 
% Neu Stopping Criteria la so luong tinh cac ham so 
%StCr='MaxFunEvals'; 
if strcmp(StCr,'MaxIter')==1 || strcmp(StCr,'MaxFunEvals')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 || strcmp(StCr,'TolFun')==1 
 GT=GTN; 
end 
options = 
optimset('Algorithm',Alg,'Display','iter',StCr,GT,'PlotFcns',@optimplotfval); 
[x,fval,exitflag,output] = fminsearch(@Objfun,x0,options) 
25 
Optimization terminated: 
 the current x satisfies the termination criteria using 
OPTIONS.TolX of 1.000000e-004 
 and F(X) satisfies the convergence criteria using 
OPTIONS.TolFun of 1.000000e-004 
26 
0 5 10 15 20 25 30 35 40 45
30
35
40
45
50
55
60
65
70
75
Iteration
F
u
n
c
ti
o
n
 v
a
lu
e
Current Function Value: 30
27 
QUY HOẠCH TUYẾN TÍNH 
 
 
minT
m n
p n
f
x c x
A x a
B x b
lb x ub
1 1 1 1 1 1
2 2 2 2 2 2
; ; ; ; ; ;
n n n n n n
c x a b l u
c x a b l u
c x a b l u
c x a b lb ub
 
 
11 12 1
21 22 2
1 2
11 12 1
21 22 2
1 2
n
n
m n
m m mn
n
n
p n
p p pn
a a a
a a a
a a a
b b b
b b b
b b b
A
B
28 
[x,fval,exitflag,output,lambda] = linprog(c,A,a,B,b,lb,ub,x0,options) 
x – xuất ra véctơ tham biến x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
• exitflag=1: có nghĩa là hàm đã hội tụ tại điểm lời giải x 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu hoặc số lượng vòng 
lặp đã đạt ngưỡng cho phép nếu Stopping Criteria = 'MaxIter' 
• exitflag=-2: Không tìm được một điểm hợp lệ nào 
• exitflag=-3: Miền ràng buộc là vô hạn 
• exitflag=-4: Giá trị NaN (không phải là số) xuất hiện trong quá trình tính 
toán 
• exitflag=-5: Cả bài toán gốc và đối ngẫu của nó đều vô nghiệm 
• exitflag=-7: Hướng tìm kiếm càng lúc càng nhỏ, không thể tìm được các 
điểm tốt hơn 
output – xuất ra các thông tin về số vòng lặp tính toán, số lần tính hàm 
số, các thuật toán tại các bước tính, v.v 
x0 – véc tơ tham biến khởi đầu 
c,A,a,B,b,lb,ub – lần lượt các véc tơ, ma trận dữ liệu bài toán 
29 
Giải bài toán quy hoạch tuyến tính sau: 
   1
2
8
2 0
; ; 10 ; ; ; ;
3 0
3
x
x
c x a b lb ub
     
3 2
1 2
2 1 ;
0 1 
A B
Đổi lại yêu cầu của mục tiêu: 1 22 3 minf x x x
Nhận diện được các véc tơ và ma trận dữ liệu đề bài: 
1) Bước 1: Phân tích bài toán 
2) Bước 2: Tạo M-file script cho bài toán, ví dụ LP1.m 
30 
clear;clc;format long;warning('off'); 
c=[-2,-3] 
A=[1,2; 
 2,1; 
 0,1] 
a=[8,10,3] 
B=[] 
b=[] 
lb=zeros(2,1) 
ub=[] 
x0=[] 
options = optimset('LargeScale', 'off', 'Simplex', 
'on','Display','iter') 
[x,fval,exitflag,output,lambda] = 
linprog(c,A,a,B,b,lb,ub,x0,options) 
31 
Giải bài toán quy hoạch tuyến tính sau: 
 
 
 
1 1
2 2
minTf 
x c x
A x a
A x a
B x b
lb x ub
1
2
3
1 2
4
5
6
1 2 7
2 5 0 2
3 7 5 10 1 2
; ; ; ; ; ; ;
4 8 1 37 15
5 8 5 4
6 1 10
x
x
x
x
x
x
c x a a b lb ub
     1 2
1 2 3 0 0 0
0 1 2 3 0 0 3 0 0 0 2 1 1 1 1 1 1 1
; ;
0 0 1 2 3 0 0 4 0 2 0 3 5 0 3 0 1 0
0 0 0 1 2 3
A A B
1) Bước 1: Phân tích bài toán 
32 2) Bước 2: Tạo M-file script cho bài toán, ví dụ LP2.m 
clear;clc;format long;warning('off'); 
c=[1,-2,3,-4,5,-6]' 
A1=[1,2,3,0,0,0; 
 0,1,2,3,0,0; 
 0,0,1,2,3,0; 
 0,0,0,1,2,3] 
A2=[3,0,0,0,-2,1; 
 0,4,0,-2,0,3] 
B=[1,1,1,1,1,1; 
 5,0,-3,0,1,0] 
a1=[5,7,8,8]' 
a2=[5,7]' 
b=[10,15]' 
A=[-A1; 
 A2] 
a=[-a1; 
 a2] 
lb=[-2,0,-1,-1,-5,1]' 
ub=[7,2,2,3,4,10]' 
x0=[] 
options = optimset('LargeScale', 'off', 'Simplex', 
'on','Display','iter') 
[x,fval,exitflag,output,lambda] = 
linprog(c,A,a,B,b,lb,ub,x0,options) 
33 
QUY HOẠCH BẬC HAI 
  
 
 
1
min
2
T T
n n
m n
p n
f
x x H x c x
A x a
B x b
lb x ub
 
 
11 12 1
21 22 2
1 2
11 12 1
21 22 2
1 2
n
n
m n
m m mn
n
n
p n
p p pn
a a a
a a a
a a a
b b b
b b b
b b b
A
B
1 1 1 1 1 1
2 2 2 2 2 2
; ; ; ; ; ;
n n n n n n
c x a b l u
c x a b l u
c x a b l u
c x a b lb ub
[H] chính là 
ma trận 
Hessian của 
hàm mục tiêu 
34 
[x,fval,exitflag,output,lambda] = quadprog(H,c,A,a,B,b,lb,ub,x0,options) 
x – xuất ra véctơ tham biến x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
• exitflag=1: có nghĩa là hàm đã hội tụ tại điểm lời giải x 
• exitflag=3: có nghĩa là sự thay đổi giá trị hàm mục tiêu đã nhỏ hơn sai số cho 
phép 
• exitflag=4: có nghĩa là cực tiểu địa phương đã được tìm thấy 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu hoặc số lượng vòng lặp đã 
đạt ngưỡng cho phép nếu Stopping Criteria = 'MaxIter' 
• exitflag=-2: Bài toán không giải được 
• exitflag=-3: Miền ràng buộc là vô hạn 
• exitflag=-4: Hướng tìm kiếm không phải là hướng tốt, không thể tìm được các 
điểm tốt hơn 
• exitflag=-7: Giá trị của hướng tìm kiếm quá nhỏ, không có tiến triển 
• output – xuất ra các thông tin về số vòng lặp tính toán, số lần tính hàm số, các 
thuật toán tại các bước tính, v.v 
x0 – véc tơ tham biến khởi đầu 
H,c,A,a,B,b,lb,ub – lần lượt các véc tơ, ma trận dữ liệu bài toán 
35 
2
1 430
2 2 2
1 2 4 1 2 1 3 2 3 2 4 3 4
0
4 8 72 5 93 6
x xx
x x x x x x x x x x x x x

 
4 4
6
12
0
0
0
7
7
4
9
9
5
8
4
5
4 8
H
Dựa vào mầu để nhận diện 
36 
Giải bài toán quy hoạch bậc 2 sau: 
1) Bước 1: Phân tích bài toán 
 
       
1
2
2 2
8
2 0
; ; 10 ; ; ; ;
3 0
3
1 2
2 0
; 2 1 ;
0 4
0 1
x
x
c x a b lb ub
H A B
2) Bước 2: Tạo M-file script cho bài toán, ví dụ QP1.m 
37 
clear;clc;format long;warning('off'); 
H=[2,0; 
 0,4]; 
c=[2,3]'; 
A=[1,2; 
 2,1; 
 0,1]; 
a=[8,10,3]'; 
B=[]; 
b=[]'; 
lb=[0,0]'; 
ub=[inf,inf]'; 
x0=[]'; 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
%Alg='interior-point'; 
%Alg='levenberg-marquardt'; 
Alg='trust-region-dogleg'; 
%Alg='lm-line-search'; 
% Chon Stopping Criteria 
GTN=1e-4;GTL=1e4; 
% Neu Stopping Criteria la sai so cua tham bien 
%StCr='TolX'; 
% Neu Stopping Criteria la sai so cua ham so 
StCr='TolFun'; 
% Neu Stopping Criteria la so luong Iterations 
%StCr='MaxIter'; 
if strcmp(StCr,'MaxIter')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 || strcmp(StCr,'TolFun')==1 
 GT=GTN; 
end 
options = optimset('Algorithm',Alg,'Display','iter',StCr,GT); 
[x,fval,exitflag,output,lambda]=quadprog(H,c,A,a,B,b,lb,ub,x0,options) 
% Ve do thi 
x1=[-4:0.1:4]; 
x2=[-4:0.1:4]; 
[X1,X2]=meshgrid(x1,x2); 
f=X1.^2 + 2*X2.^2 + 2*X1 + 3*X2; 
meshc(X1,X2,f); 
hold on 
plot(x(1),x(2),'r*') 
38 Giải bài toán quy hoạch bậc 2 sau: 
1) Bước 1: Phân tích bài toán 
 
       
1
2
3
3 3
4 0
6
6 ; ; ; ; 0 ; ;
2
12 0
2 1 0
1 1 1
1 4 2 ; ;
1 1 2
0 2 4
x
x
x
c x a b lb ub
H A B
2) Bước 2: Tạo M-file script cho bài toán, ví dụ QP2.m 
39 
clear; 
clc; 
format long; 
H=[2,1,0; 
 1,4,2; 
 0,2,4]; 
c=[4,6,12]'; 
A=[1,1,1; 
 -1,-1,2]; 
a=[6,2]'; 
B=[]; 
b=[]'; 
lb=[0,0,0]'; 
ub=[inf,inf,inf]'; 
x0=[]'; 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
%Alg='interior-point'; 
%Alg='levenberg-marquardt'; 
Alg='trust-region-dogleg'; 
%Alg='lm-line-search'; 
% Chon Stopping Criteria 
GTN=1e-4;GTL=1e4; 
% Neu Stopping Criteria la sai so cua tham bien 
%StCr='TolX'; 
% Neu Stopping Criteria la sai so cua ham so 
StCr='TolFun'; 
% Neu Stopping Criteria la so luong Iterations 
%StCr='MaxIter'; 
if strcmp(StCr,'MaxIter')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 || strcmp(StCr,'TolFun')==1 
 GT=GTN; 
end 
options = optimset('Algorithm',Alg,'Display','iter',StCr,GT); 
[x,fval,exitflag,output,lambda]=quadprog(H,c,A,a,B,b,lb,ub,x0,options) 
%Neu tim max f: doi cua f se chinh la gia tri cuc dai 
%[x,fval,exitflag,output,lambda]=QUADPROG(H,c,A,a,B,b,lb,ub,x0,options) 
40 
Cực tiểu hóa hàm phi tuyến 
với các ràng buộc 
 
 
 1 2
min
0; 1..
0; 1..
j
k
m n
p n
T
n
f
g j r
h k s
x x x
x
x
x
A x a
B x b
lb x ub
x
1 1 1 1 1
2 2 2 2 2
; ; ; ;
n n n n n
x a b l u
x a b l u
x a b l u
x a b lb ub
 
 
11 12 1
21 22 2
1 2
11 12 1
21 22 2
1 2
n
n
m n
m m mn
n
n
p n
p p pn
a a a
a a a
a a a
b b b
b b b
b b b
A
B
41 
[x,fval,exitflag,output,lambda,grad,hessian]= 
fmincon(@Objective,x0,A,a,B,b,lb,ub,@NonLinCon,options) 
x – xuất ra véctơ tham biến x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
a) Đối với tất cả các thuật toán: 
• exitflag=1: Thước đo tối ưu bậc nhất (First-order optimality measure) nhỏ hơn 
giá trị thông số TolFun và vi phạm ràng buộc lớn nhất (maximum constraint 
violation) nhỏ hơn giá trị thông số TolCon 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu đã đạt ngưỡng cho phép nếu 
Stopping Criteria = FunEvals hoặc số lượng vòng lặp đã đạt ngưỡng cho phép 
nếu Stopping Criteria = 'MaxIter' 
• exitflag=-1: Hàm đầu ra kết thúc thuật toán 
• exitflag=-2: Không tìm được lời giải hợp lệ nào 
b) Đối với 2 thuật toán là ‘trust-region-reflective’ và ‘interior-point’: 
• exitflag=2: Sự biến thiên của x nhỏ hơn giá trị thông số TolX và vi phạm ràng 
buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
c) Đối với riêng thuật toán ‘trust-region-reflective’: 
• exitflag=3: có nghĩa là sự thay đổi giá trị hàm mục tiêu đã nhỏ hơn sai số cho 
phép TolFun và vi phạm ràng buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
42 d) Đối với riêng thuật toán ‘active-set’: 
• exitflag=4: giá trị của hướng tìm kiếm (search direction) đã nhỏ hơn 2 lần giá trị 
cho phép TolX và vi phạm ràng buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
• exitflag=5: giá trị của hướng tìm kiếm (search direction) đã nhỏ hơn 2 lần giá trị 
cho phép TolFun và vi phạm ràng buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
e) Đối với riêng thuật toán ‘interior-point’: 
• exitflag=-3: Điểm x đang xét đi xuống dưới giá trị ObjectiveLimit và vi phạm 
ràng buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
lambda – Giá trị toán tử Lagrange tại điểm lời giải 
grad – xuất ra véctơ Gradient của hàm số tại điểm lời giải 
hessian – xuất ra ma trận Hessian của hàm số tại điểm lời giải 
output – xuất ra các thông tin về số vòng lặp tính toán, số lần tính hàm số, các 
thuật toán tại các bước tính, v.v 
x0 – véc tơ tham biến khởi đầu 
A,a,B,b,lb,ub – lần lượt các véc tơ, ma trận dữ liệu bài toán 
Objective - tên của M-file xác định hàm mục tiêu 
NonLinCon - tên của M-file xác định hàm ràng buộc phi tuyến đẳng thức và bất 
đẳng thức 
43 Giải bài toán tối ưu hóa phi tuyến với các ràng buộc sau 
1 2 1 2 1
1 2
2 2 2
1 2
1 2
3
1 2 2
2
1 2
1
1 2
1 2
1 2
1 2
, 9.82 2 min
2500
500 0
2500
0
0.5882
0.3 48.59 0
10
1.056 0
3 18 12
3 8
10 2.531
2 14;0.2 0.8
f f x x x x x
x x
x x
x x
x x x
x
x x
x
x x
x x
x x
x x
x
 
   
 
2 2
1 2
1
2
3 18
1 3
1 10
12
; ; 2.531 ;
8
2 14
;
0.2 0.8
x
x
A
B
x a b
lb ub
44 1) Bước 1: Tạo 1 thư mục cho bài toán, ví dụ NonlinCon1 
2) Bước 2: Tạo M-file hàm mục tiêu, ví dụ Objective.m 
function f = Objective( x ) 
%OBJECTIVE Summary of this function goes here 
% Detailed explanation goes here 
f= 9.82*x(1)*x(2)+2*x(1); 
end 
3) Bước 3: Tạo M-file hàm các hàm ràng buộc phi tuyến đẳng thức và 
bất đẳng thức, ví dụ NonLinCon.m 
function [ g, h ] = NonLinCon( x ) 
%NONLINCON Summary of this function goes here 
% Detailed explanation goes here 
% Cac rang buoc bat dang thuc 
g = [2500/(pi*x(1)*x(2))- 500; 
 2500/(pi*x(1)*x(2))- (pi^2*(x(1)^2+x(2)^2))/0.5882]; 
% Cac rang buoc dang thuc 
h = [0.3*x(1)^3*x(2)-48.59*x(2); 
 10*x(2)/x(1) - x(1)*x(2) + 1.056]; 
end 
4) Bước 4: Tạo M-file Script lời giải, ví dụ Solve.m 
45 clear;clc;format long;warning off; 
A=[3,-18; 
 1,3]; 
a=[12,8]'; 
B=[1,-10]; 
b=[2.531]'; 
lb=[2,0.2]'; 
ub=[14,0.8]'; 
x0 = [lb(:)+(ub(:)-lb(:))*rand(1,1)]' % Diem khoi dau 
% Kiem tra gia tri cua ham muc tieu cung nhu cac rang buoc phi tuyen tai diem khoi dau 
fprintf ('Gia tri ham muc tieu va rang buoc phi tuyen tai diem khoi dau\n'); 
f=Objective(x0) 
[g, h] = NonLinCon(x0) 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
Alg='interior-point'; 
% Chon Stopping Criteria 
GTN=1e-6;GTL=1e6; 
% Neu Stopping Criteria la sai so cua tham bien 
%StCr='TolX'; 
% Neu Stopping Criteria la sai so cua ham so 
%StCr='TolFun'; 
% Neu Stopping Criteria la so luong Iterations 
StCr='MaxIter'; 
% Neu Stopping Criteria la so luong tinh cac ham so 
%StCr='MaxFunEvals'; 
if strcmp(StCr,'MaxIter')==1 || strcmp(StCr,'MaxFunEvals')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 || strcmp(StCr,'TolFun')==1 
 GT=GTN; 
end 
options = optimset('Algorithm',Alg,'Display','iter-detailed',StCr,GT,'PlotFcns',@optimplotfval); 
[x,fval,exitflag,output,lambda,grad,hessian]=fmincon(@Objective,x0,A,a,B,b,lb,ub,@NonLinCon,options
) 
fprintf('Gia tri rang buoc phi tuyen tai diem cuc tri\n'); 
[g, h] = NonLinCon(x) 
46 CỰC TIỂU HÓA NỬA VÔ HẠN 
(Semi-Infinite Minimization) 
 
 
 
1 2
1 1
2 2
min
0; 1..
0; 1..
, 0
, 0
, 0
j
k
m n
p n
T
n
q q
f
g j r
h k s
x x x
K w
K w
K w
x
x
x
A x a
B x b
lb x ub
x
x
x
x
 
 
11 12 1
21 22 2
1 2
11 12 1
21 22 2
1 2
n
n
m n
m m mn
n
n
p n
p p pn
a a a
a a a
a a a
b b b
b b b
b b b
A
B
1 1 1 1 1
2 2 2 2 2
; ; ; ;
n n n n n
x a b l u
x a b l u
x a b l u
x a b lb ub
47 [x,fval,exitflag,output,lambda]=fseminf(@Objective,x0,q,
myinfcon,A,a,B,b,lb,ub,options) 
x – xuất ra véctơ tham biến x làm cho hàm mục tiêu đạt cực tiểu 
fval – xuất ra giá trị hàm mục tiêu tại điểm cực tiểu x 
exitflag – xuất ra giá trị để xác định điều kiện dừng tính toán, cụ thể là: 
• exitflag=1: Hàm hội tụ tại điểm lời giải x 
• exitflag=4: giá trị của hướng tìm kiếm (search direction) đã nhỏ hơn giá trị cho 
phép và vi phạm ràng buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
• exitflag=5: giá trị của đạo hàm định hướng (directional derivative) đã nhỏ giá trị 
cho phép và vi phạm ràng buộc lớn nhất nhỏ hơn giá trị thông số TolCon 
• exitflag=0: có nghĩa là số lượng tính hàm mục tiêu đã đạt ngưỡng cho phép nếu 
Stopping Criteria = FunEvals hoặc số lượng vòng lặp đã đạt ngưỡng cho phép 
nếu Stopping Criteria = 'MaxIter' 
• exitflag=-1: Hàm đầu ra kết thúc thuật toán 
• exitflag=-2: Không tìm được lời giải hợp lệ nào 
x0 – véc tơ tham biến khởi đầu 
A,a,B,b,lb,ub – lần lượt các véc tơ, ma trận dữ liệu bài toán 
Objective - tên của M-file xác định hàm mục tiêu 
myinfcon - tên của M-file xác định hàm ràng buộc phi tuyến đẳng thức, bất đẳng 
thức và ràng buộc nửa vô hạn 
48 myinfcon
function [g,h,K1,K2,..,Kq,S] = myinfcon(x,S) 
% Initial sampling interval 
if isnan(S(1,1)), 
 S = [s1 0; 
 s2 0; 
 ... 
 sq 0]; % S co q hang and 2 columns 
end 
% Sample set 
w1 = w1_begin:S(1,1):w1_end; % Compute sample set 
w2 = w2_begin:S(2,1):w2_end; % Compute sample set 
... 
wq = wq_begin:S(q,1):wq_end; % Compute sample set 
% Semi-infinite constraints 
K1 = ... % 1st semi-infinite constraint at x and w 
K2 = ... % 2nd semi-infinite constraint at x and w 
... 
Kq = ...% Last semi-infinite constraint at x and w 
% Cac rang buoc bat dang thuc phi tuyen 
g = [ham_g1; 
 ham_g2; 
 ... 
 ham_gr]; % Compute nonlinear inequalities at x 
% Cac rang buoc dang thuc phi tuyen 
h = [ham_h1; 
 ham_h2; 
 ... 
 ham_hs]; % Compute nonlinear equalities at x 
% Plot a graph of semi-infinite constraints 
% plot(w1,K1,'-',w2,K2,':'),title('Semi-infinite constraints') 
% drawnow 
49 Giải bài toán tối ưu hóa nửa vô hạn sau: 
min 
1) Bước 1: Tạo 1 thư mục cho bài toán, ví dụ Seminf1 
2) Bước 2: Tạo M-file hàm mục tiêu, ví dụ Objective.m 
function f = Objective(x,s) 
%OBJECTIVE Summary of this function goes here 
% Detailed explanation goes here 
f=sum((x-0.5).^2); 
%f=(x(1)-0.5)^2 + (x(2)-0.5)^2 + (x(3)-0.5)^2; 
end 
3) Bước 3: Tạo M-file hàm các hàm ràng buộc phi tuyến đẳng thức, 
bất đẳng thức và nửa vô hạn, ví dụ myinfcon.m 
50 function [g,h,K1,K2,s] = myinfcon(x,s) 
% Initial sampling interval 
if isnan(S(1,1)), 
 s = [0.2 0; 
 0.2 0]; % S co q hang and 2 columns 
end 
% Sample set 
w1 = 0:s(1,1):100; % Compute sample set 
w2 = 0:s(2,1):100; % Compute sample set 
% Semi-infinite constraints 
K1 = sin(w1*x(1)).*cos(w1*x(2)) - 1/1000*(w1-50).^2 - sin(w1*x(3))-x(3)-1; % 1st 
semi-infinite constraint at x and w 
K2 = sin(w2*x(2)).*cos(w2*x(1)) - 1/1000*(w2-50).^2 - sin(w2*x(3))-x(3)-1; % 2nd 
semi-infinite constraint at x and w 
% Cac rang buoc bat dang thuc phi tuyen 
g = []; %Trong bai toan nay khong co 1 rang buoc phi tuyen bat dang thuc nao 
% Cac rang buoc dang thuc phi tuyen 
h = []; %Trong bai toan nay khong co 1 rang buoc phi tuyen dang thuc nao 
% Plot a graph of semi-infinite constraints 
% plot(w1,K1,'-',w2,K2,':'),title('Semi-infinite constraints') 
% drawnow 
4) Bước 4: Tạo M-file Script lời giải, ví dụ Solve.m 
51 clear;clc;format long;warning off; 
A=[]; %Trong bai tap nay khong co rang buoc tuyen tinh bat dang thuc 
a=[]'; %Trong bai tap nay khong co rang buoc tuyen tinh bat dang thuc 
B=[]; %Trong bai tap nay khong co rang buoc tuyen tinh dang thuc 
b=[]'; %Trong bai tap nay khong co rang buoc tuyen tinh dang thuc 
lb=[]';%Trong bai tap nay khong co rang buoc phia truoc cua tham bien 
ub=[]';%Trong bai tap nay khong co rang buoc phia sau cua tham bien 
x0 = [0.5;0.2;0.3]; % Diem khoi dau 
% Kiem tra gia tri cua ham muc tieu cung nhu cac rang buoc phi tuyen tai diem khoi dau 
fprintf ('Gia tri ham muc tieu va rang buoc phi tuyen tai diem khoi dau\n'); 
f=Objective(x0,NaN) 
[g,h] = myinfcon(x0,NaN) 
% Chon Algorithm (Chon mot trong so duoi day) 
%Alg='active-set'; 
%Alg='trust-region-reflective'; 
%Alg='interior-point'; 
%Alg='levenberg-marquardt'; 
%Alg='trust-region-dogleg'; 
Alg='lm-line-search'; 
% Chon Stopping Criteria 
GTN=1e-6;GTL=1e6; 
% Neu Stopping Criteria la sai so cua tham bien 
%StCr='TolX'; 
% Neu Stopping Criteria la sai so cua ham so 
%StCr='TolFun'; 
% Neu Stopping Criteria la so luong Iterations 
StCr='MaxIter'; 
% Neu Stopping Criteria la so luong tinh cac ham so 
%StCr='MaxFunEvals'; 
if strcmp(StCr,'MaxIter')==1 || strcmp(StCr,'MaxFunEvals')==1 
 GT=GTL; 
elseif strcmp(StCr,'TolX')==1 || strcmp(StCr,'TolFun')==1 
 GT=GTN; 
end 
options = optimset('Algorithm',Alg,'Display','iter-detailed',StCr,GT,'PlotFcns',@optimplotfval); 
[x,fval,exitflag,output,lambda]=fseminf(@Objective,x0,2,@myinfcon,A,a,B,b,lb,ub,options) 
fprintf('Gia tri rang buoc phi tuyen tai diem cuc tri\n'); 
[g, h] = myinfcon(x,NaN) 

File đính kèm:

  • pdfbai_giang_toi_uu_hoa_trong_thiet_ke_co_khi_chuong_11_ung_dun.pdf