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
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="">0:>
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
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
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:
bai_giang_toi_uu_hoa_trong_thiet_ke_co_khi_chuong_11_ung_dun.pdf

