Áp dụng các mô hình thủy lực để mô phỏng và dự báo lũ<br />
trên hệ thống sông Hồng.<br />
<br />
Nguyễn Tuấn Anh, Hoàng Văn Lai, Nguyễn Hồng Phong<br />
Viện Cơ học, Viện Khoa học và Công nghệ Việt Nam<br />
<br />
TÓM TẮT: Mô hình thủy lực đã và đang được áp dụng để mô phỏng và dự báo lũ trên<br />
các hệ thống sông. Trong các phần mềm thủy lực hiện đại, thí dụ như MIKE và HEC, thường sử<br />
dụng các hệ phương trình mô tả: sóng động lực học (dynamic wave model), mô hình sóng khuếch<br />
tán (diffusive wave model), mô hình sóng động học (kinematic wave model). Sóng khuếch tán và<br />
sóng động học là các mô hình được đơn giản hóa từ mô hình sóng động lực học và phù hợp hơn<br />
trong trường hợp số liệu không đầy đủ và chi tiết.<br />
Báo cáo này trình bày một số kết quả trong nghiên cứu và áp dụng các mô hình thủy lực<br />
trên để mô phỏng và dự báo lũ trên hệ thống sông Hồng.<br />
Báo cáo gồm 3 phần.<br />
Phần 1 giới thiệu tổng quan về các mô hình thủy lực, phương pháp giải số và các mô đun<br />
tính toán tương ứng.<br />
Phần 2 trình bày một số kết quả áp dụng các mô hình cho một số bài toán mẫu.<br />
Phần 3 giới thiệu việc thu thập và xử lý số liệu trên hệ thống sông Hồng phục vụ cho các<br />
mô hình và kết quả áp dụng để mô phỏng lũ năm 2008 và dự báo lũ năm 2009.<br />
<br />
Mở đầu.<br />
Sử dụng hợp lý tài nguyên nước nói chung, kiểm soát lũ lụt nói riêng, trên một lưu vực<br />
sông phức tạp là một vấn đề lớn phục vụ trực tiếp cho việc phát triển kinh tế - xã hội và bảo vệ<br />
môi trường. Mô hình hóa là một công cụ quan trọng trong sử dụng tài nguyên nước. Trong hệ<br />
thống mô hình cho một lưu vực sông, luôn luôn có mặt 2 loại mô hình là mô hình thủy văn và mô<br />
hình thủy lực. Trong trường hợp dòng chảy trên lưu vực sông chủ yếu được hạn chế trong lòng<br />
sông do có đê bao (thí dụ như lưu vực hệ thống sông Hồng - sông Thái Bình), thì cùng với mô<br />
hình thủy văn, mô hình thủy lực 1 chiều là các mô hình quan trọng nhất. Vì vậy, trong báo cáo<br />
này chúng tôi xin trình bày một số vấn đề liên quan đến các mô hình thủy lực 1 chiều trong mô<br />
phỏng và dự báo lũ cho hệ thống sông Hồng - sông Thái Bình.<br />
<br />
1. Giới thiệu tổng quan về các mô hình thủy lực 1 chiều<br />
Mô hình thủy lực 1 chiều, mô tả chuyển động của nước trông sông hoặc kênh hở, được<br />
xây dựng dựa trên hệ phương trình Saint Venant ([1]) đầy đủ (mô hình sóng động lực học SDLH<br />
- dynamic wave model). Một trong những dạng của hệ phương trình này là:<br />
Ac Q<br />
q (1)<br />
t x<br />
Q Q 2 Z <br />
gA Sf 0 (2)<br />
t x A x <br />
Trong hệ phương trình (1), (2): Q=Q(x,t) – lưu lượng của dòng chảy; Z=Z(x,t) – mực nước; q -<br />
lưu lượng phụ; Ac - Diện tích mặt cắt (kể cả vùng chứa); A - Diện tích chảy; Sf - Sức cản đáy<br />
<br />
n 2Q Q<br />
Sf (3)<br />
A 2R 4 / 3<br />
R là bán kính thủy lực; n - Hệ số Manning.<br />
Phương trình (1) là phương trình bảo toàn khối lượng, phương trình (2) là phương trình<br />
bảo toàn động lượng.<br />
Nếu bỏ đi 2 số hạng đầu trong phương trình bảo toàn động lượng (2), ta thu được phương<br />
trình:<br />
Z<br />
Sf 0 (4)<br />
x<br />
Mô hình được xây dựng từ các phương trình (1) và (4) là mô hình sóng khuếch tán SKT<br />
(diffusive wave model, [2], [3]). Mô hình này có thể áp dụng để mô phỏng sóng chậm như sóng<br />
lũ trong tự nhiên ([4],[5]).<br />
h Z h<br />
Nếu bỏ đi tiếp thành phần trong biểu thức = - S0 , h là độ sâu, S0 là độ dốc<br />
x x x<br />
của đáy sông, ta thu được phương trình:<br />
S0 S f (5)<br />
Mô hình được xây dựng từ các phương trình (1) và (5) là mô hình sóng động học SDH<br />
(kinematic wave model, [2], [3]). Mô hình này có thể áp dụng để mô phỏng sóng lũ chậm trên<br />
sông có độ dốc lớn ([5]). Khác với mô hình SDLH và SKT, trong mô hình SDH chi cần cho điều<br />
kiện biên lưu lượng ở thượng du. Mô hình SDH không mô tả được ảnh hưởng của mực nước ở hạ<br />
du đến dòng chảy.<br />
Trong áp dụng mô hình thủy lực 1 chiều, cùng với 3 mô hình SDLH, SKT, SDH còn 2<br />
cách tiếp cận (approach) sau được nhiều người sử dụng, đó là cách tiếp cận kiểu Muskingum và<br />
Muskingum-Cunge.<br />
Cách tiếp cận kiểu Muskingum (MUS) là cách đơn giản nhất trong áp dụng mô hình thủy<br />
lực 1 chiều. Trong cách tiếp cận này, người ta chỉ sử dụng phương trình bảo toàn khối lương (1)<br />
và để giải một phương trình gồm 2 ẩn (A, Q), MUS đã sử dụng một giả thiết quan trọng là lưu<br />
lượng Q có quan hệ với mặt cắt A:<br />
Q = f(A) (6)<br />
Đặt c = dQ/dA (c có tên gọi là vận tốc truyền sóng) và ta có thể viết phương trình (1) trong<br />
trường hợp không có lưu lượng phụ về dạng sau:<br />
dA Q Q 1 Q Q<br />
0<br />
dQ t x c t x<br />
Từ phương trình này ta thu được phương trình<br />
Q Q<br />
c 0 (7)<br />
t x<br />
Trong MUS, lưu lượng Q trong phương trình (7) được tính tương đối đơn giản qua việc lựa chọn<br />
các hệ số cho từng đoạn sông. Vì vậy, chúng tôi không áp dụng cách tiệm cận này trong mô<br />
phỏng lũ trên hệ thống sông Hồng - sông Thái Bình.<br />
Cách tiếp cận kiểu Muskingum-Cunge (MUSC) được phát triển từ MUS qua phương trình<br />
khuếch tán sau ([6])<br />
Q Q 1 2Q<br />
c 0 (8)<br />
t x 2 BS 0 x 2<br />
trong đó B là chiều rộng của sông. Phương trình (8) có thể thu nhận được từ các phương trình (1),<br />
(4) và (5). Tương tự như mô hình SDH, để tính lưu lượng Q theo MUSC, cũng chỉ cần cho điều<br />
kiên biên lưu lượng ở thượng du và MUSC cũng không mô tả được ảnh hưởng của mực nước ở<br />
hạ du đến dòng chảy.<br />
Phương pháp sai phân được sử dụng để tìm nghiệm số trị trong các mô hình thủy lực.<br />
<br />
2. Kết quả kiểm định các mô hình thủy lực cho một số bài toán mẫu<br />
2.1. Mô tả các bài toán mẫu. Có nhiều bài toán mẫu để kiểm tra khả năng mô phỏng của<br />
mô hình thủy lực 1 chiều. Trong trường hợp mô phỏng sóng lũ trong tự nhiên chúng tôi đã chọn<br />
các bài toán mẫu số 2 (dòng chảy dừng), số 6 (sóng khuếch tán) và số 7 (sóng động học) trong<br />
[5]. Bài toán mẫu số 2 có nghiệm giải tích (xem [5]) và được chọn để kiểm tra độ chính xác của<br />
nghiệm số trị thu được từ các loại mô hình. Bài toán mẫu số 6 là bài toán mẫu đặc trưng cho<br />
SKT. Bài toán mẫu số 7 là bài toán mẫu đặc trưng cho SDH.<br />
Bài toán mẫu số 2 - Dòng chảy dừng. Bài toán mẫu mô tả dòng chảy dừng trong kênh chữ<br />
nhật dài 10,000 m, có độ dốc dáy, biên lưu lượng trên thượng du, biên mực nước dưới hạ du<br />
không thay đổi theo thời gian.<br />
Bài toán mẫu số 6 - Sóng khuếch tán. Bài toán mô tả sóng lũ trên một kênh chữ nhật dài<br />
50,000 m, có độ dốc là 0.0005 với lưu lượng vào đầu kênh biến đổi từ Q0 = 750 m3/giây lên Qm<br />
=2,000 m3/giây trong khoảng thời gian là T = 1,800 giây, sau đó lưu lượng lại giảm trở lại Q0<br />
cũng trong khoảng thời gian là T.<br />
Bài toán mẫu số 7 - Sóng động học. Tương tự như sóng khuếch tán, sóng động học mô tả<br />
sóng lũ trên một kênh chữ nhật dài 100,000 m, có độ dốc lớn là 0.005 với lưu lương vào đầu kênh<br />
biến đổi từ Q0 =1,000 m3/giây lên Qm = 2,000 m3/giây trong khoảng thời gian tương đối dài T =<br />
43,200 giây (12 giờ), sau đó lưu lượng lại giảm trở lại Q0 cũng trong khoảng thời gian là 12 giờ.<br />
Tài liệu [5] chỉ mô tả nghiệm chính xác của một số bài toán mẫu (thí dụ số 1, số 2, số 9).<br />
Tài liệu [5] không mô tả nghiệm chính xác của các bài toán mẫu số 6 và số 7.<br />
Trong báo cáo này, chúng tôi xin mô tả nghiệm chính xác của bài toán mẫu số 7. Từ công<br />
thức (3) tính Sf, phương trình (5) và giả thiết rằng chu vi ước có thể lấy bằng chiều rộng B của<br />
sông, ta thu được các quan hệ sau:<br />
5/3 3/ 5<br />
2/31 A S0 5/3 S0 Q<br />
Q AR S0 k 1 A , k1 , A , (9)<br />
n nB 2 / 3 nB 2 / 3 k1 <br />
dQ 5 5 3/ 5<br />
c k1 A2 / 3 k 2Q 2 / 5 f (Q), k 2 k1 . (10)<br />
dA 3 3<br />
Từ phương trình (7) và đẳng thức (10) ta thu được phương trình tải phi tuyến sau:<br />
Q Q<br />
f (Q) 0 (11)<br />
t x<br />
Nghiệm Q của phương trình (11) cần thỏa mãn điều kiện ban đầu và điều kiện biên của bài toán<br />
mẫu số 7. Sử dụng cách biều diễn nghiệm cho phương trình đạo hàm riêng phi tuyến trong [7], ta<br />
có thể biểu diễn nghiệm Q như sau:<br />
x<br />
Q = Q(x,t) = Q0 , nếu t ≤ x1 , x1 = x1(x,Q) = (12.1)<br />
f (Q)<br />
= Q0 + (Qm - Q0)(-x1+t)/T , nếu x1≤ t ≤ x1+T (12.2)<br />
= Qm - (Qm - Q0)[(-x1+t)/T-1] , nếu x1+T ≤ t ≤ x1+2T (12.3)<br />
= Q0 , nếu x1+2T ≤ t (12.4)<br />
Áp dụng công thức tính đạo hàm của<br />
hàm ẩn ta có thể kiểm tra được ngay Q là<br />
nghiệm chính xác của bài toán mẫu số 7.<br />
Xét điểm giữa x2 = 50,000 của kênh.<br />
Tại điểm x2 này nghiệm Q bắt đầu tăng khi t =<br />
x1/f(Q0) = x1(x2,Q0)/f(Q0) ≈ 2.072 giờ. Kết quả<br />
này phù hợp với kết quả đã công bố trong [5]<br />
(hình 1). Từ thời điểm này đến x1+2T, x1 trong<br />
(12.2) và (12.3) phụ thuộc phi tuyến vào<br />
nghiệm Q. Để tìm nghiệm Q cần giải phương<br />
trình phi tuyến (12.2) và (12.3).<br />
Có thể xấp xỉ nghiệm Q này bằng biểu<br />
thức tuyến tính sau<br />
x<br />
Hình 1. NghiệmQ ([5]) tại x=0 và x=50,000 m Q’ = Q0 + (Qm - Q0)(- +t)/T (13)<br />
f (Q0 )<br />
Sai số của phép xấp xỉ này là :<br />
x x (Qm Q0 ) x f (Q0 ) f (Q)<br />
Q’ - Q = (Qm - Q0)( - )/T = *<br />
f (Q) f (Q0 ) T f (Q) f (Q0 )<br />
Với các giá trị đã cho trong bài toán mẫu số 7, ta có<br />
6<br />
thể chứng minh được │Q’ - Q│≤ 0.07 . Tương tự ta<br />
5<br />
có thể xấp xỉ nghiệm Q trong (12.3) bằng biểu thức<br />
4 tuyến tính sau<br />
H_TT<br />
3<br />
H_CX x<br />
Q’ = Q m - (Q m - Q 0 )[(- +t)/T – 1] (14)<br />
2<br />
f (Q0 )<br />
1<br />
Nghiệm xấp xỉ (13) và (14) là các tuyến tính theo t<br />
0<br />
1 8 15 22 29 36 43 50 57 64 71 78 85 92 99<br />
và phù hợp với kết quả đã công bố trong [5] (hình<br />
1).<br />
Hình 2. Kết quả kiểm định mô hình SDLH 2.2. Kết quả kiểm định các bài toán mẫu.<br />
Trong bài toán mẫu số 2 SDLH cho kết quả rất chính xác (hình 2) với sai số tuyệt đối nhỏ<br />
hơn 0.015 m. SKT cho kết quả tương đối chính xác với sai số tuyệt đối nhỏ hơn 0.15 m. MUSC<br />
và SDH trong trường hợp này cho kết quả sai (vì không mô tả được ảnh hưởng của mực nước ở<br />
hạ du đến dòng chảy)<br />
Trong bài toán mẫu số 6 2100<br />
<br />
SDLH và SKT cho các kết quả<br />
1900<br />
trùng với kết quả đã được công bố<br />
trong [5]. MUSC và SDH cho kết<br />
L u u lu o n g (m 3/s )<br />
<br />
<br />
<br />
<br />
1700<br />
quả xấu hơn 2 mô hình trên.<br />
1500<br />
Trong bài toán mẫu số 7. x =0<br />
x =50,000<br />
Mô hình SDH (hình 3) và MUSC 1300<br />
cho kết quả trùng với kết quả đã<br />
1100<br />
được công bố trong [5] trong hình<br />
1. Mô hình SDLH và mô hình SKT 900<br />
cho các kết quả giống nhau và xấu 0 21600 43200 64800 86400 108000 129600 151200<br />
Thoi g ian (g ia y)<br />
hơn kết quả của mô hình SDH và<br />
MUSC. Hình 3. Kết quả kiểm định mô hình SDH<br />
3. Kết quả áp dụng các mô hình thủy lực để mô phỏng và dự báo lũ trên hệ thống<br />
sông Hồng<br />
Theo SDLH phần mềm tính toán thủy<br />
1000 lực 1 chiều IMECH_1D cho mạng sông<br />
900<br />
800<br />
Hồng - sông Thái Bình đã được chúng tôi<br />
700 xây dựng và áp dụng để dự báo lũ từ năm<br />
Muc nuoc (cm)<br />
<br />
<br />
<br />
<br />
600 2002 đến nay. Kết quả dự báo được đánh giá<br />
Thuc do<br />
500<br />
24 gio<br />
là tương đối tốt. Hình 4 là kết quả so sánh<br />
400<br />
300<br />
48 gio mực nước thực đo với mực nước dự báo<br />
200 trước 24 giờ và trước 48 giờ tại trạm thủy<br />
100 văn Hà Nội.<br />
0<br />
0 240 480 720 960 1200 1440 1680<br />
Trên thượng du của hệ thống sông<br />
Thoi gian (gio) Hồng, việc đo đạc mặt cắt các sông hiện còn<br />
gặp một số khó khăn. Vì vậy, chúng tôi chỉ<br />
thu thập được thông tin về 2 đợt đo đạc trên thượng du sông Lô (hình 5) và các mặt cắt được trích<br />
từ bản đồ trên thượng du sông Đà (hình 6).<br />
<br />
<br />
<br />
<br />
Hình 5. Cao độ đáy sông Lô qua 2 lần đo Hình 6. Mặt cắt trên thượng du sông Đà<br />
<br />
Với các số liệu này, chúng tôi đã áp dụng SDH để mô phỏng lũ trên thượng du sông Lô<br />
năm 2008 (hình 7). Chúng tôi cũng đã kết nối mô hình thủy văn MARINE với MUSC và áp dụng<br />
áp dụng thử nghiệm để dự báo lưu lượng vào hồ Hòa Bình 2009 (hình 8).<br />
<br />
<br />
<br />
<br />
H. 7. Mô phỏng lũ trên thượng du sông Lô. H. 8. Dự báo lưu lượng vào hồ Hòa Bình, 5/7/09<br />
Kết luận.<br />
SDT và MUSC là các mô hình thủy lực đơn giản hơn so với các mô hình SDLH và SKT.<br />
Các mô hình này không mô tả được ảnh hưởng của mực nước ở hạ du đến dòng chảy. Vì vậy<br />
SDT và MUSC chỉ nên áp dụng cho các sông có độ dốc đáy tương đối lớn và số liệu địa hình<br />
không đủ chi tiết như thượng nguồn sông Lô hoặc thượng nguồn sông Đà<br />
Mô hình SDLH mô tả tương đối tốt các đặc trưng thủy lực của dòng chảy 1 chiều khi số<br />
liệu địa hình đủ chi tiết như ở hạ du hệ thống sông Hồng - sông Thái Bình.<br />
SKT đang được chúng tôi phát triển để có thể sử dụng trong dự báo và kiểm soát lũ lụt<br />
trên các hệ thống sông trong trường hợp cần dự báo nhanh với độ chính xác vừa phải.<br />
<br />
Lời cảm ơn<br />
Các tác giả của bài báo chân thành cám ơn các cán bộ tham gia thực hiện đề tài cấp nhà<br />
nước KC.08.17/06-10 ([8]) đã thảo luận và đóng góp nhiều ý kiến bổ ích về các nội dung của bài<br />
báo này.<br />
<br />
Tài liệu tham khảo<br />
[1] Cung J.A., Holly F.M. Verwey A. (1980) Practical Aspects of Computational River<br />
Hydraulics. Pitman Advanced Publishing Program.<br />
[2] MIKE11 (2004). Reference Manual<br />
[3] Hydrologic Modeling System HEC-HMS, (2000), Technical Reference Manual.<br />
[4] Hồ Trọng Tiến. (2009) “Nghiên cứu mô hình dự báo lũ sông Cửu Long dựa trên đánh giá vai<br />
trò tác động các nguồn nước”. Luận án tiến sĩ kỹ thuật. Viện khoa học thủy lợi miền Nam.<br />
[5] Lebosse A., (1993) Codes de Calcul d’ecoulement a surface libre filaire “Lido” , “Sara”, et<br />
“Zero” (version 2). Note de Validation. HE-43/92-95, EDF.<br />
[6] E. Todini, (2007) A mass conservative and water storage consistent variable parameter<br />
Muskingum-Cunge approach, Hydrol. Earth Syst. Sci., 11, 1645–1659.<br />
[7] Peter D. Lax, (1972) The Formation and Decay of Shock Wave, The American Mathematical<br />
Monthly, vol. 79, pp. 227-241<br />
[8] §Ò tµi cÊp Nhµ níc KC-08-17/06-10. (2007-2010). Hoàn thiện công nghệ dự báo lũ cho hệ<br />
thống sông Hồng, sông Thái Bình (Phát triển kết quả đề tài KC-08-13).<br />