Ứng dụng phương pháp Hiber-Hughes-Taylor-α phân tích kết cấu tháp điều áp chịu tải trọng động

KHOA HỌC CÔNG NGHỆ 
 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 1
ỨNG DỤNG PHƯƠNG PHÁP HIBER-HUGHES-TAYLOR-α 
PHÂN TÍCH KẾT CẤU THÁP ĐIỀU ÁP CHỊU TẢI TRỌNG ĐỘNG 
Nguyễn Quang Hùng 
Trường Đại học Thủy lợi 
Tóm tắt: Tháp điều áp là một trong những hạng mục quan trọng trong hệ thống nhà máy thủy 
điện và nhất là đối với công trình thủy điện nhỏ chủ yếu dùng hình thức khai thác kiểu đường 
dẫn. Kết cấu tháp điều áp đặt ngầm và có mực nước trong tháp thường xuyên thay đổi nhằm 
đảm bảo điều kiện hoạt động ổn định của tuyến áp lực nên điều kiện chịu lực tương đối phức 
tạp. Nội dung nghiên cứu chủ yếu của bài báo này ứng dụng phương pháp HHt-α tiến hành 
nghiên cứu tác động của dao dộng mực nước tới tháp điều áp. Thông qua đó xác định trường 
ứng suất trong kết cấu tháp điều áp thay đổi theo thời gian ứng với diễn biến mực nước thay đổi 
trong tháp theo điều kiện thực tế nhằm có những ứng xử về công trình thích hợp đảm bảo điều 
kiện kinh tế kỹ thuật. 
Từ khóa: Tháp điều áp, Chu kỳ, Ứng suất, Phần tử hữu hạn, HHT-α 
Summary: Pressurized tower is one of the important items in the hydropower plant system and 
especially for small hydropower projects, mainly using the type of path exploitation. The 
structure of the imposing control tower and the water level in the tower is constantly changing to 
ensure stable operating conditions of the pressure line, so the bearing conditions are relatively 
complicated. The main research contents of this paper apply HHt-α method to study the impact 
of water level fluctuation on the pressure regulator tower. Thereby determining the stress field in 
the structure of the pressurized tower structure changes over time in response to changes in the 
water level in the tower according to the actual conditions in order to have appropriate 
construction behavior to ensure economic conditions and Technology. 
Keywords: Pressure regulator, Cycle, Stress, Finite element, HHT-α 
1. ĐẶT VẤN ĐỀ * 
Việt Nam là một quốc gia có tiềm năng thủy 
điện khá lớn, tập trung ở vùng Tây Bắc, miền 
Trung và Tây Nguyên. Việc khai thác nguồn 
tài nguyên thiên nhiên quý giá này sẽ góp phần 
đảm bảo an ninh, an toàn năng lượng quốc gia, 
nhất là trong bối cảnh các nguồn năng lượng 
hóa thạch, năng lượng không tái tạo ngày càng 
suy giảm, cạn kiệt, qua đó cũng góp phần 
chuyển dịch cơ cấu kinh tế theo hướng tăng tỷ 
trọng công nghiệp của tỉnh trong khu vực này 
Ngày nhận bài: 03/4/2019 
Ngày thông qua phản biện: 22/4/2019 
Ngày duyệt đăng: 26/4/2019 
[1]. Các hình thức khai thác thủy năng có thể 
được kể đến như: Xây dựng đập tạo cột nước, 
tập trung cột nước bằng đường dẫn và xây 
dựng đập kết hợp đường dẫn. Đối với loại hình 
tập trung cột nước bằng đường dẫn: cột nước 
chủ yếu do đường dẫn tạo nên. Với loại hình 
này trên tuyến áp lực phải bố trí tháp điều áp 
nhằm điều hòa áp lực nhằm giảm áp lực nước 
va khi tăng, giảm phụ tải. Hạng mục công 
trình này giúp cho tua bin làm việc ổn định, an 
toàn Sơ đồ bố trí điển hình được thể hiện trên 
hình 1 và quá trình dao động của mực nước 
trong tháp được thể hiện trên hình 2. 
KHOA HỌC CÔNG NGHỆ 
TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 2
Hình 1. Sơ đồ đặt tháp điều áp 
1- Tháp điều áp phía thượng lưu; 2- Tháp điều 
áp phía hạ lưu; 3- Nhà máy thuỷ điện; 4- 
Đường hầm dẫn nước; 5- Đường ống áp lực 
dẫn nước vào turbin 
Hình 2. Sơ đồ dao động mực nước trong 
tháp điều áp 
Hình 2 thể hiện quá trình giảm tải đột ngột 
turbin từ Q0 xuống Q1. Do quán tính của dòng 
chảy, lưu lượng vào đường hầm dẫn nước vẫn 
là Q0, như vậy sẽ có một trị số lưu lượng Q = 
Q0 - Q1 chảy vào tháp, làm cho mực nước 
trong tháp dâng lên dần, từ đó độ chênh lệch 
mực nước giữa thượng lưu (trong hồ chứa) và 
trong tháp giảm dần, dẫn đến vận tốc dòng 
chảy giảm dần, do đó lưu lượng trong đường 
hầm giảm dần. Nhưng cũng do quán tính của 
dòng chảy, mực nước trong tháp không dừng ở 
mực nước tương ứng với lưu lượng Q1 trong 
đường hầm mà vẫn tiếp tục dâng lên thậm chí 
cao hơn cả mực nước thượng lưu. Sau đó, để 
cân bằng thuỷ lực nước phải chảy ngược trở 
lại về thượng lưu, mực nước trong tháp hạ 
xuống. Nhưng cũng do lực quán tính nó lại hạ 
xuống quá mức nước cân bằng và dòng chảy 
lại chảy vào tháp. Cứ như vậy, mực nước trong 
tháp dao động theo chu kỳ và tắt dần do ma 
sát. Cuối cùng mực nước trong tháp dừng ở 
mực nước ổn định mới ứng với lưu lượng Q1. 
Chính sự dao động này dẫn tới phát sinh tải 
trọng động tác động lên kết cấu tháp điều 
áp.[2,3,4] 
2. PHƯƠNG PHÁP NGHIÊN CỨU 
2.1 Phương trình động lực cơ bản 
Phương trình cơ bản của bài toán phần tử hữu 
hạn dưới tác dụng của tải trọng động [5] viết 
dưới dạng: 
        
  
nn
nddndndd
FP
uKuCuM
 
(1) 
Trong đó: 
[M] and [K] lần lượt là ma trận khối lượng 
và ma trận cứng phần tử; [C] là ma trận cản 
Rayleigh là hàm số tuyến tính của ma trận 
khối lượng [M] và ma trận cứng [K], {∆P} là 
véc tơ gia tăng nội lực, {∆F} là véc tơ gia 
tăng ngoại lực, {∆u} là véc tơ gia tăng 
chuyển vịare the incremental displacement 
vectors; Chỉ số n biểu thị bước thời gian. 
Phương pháp Newmark β sử dụng khai triển 
taylor thu được [6]: 
    
     tuuu
t
u
uu
t
u
t
u
nnnn
nnnn









2
2
2
111
2
(2) 
Trong đó: 
∆t là số gia thời gian, và β là các tham số của 
phương pháp. Phương trình (1) có thể được 
giải bằng phương pháp Newmark, phương 
trình (1) được viết lại thành 
KHOA HỌC CÔNG NGHỆ 
 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 3
      
       
    ndd
nddnf
nddd
uCM
t
utCMFD
uKC
t
M
t











1
2
2
2
1
1
(3) 
Trong đó: 
∆t là số gia thời gian, δ và β là các tham số của 
phương pháp và chỉ số t biểu thị một bước 
thời gian. 
2.2. Phương pháp giải. 
Khi giải phương trình động lực (1), thường 
dùng phương pháp Newmark-β với khai triển 
taylor được quy định bởi công thức (2). Tuy 
nhiên khi giải phương trình (3) phương pháp 
này vẫn còn hạn chế do không thể đồng thời 
đạt được cả độ chính xác về thành phần cao 
tần của tổn hao năng lượng và độ chính xác 
của nghiệm giải phương trình. Hiber và đồng 
nghiệp (1977) đã đưa ra phương pháp Hilber-
Hughes-Taylor-α [7,8] để giải quyết vấn đề 
này. Phương pháp này dùng 3 tham số α、β、
γ, với điều kiện 0
3
1
 , ）（  21
2
1
 ,
21
4
1
）（  , đối với bài toán tuyến tính, 
không những ổn định vô điều kiện mà đối với 
thành phần cao tần của lực cản còn đạt được 
độ chính xác bậc 2, đối với thành phần thấp 
tần thì những hiên tượng dị thường đương đối 
ít hơn so với các phương pháp khác. 
Dựa vào phương pháp Hilber-Hughes-
Taylor-α (HHT-α)，phương trình (1) được 
biến đổi thành: 
 {F}-){F}+(1={P}n -
){P}++(1[Kp]{u}-}u)[Kp]{+(1
 + }[C]{-}u)[C]{-(1+[M]{u}
n1+n
1+nn1+n
n1+n1+n

 u
(4) 
Đối với mỗi bước thời gian đều tiến hành tính 
toán theo các bước dưới đây để tìm ra được số 
gia gia tốc 1+i}u{  . Từ đó có thể tiến hành tính 
toán được số gia chuyển vị, biến dạng, ứng 
suất cũng như áp lực lỗ rỗng tại bước thời gian 
tính toán: 
 (i)=}u{[M] 1+iieff  (5) 
Trong đó i là bước tính lặp thứ i trong bước 
tính toán n. 
n
i
nnp
i
np
i
n
i
nnn
i
PPuKuK
uCuMFF
}{}){1(}]{[}]{)[1(
}]{)[1(}]{[}{}){1(
)(
1
)(
1
)(
1
)(
11
)(
 
  (6) 
)(2
2)(
][)1(
][)1(][)1(][][
i
T
p
i
eff
Kt
KtCtMM
 
  (7) 
Trong đó: 
[KT] Ma trận tiếp tuyến cứng tổng thể. 
Từ phương trình (5) giải ra đượcsố gia gia tốc 
{Δa}i+1. Tiến hành điều chỉnh lại chuyển vị, 
vận tốc và gia tốc tại bước lặp tiếp theo: 
2)1(
1
1
1
)1(
1
1
1
)1()(
1
1
1
}{}{}{
}{}{}{
}{}{}{
tuuu
tuuu
uuu
i
n
i
n
i
n
i
n
ii
n
i
n




 
(8) 
Dùng kết quả 11
1
1
1
1 }{;}{;}{
i
n
i
n
i
n uuu  tính tóan 
được từ phương trình (8) thay vào phương 
trình (6) tính được giá trị ψ(i+1) tại bước 
tiếp theo: 
n
i
nnp
i
np
i
n
i
nnn
i
PPuKuK
uCaMFF
}{}){1(}]{[}]{)[1(
}]{)[1(}]{[}{}){1(
)1(
1
)1(
1
)1(
1
)1(
11
)1(
 
  (9) 
Điều kiện hội tụ: 
2)1(
1
)1(
1)0(
)1(
}{
}{
; 


i
n
ii
u
u


 (10) 
Tiến hành lặp từ phương trình (6)-(9) đến 
khi nào thỏa mãn điều kiên hội tụ (10) đồng 
nghĩa với việc giải xong bài toán tại bước 
thời gian n. 
3. KẾT QUẢ VÀ THẢO LUẬN 
Nghiên cứu được tiến hành với một tháp điều 
áp cụ thể có các đặc trưng như sau: 
KHOA HỌC CÔNG NGHỆ 
TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 4
Hình 3. Cắt dọc tuyến áp lực 
Hình 4. Diễn biến mực nước dao động 
trong tháp điều áp 
F11 F22 
M11 M22 
Hình 5. Kết quả nội lực trường hợp 
mực nước min 105.95 m 
F11 F22 
M11 M22 
Hình 6. Kết quả nội lực trường hợp 
mực nước min 121.11 m 
Bảng 1: Kết quả nội lực thp điều áp ứng với mực nước tĩnh 
Nội lực 
MN=121.11 MN=105.95 
Min Max Min Max 
F11(T/m) -36.22 116.46 -48.48 42.25 
F22(T/m) -129.04 77.1 -93.0 14.53 
M11(T.m/m) -69.98 40.29 -48.59 34.37 
M22 (T.m/m) -72.47 40.29 -50.3 34.37 
V13 (T/m) -82.23 82.23 -56.53 56.53 
V23 (T/m) -82.23 82.23 -56.53 -56.53 
Bảng2: Bảng tổng hợp kết quả nội lực tháp điều áp khi mực nước dao động 
Giá trị 
/Vị trí 
F11 F22 M11 M22 
T/m T/m T-m/m T-m/m 
Max 
value 79.88 26.61 63.61 63.62 
KHOA HỌC CÔNG NGHỆ 
 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 5
Giá trị 
/Vị trí 
F11 F22 M11 M22 
T/m T/m T-m/m T-m/m 
Area 
/Joint 
2691 
/1297 
2160 
/3 
2161 
/5 
2187 
/765 
Min 
 value -96.55 -121.31 -88.64 -91.81 
Area 
/Joint 
2076 
/665 
1615 
/217 
2160 
/3 
2082 
/705 
Giá trị 
/Vị trí 
V13 V23 F12 M12 
T/m T/m T-m/m T-m/m 
Max 
value 101.72 103.08 37.05 32.03 
Area 
/Joint 
2142 
/799 
2152 
/808 
2156 
/707 2183/757 
Min 
value -103.05 -101.71 -35.02 -32.39 
Area 
/Joint 2160/32 2134/790 
2147 
/699 2192/773 
F11 max thay đổi theo thời gian F11 min thay đổi theo thời gian 
F22 max thay đổi theo thời gian F22 min thay đổi theo thời gian 
KHOA HỌC CÔNG NGHỆ 
TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 6
M11 max thay đổi theo thời gian M11 min thay đổi theo thời gian 
M22 max thay đổi theo thời gian M22 min thay đổi theo thời gian 
Hình 7. Thay đổi nội lực tháp điều áp theo thời gian 
Qua phân tích kết quả tính toán các thành 
phần nội lực trong tháp điều áp nhận thấy 
rằng : các thành phần nội lực theo phương dọ 
trục tháp F11 và thành phần nội lực theo 
phương pháp tuyến F22 trong trường hợp 
mực nước dao động không xảy ra nguy hiểm 
so với trường hợp mực nước tĩnh, Do vậy 
nên nghiên cứu chỉ tập trung phân tích đánh 
giá sự thay đổi của các thành phần nội lực 
khác như bảng. 
Bảng 3: Chênh lệch nội lực khi xét đến dao 
động mực nước trong tháp điều áp 
K=Động/Tĩnh Max value Min value 
V13 1.24 1.25 
V23 1.25 1.24 
M11 1.58 1.27 
M22 1.57 1.27 
M12 1.49 1.51 
4. KẾT LUẬN 
Thông qua việc ứng dụng phương pháp HHT-
α phân tích ứng trạng thái ứng suất trong tháp 
điều áp khi chịu tải trọng động do mực nước 
trong tháp dao động đã nâng cao độ tin cậy của 
lời giải bài phân tích động phi tuyến, giảm 
thiểu sai số đối với các thành phần cao tần 
trong phương pháp giải bài toán dao động. Từ 
những kết quả phân tích về các giá trị cực trị 
ứng suất xuất hiện trong tháp điều áp chịu tác 
động của mực nước dao động đã có thể thấy 
rằng hệ số lệch tải đối với bài toán động trong 
nghiên cứu này nằm trong khoảng 1.2÷1.6. 
Kết quả nghiên cứu cũng cho thấy rõ diễn biến 
trạng thái ứng suất trong tháp theo thời gian 
trong quá trình vận hành góp phần minh chứng 
rõ ràng hơn nữa về ứng xử của tháp điều áp 
chịu tác động của tải trọng động. 
Những kết quả nghiên cứu này có cũng là tiền 
đề để nghiên cứu điều chỉnh hệ số lệch tải đối 
KHOA HỌC CÔNG NGHỆ 
 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 53 - 2019 7
trong bài toán phân tích kết cấu tháp điều áp 
nhằm giảm thiểu khối lượng tính toán bài toán 
động mà cẫn đảm bảo yêu cầu về kỹ thuật 
trong bài toán thiết kế. 
TÀI LIỆU THAM KHẢO 
[1] Bộ công thương, tạp chí công nghiệp - vụ năng lượng, Thủy điện Việt Nam tiềm năng và 
triển vọng phát triển. Nhà xuất bản công thương 2010. 
[2] Nguyễn Duy Thiện, Thiết kế và thi công Trạm thủy điện nhỏ. NXB xây dựng 2010. 
[3] Hoàng Đình Dũng - Phạm Hồng Nhật - Vũ Hữu Hải - Nguyễn Thượng Bằng (1991), Trạm 
thủy điện - Các công trình trên tuyến năng lượng. NXB Giao thông vận tải. 
[4] PGS. TS. Hồ Sỹ Dự - PGS. TS. Nguyễn Duy Hạnh - TS. Huỳnh Tấn Lượng - PGS. TS. 
Phan Kỳ Nam Công trình trạm thuỷ điện, NXB Xây dựng 2009. 
[5] Thomas J. R. Hughes. The Finite Element Method: Linear Static and Dynamic Finite 
Element Analysis. Prentice-Hall, Inc., Englewood Cliffs, New Jersey, 1987 edition. 
[6] Newmark, N. M. (1959) A method of computation for structural dynamics. Journal of 
Engineering Mechanics, ASCE, 85 (EM3) 67-94. 
[7] HILBER, H. M., HUGHES, T. J. R., AND TAYLOR, R. L. Improved numerical 
dissipation for time integration algorithms in structural dynamics. 
Earthquake Engineering and Structural Dynamics 5 (1977), 283-292. 
[8] Erlicher, Silvano, Luca Bonaventura, and Oreste S. Bursi. "The analysis of the generalized-α 
method for non-linear dynamic problems." Computational Mechanics 28.2 (2002): 83-104. 
[9] CiVilax.SAP2000 Manual 2014