Điều khiển dự báo dựa trên mô hình cho hệ phi tuyến với tầm dự báo bằng 1

Lê Thị Huyền Linh và Đtg

Tạp chí KHOA HỌC & CÔNG NGHỆ

112(12)/2: 55 - 61

ĐIỀU KHIỂN DỰ BÁO DỰA TRÊN MÔ HÌNH CHO
HỆ PHI TUYẾN VỚI TẦM DỰ BÁO BẰNG 1
Lê Thị Huyền Linh1*, Lại Khắc Lãi2,
Nguyễn Thị Mai Hương1
1

Trường Đại học Kỹ thuật Công nghiệp – ĐH Thái Nguyên
2
Đại học Thái Nguyên

TÓM TẮT
Điều khiển dự báo theo mô hình (MPC Model Predictive Control) đã được tác giả đề cập và
nghiên cứu trong một số công trình [1], [2]. Phương pháp điều khiển dự báo đã cải thiện chất
lượng điều khiển một cách đáng kể so với các phương pháp khác. MPC đã được nghiên cứu và
ứng dụng rộng rãi trong công nghiệp, đặc biệt là đối với hệ tuyến tính biến đổi chậm[3], [4]. Bài
báo này sẽ đi sâu phân tích xây dựng mô hình dự báo cho hệ phi tuyến với thuật toán xác định
phiếm hàm mục tiêu và xây dựng phương pháp điều khiển dự báo cho hệ phi tuyến với tầm dự báo
bằng 1.
Từ khoá: Điều khiển dự báo, mô hình dự báo cho hệ phi tuyến với tầm dự báo bằng 1.

GIỚI THIỆU CHUNG*
Điều khiển dự báo dựa trên mô hình là sự kết
hợp của một số lĩnh vực đã được phát triển
trong lý thuyết điều khiển hiện đại, điển hình
đó là hai lĩnh vực điều khiển tối ưu và nhận
dạng hệ thống. Ngay như tên của nó “điều
khiển dự báo dựa trên mô hình” có nghĩa là
trong đó cần phải sử dụng một mô hình dự
báo để ước lượng (dự báo) các giá trị của đầu
ra trong tương lai để phục vụ cho bài toán
điều khiển. Điều khiển dự báo dựa trên mô
hình có thể kết hợp chặt chẽ hay đưa được các
điều kiện ràng buộc về mặt vật lý của quá
trình (như độ mở van, các hạn chế của cơ cấu
chấp hành, các giới hạn của tín hiệu điều
khiển v.v) trong thiết kế bộ điều khiển và
chuyển hóa bài toán thiết kế bộ điều khiển
thành một bài toán tối ưu. Hiện nay MPC đã
trở thành một sách lược điều khiển cao cấp
được chấp nhận khá rộng rãi trong một số lĩnh
vực công nghiệp. Đã có hơn 3000 ứng dụng
của MPC đã được thương mại hóa trong các
lĩnh vực khác nhau bao gồm: công nghệ lọc
hóa dầu, công nghệ xử lý thực phẩm, công
nghệ ô tô, công nghệ không gian, công nghệ
bột giấy và giấy v.v [4].
*

Tel: 0918127781; Email: lethihuyenlinh@gmail.com

XÂY DỰNG MÔ HÌNH DỰ BÁO CHO HỆ
PHI TUYẾN
Giả thiết cho hệ phi tuyến được biểu diễn tại
t ,k = 0,1,...,∞

thời điểm rời rạc k
như sau:
x = f (x ) + h (x )u + υ
y = Cx

với các ràng buộc như sau:
k +1

k

k

u

y

≤ u

m in

∆u

x

k

u ∈ℝ
hình, k

≤

y

∈ ℝ

m

k

f (x ), h (x )
k
k

chiều phù hợp,

u

≤∆u
≤

k

k

m ax

m ax

y

(1)

,

,

m ax

(2)

n

là biến trạng thái của mô

υ ∈ℝ
là véc tơ đầu vào, k

y ∈ℝ
là véc tơ nhiễu, k

hệ.

k

≤

k

≤∆u

m in

m in

trong đó:

k

k

n

p

là véc tơ đầu ra của

là các hàm phi tuyến với số

C ∈ ℝ

p ×n

,

u
≤u ≤u
, ∆u
≤ ∆u ≤ ∆u
,y
≤y ≤y
min
k
max min
k
max min
k
max

là các véc tơ chặn dưới và chặn trên.
xˆ
Giả thiết rằng k + j |k là véc tơ giá trị dự báo

của biến trạng thái tại thời điểm

t
k+j

được
55

Lê Thị Huyền Linh và Đtg

Tạp chí KHOA HỌC & CÔNG NGHỆ

t , ∆u = u − u
k
k
k −1
ước lượng tại thời điểm k

∆uˆ
= uˆ
− uˆ
k + j |k
k + j |k
k + j −1|k là giá trị dự
và
báo lượng số gia đầu vào (lượng gia tăng tín
u

hiệu điều khiển) tại thời điểm k + j được ước
lượng (dự báo) tại thời điểm k , vậy hàm mục
tiêu có thể được viết như sau:

 p −1 

J =F  xˆ
,∆u
+ ∑ G  xˆ
k
k
+
p
|
k
k + j |k 

 j =0  k + j |k

(3)

F (.),G (.)

với các hàm
được gọi là các hàm
phạt và hàm giá tại trạng thái cuối (kết thúc),
p được gọi là tầm dự báo.

J
Hàm mục tiêu

k thường có dạng toàn

phương, giả thiết rằng

w
k + j |k

là giá trị đặt của

x
k +j

tại thời điểm k . Gọi ma trận xác định
bán dương Q và ma trận xác định dương R là
các ma trận trọng số thì hàm mục tiêu (3) có
thể viết lại dưới dạng như sau:
2
2
p
p −1
J = ∑ xˆ
−w
+ ∑ ∆u
k
k + j |k
k + j |k
k + j |k R
j =1
Q j =0

(4)
Vậy từ công thức (1) và (4), bài toán điều
khiển dự báo cho hệ phi tuyến tại từng thời
điểm trích mẫu trở thành bài toán cực tiểu hóa
J

hàm k và xác định lượng số gia tương ứng
cho
tín
hiệu
điều
khiển


∆u
∆u
∆u
...∆u
k +1|k
k + 2|k
k + p −1|k 
 k |k

với các ràng
buộc (2). Để đơn giản, từ công thức (4), ta
xˆ

thấy k + j |k có thể xác định được thông qua dự
báo đầu ra nếu ma trận C là ma trận tuyến
tính hằng, như vậy hàm mục tiêu (4) cho bài
toán điều khiển dự báo có thể chuyển thành:
Jk =
=

p

∑

j =1

p

∑

j =1

2

C xˆk + j |k − w k + j |k
2

yˆk + j |k − w k + j |k

+
Q

+
Q
p −1

∑

j =0

p −1

∑

j =0

∆ u k + j |k

∆ u k + j |k

2

u
k + j |k

thỏa mãn yêu cầu thực tế của các bài
toán điều khiển.
ĐIỀU KHIỂN DỰ BÁO CHO MÔ HÌNH HỆ
PHI TUYẾN VỚI TẦM DỰ BÁO BẰNG 1
Ngoại trừ các hệ phi tuyến đặc biệt như mô
hình phi tuyến Hammerstein thì điều khiển dự
báo dựa trên mô hình cho hệ phi tuyến với
tầm dự báo lớn hơn 1 là rất khó khăn và hầu
hết không thỏa mãn điều kiện về lời giải, đặc
biệt là trong trường hợp có kể đến nhiễu. Tuy
nhiên, điều này có thể khắc phục được nếu hệ
là khả nghịch, khi đó điều khiển dự báo cho
hệ phi tuyến luôn luôn có nghiệm [5]. Vì vậy
bài báo này chỉ giới hạn ở việc nghiên cứu
điều khiển dự báo mô hình hệ phi tuyến với
tầm dự báo bằng 1.
Quay trở lại với hệ (1), điều khiển dự báo với
tầm dự báo bằng 1 có thể được suy ra trực
∆u

=u

 
 
xˆ
= f  x  + g x  u
k +1|k
 k
 k  k |k
 
 
 
= f x  + g  x  u
+ g  x  ∆u
k |k
 k
 k  k −1
 k
1
= x
+ g (x )∆u
k +1|k
k
k |k

(6)

x1

Trong biểu thức trên, k +1|k là ký hiệu của
phần trong đó chứa các dữ liệu đã biết

(xk ,uk −1 )

tại thời điểm k , và

( k )∆uk|k

g x

là

xˆ
.
k +1|k

phần chưa biết của trạng thái dự báo
Nếu mô hình không chính xác do ảnh hưởng
của nhiễu và các sai lệch của mô hình thì sai
lệch dự báo của (6) sẽ có dạng như sau:
xɶ
=x
−xˆ
=ζ
k +1|k k +1 k +1|k k +1

Trong đó

ζ

k +1

(7)

là nhiễu do sai lệch mô hình và
ζ

(5)

không và có phương sai

∆u
k + j |k

Tại các thời điểm
trong phiếm hàm
mục tiêu (5) có thể được thay đổi sao cho
56

−u

tiếp k |k k|k k −1 với một dữ liệu chưa biết
tại thời điểm k như sau:

nhiễu của hệ (1), nếu k là nhiễu ngẫu nhiên
ngẫu nhiên dừng, có kỳ vọng toán bằng

R

2
R

112(12)/2: 55 - 61

thể thấy rằng

E ζ  =δ 2 ,
 k



E xɶ
 = 0
k
+
|
k
1



và

thì ta có

(5)

Lê Thị Huyền Linh và Đtg

Tạp chí KHOA HỌC & CÔNG NGHỆ



T


 

 
2
ɶ
ɶ
E  xɶ
−E xɶ
x
−
E
x



 k +1|k   =nδ
 k +1|k
 k +1|k    k +1|k




hay nói cách khác cả kỳ vọng toán và phương
sai của sai lệch dự báo là cực tiểu, do vậy bài
toán dự báo (6) là bài toán dự báo tối ưu.
x

Nếu giá trị đặt là sp , đủ trơn và là đường
cong mong muốn của trạng thái trong tương
lai thì giá trị trạng thái mong muốn tại thời
w

điểm k +1 được chọn là k +1|k

=α x +(1−α )x
k
sp

α ∈0,1

  được gọi là hệ số trơn hóa,
trong đó
do vậy hàm mục tiêu của bài toán điều khiển
dự báo cho hệ phi tuyến có nhiễu với tầm dự
báo bằng 1 có thể được viết như sau:

2

2

J = xˆ
−w
+ ∆u
k
k + j |k k + j |k Q
k + j |k R

(8)
Để cực tiểu hóa (8) với các ràng buộc, ta phải
∂ 2J

∂J

có

k =0
∂∆u
k |k

và

k >0
∂∆u 2
k |k

Trước hết, giả sử rằng các ràng buộc có thể
aT ∆u ≤b ,i =1, 2,...,q
i
k |k i

viết lại được dưới dạng
khi đó các ràng buộc có thể được biểu diễn
dưới dạng ma trận như sau:
A∆u ≤B
k |k
(11)
Trong đó:

A= aT
1

aT
2


T 
 
L  λ  = J + λ  aT ∆u −b  , i = 1, 2, ..., q ,
k  i
k
i  i
k |k i 

đặt

∂L
= H ∆u
+F +a λ = 0
k |k
i i
∂∆u
k |k

∂L
T
= a ∆u
− b = 0,
i
k |k
i
∂∆λ
i
và

∆u

k |k

= −H

−1 

 F +a λ 

i i

(12)

aT H −1F +b
i
λ = − i
i
aT H −1a
i
i

do vậy:
(9)

Ký hiệu
H = g ( x )T Qg ( x )+R
k
k



F =g ( x )Q  xˆ
−w
k  k +1|k k +1|k 

và
thì lượng gia tăng tín hiệu điều khiển tại bước
k là:
∆u =− H −1F
k |k

(10)
Nhưng trong thực tế điều khiển, các tín hiệu
đầu vào và đầu ra luôn có những giới hạn của
nó, do vậy có những kết quả thu được từ biểu
thức (10) trên thường không phù hợp. Để thỏa
mãn các ràng buộc, chúng ta phải đưa các
giới hạn mang tính logic vào trong các giá trị
u ,x
k k . Để đơn giản chúng ta sẽ sử dụng

phương pháp Lagrange [5] .

T
T

⋯ aT  ,B = a a ⋯ a 
q 
q
1 2
.

Chọn hàm Lagrange là

thì

−1

 


∆u =− g (x )T Qg (x )+R  g (x )Q xˆ
−w

k|k  k
k
  k  k +1|k k +1|k  

112(12)/2: 55 - 61

(13)

λ

Nếu i thu được từ (13) nhỏ hơn hoặc bằng
0, có nghĩa là điều kiện ràng buộc tương ứng
sẽ không ảnh hưởng đến
ta có thể chọn

λ = 0,
i

∆u
k |k

, do vậy chúng

nhưng nếu

λ >0
i

ràng buộc tương ứng sẽ ảnh hưởng đến

, các
∆u
k

λ =λ ,

do đó chúng ta phải chọn i i cuối cùng
lượng số gia tín hiệu điều khiển cho bài toán
điều khiển dự báo cho hệ phi tuyến với tầm
dự báo bằng 1 là:


∆u =−H −1 F +AT Λ 
k |k



trong đó

(14)

T
Λ= λ λ ⋯ λ 
1
2
q



Tóm lại ta có thuật toán tại từng thời điểm
trích mẫu để xác định lượng số gia tín hiệu
điều khiển như sau:
57

Lê Thị Huyền Linh và Đtg

Tạp chí KHOA HỌC & CÔNG NGHỆ

Bước 1: Xác định ma trận Q và R từ hàm
mục tiêu của bài toán, xác định hệ số trơn
hóa α .
ω
k +1|k
Bước 2: Xác định tín hiệu đặt
cho
bước tiếp theo
Bước 3: Xác định các ma trận H và F
Bước 4: Xác định lượng gia tăng tín hiệu
điều khiển theo công thức (12)
Bước 5: Tính tín hiệu điều khiển tại thời
điểm k theo công thức

+∆u
u =u
k |k k −1
k |k
y ,

Bước 6: Đo đầu ra của hệ k xác định
y

x

k

= Cx

k
bằng công thức k
Bước 7: Xác định tín hiệu trạng thái của hệ
tại thời điểm k+1 tiếp theo bằng mô hình dự
báo (6)
Bước 8: Thời điểm k+1 tiếp theo, quay lại
bước 1
MINH HỌA BẰNG VÍ DỤ MÔ PHỎNG
Để kiểm chứng các kết quả trên, ta áp dụng
thuật toán đã trình bày ở phần 3. Cho hệ thống
bồn nước [8] như mô tả trên hình vẽ sau:

γ

biến trạng thái k . Trong bài toán điều khiển
bồn nước, chúng ta chọn trạng thái hệ thống
y = x
k hay hàm hệ
là đầu ra, có nghĩa là k

thống

x

k +1

= x − 0.2021 x + 0.01923u + γ
k
k
k
k
x ∈0%,100%
k
u ∈ 0%,100%
k

Trong đó

xk

và

Để thay đổi chiều cao của mức nước trong
bồn, chúng ta thay đổi lưu lượng dòng nước
đổ vào bồn bằng cách điều chỉnh van V1 và
quan hệ thông thường giữa độ mở van và lưu
lượng dòng chảy đổ vào bồn được cho trong
hình sau:

Hình 2. Quan hệ giữa độ mở van và lưu lượng
dòng chảy vào trong bồn

Sau đây ta sẽ áp dụng thuật toán điều khiển
dự báo với tầm dự báo bằng 1 cho bài toán ổn
định mức trong bồn với giả thiết van V2 luôn
được mở, đây là điểm đặc trưng cho nhiễu tác
động lên hệ trong mô hình (1).
Trước hết chọn hàm mục tiêu là:
k

x
sp

Giả thiết các biến trong hệ thống ở điều kiện
thông thường và thời gian trích mẫu là 1s, ta
có mô hình dạng (1) như sau:

là

 
f  x  = x − 0.2021 x
k
 k
k

g (x ) = 0.01923
k

J
Hình 1. Mô hình hệ thống bồn nước

112(12)/2: 55 - 61

2
2
= (xˆ
−w
) + 0.001∆u
k +1|k
k +1|k
k |k
= 30%

(16)
và chọn hệ số trơn hóa α = 0.975 .
Giả thiết rằng mô hình không có sai lệch, tức
υ

= 0.

là trong mô hình (1) có k
Kết quả mô phỏng được thể hiện như trên các
hình vẽ sau rõ ràng đã đáp ứng yêu cầu về
mục tiêu điều khiển.

(15)

là độ cao của mức nước trong

u

bồn, k lưu lượng dòng nước đổ vào trong
bồn từ bơm P1 và van V1, trong khi van V2
luôn được mở ở một góc bất kỳ nào đó, đặc
trưng cho nhiễu tác động lên hệ thông qua
58

Hình 3. Đáp ứng mức nước x và lưu lượng dòng chảy
u (không có sai lệch mô hình)

Lê Thị Huyền Linh và Đtg

Tạp chí KHOA HỌC & CÔNG NGHỆ

Để nghiên cứu ảnh hưởng của sai lệch mô
hình, ta sẽ sử dụng mô hình của bồn nước có
sự sai lệch so với mô hình (1) như sau:
x

k +1

= x

k

− 110% × 0.2021 x

+90% × 0.01923u

+ γ

k

k

x

k +1

= x

k

− 0.2021 x

k

e = x −x
s
sp giữa mô
Bảng 1 so sánh sai lệch
phỏng và phân tích lý thuyết với các giả
thuyết:
+ Mô hình mô phỏng:

x

k

(17)
trong khi đó ta vẫn sử dụng mô hình dự báo là:
+ 0.01923u

(18)
Sau khi thực hiện mô phỏng ta nhận được kết
quả như sau:

k +1

=x −110%×0.2021 x + 90%×0.01923u
k
k
k và

+ Mô hình dự báo:
x

k

112(12)/2: 55 - 61

k +1

=x − 0.2021 x + 0.01923u
k
k
k

Từ bảng 1 chúng ta thấy rằng, chúng ta
không thể bỏ qua sai lệch tĩnh bằng cách điều
chỉnh α , do vậy chúng ta có thể sử dụng
mạch bù phản hồi để tạo ra đại lượng bù bù
lại sai lệch này, giải thích sai lệch dự báo tại
e
thời điểm k là k như sau:



=x − xˆ
+g (x
e =x −xˆ
)∆
k k k|k −1 k  k|k −1
k −1 k −1 

Trong đó:
Hình 4. Đáp ứng mức nước x và lưu lượng dòng
chảy u (có sai lệch mô hình)

Ta nhận thấy rằng kết quả điều khiển không
đáp ứng được yêu cầu vì đã tồn tại sai lệch
điều khiển, điều này chính là do ảnh hưởng
của sai lệch mô hình. Để giải quyết vấn đề
này, ta sử dụng phương pháp bù sai lệch sử
dụng phản hồi.
Khi mô hình có sai lệch thì khi đó có sai lệch
tĩnh, sai lệch đó không phụ thuộc vào ma trận
Q mà phụ thuộc và hệ số trơn hóa α .
e = x −x
s
sp giữa mô phỏng
Bảng 1: So sánh
và phân tích lý thuyết
α

Q

e = x −x
s
sp

e = x −x
s
sp

Mô phỏng %

Giá trị của
15%

0.975 0
0.001
0.01

-8.3489
-8.3489
-8.3489

-8.3489
-8.3489
-8.3489

0.95

-4.5279
-4.5279
-4.5279

-4.5279
-4.5279
-4.5279

0
0.001
0.01

x

k

(19)

nhận được bởi phản hồi hệ

xˆ
thống tại thời điểm k và k |k −1 là giá trị dự
x
báo của k tại thời điểm k-1.
e

x

Sau đó thêm k vào giá trị dự báo của k +1
tại thời điểm k một cách trực tiếp, qua đó (6)
có thể được viết lại như sau:
xˆ
= f (x ) + g(x )u
+ g(x )∆u + e
k+1|k
k
k k −1
k k|k
k
= xˆ
+ g(x )∆u + e
k +1|k
k k|k
k

(20)
Sử dụng giá trị dự báo mới để tiến hành thuật
toán NMPC, kết quả mô phỏng như ở hình 5.
Có thể thấy rằng chúng ta đã thay đổi được
tính bền vững của nó dưới sai lệch mô hình,
với phương pháp bù phản hồi này chúng ta đã
hoàn toàn loại bỏ được sai lệch tĩnh.

Hình 5. Đáp ứng mức nước x và lưu lượng dòng chảy
u (có sai lệch mô hình) và có bù phản hồi trực tiếp

59