Study of flow in Rua island on Dong Nai river using 3D numerical method

Abstract

Erosion and bank erosion are processes occurring to the riverbed and are issues closely related to the people and society. Rua Island is located in Thanh Hoi commune, Tan Uyen town, Binh Duong province which is one of serious bank erosion areas in the downstream of the Dong Nai river. Riverbank erosion in Rua islet in general has been going on for a long time and the rate of bank erosion has increased after the Tri An hydropower project came into operation and river-bed sand mining activities accelerated. At present, the Rua Islet is in danger of being cut in two at the position of Rua neck. Therefore, researching to determine the mechanism and cause of landslide of Rua Islet in order to point out effective measures to prevent landslide is urgent. A three-dimensional numerical modeling has been developed for the river section in this island to study the flow structure. This model was calibrated with observation data in 2017. Using this model flow in the entire river section was simulated. The calculation results clearly showed the three-dimensional structure of the flow at this river section, especially at the curved sections where bank erosion often occurs. Calculation results also show the distribution of bottom stresse in the river section, it is a parameter that is crucial to the process of changing of the river bed.

GIỚI THIỆU

Cù Lao Rùa (cù lao Thanh Hội) là một địa danh thuộc xã Thạnh Hội, thị xã Tân Uyên, tỉnh Bình Dương được bao bọc bởi dòng chảy chính của sông Đồng Nai và dòng chảy phụ tẻ nhánh bao trọn cù lao nhập vào dòng chính và chảy xuôi về Sài Gòn. Cù lao này là một trong những khu vực sạt lở nghiêm trọng ở vùng hạ du sông Đồng Nai. Một trong các nguyên nhân gây ra sạt lở là đoạn sông này uốn khúc quanh co làm dòng chảy có nhiều chỗ ép sát bờ và gây xói hàm ếch.

Hiện nay trên toàn bộ Cù Lao Rùa có khoảng 7 điểm sạt lở với tổng chiều dài 2.5 km, trong đó đoạn hẹp nhất của cù lao Rùa gọi là "cổ rùa" bị xói lở từ hai phía, chiều rộng chỉ còn khoảng 85m. Sạt lở không những làm mất đất, phá vỡ các cơ sở hạ tầng như đường giao thông, nhà cửa mà còn ảnh hưởng lớn đến tâm lý, đời sống, thậm chí là tính mạng của nhiều hộ dân 1 .

Để nghiên cứu dòng chảy, người ta thường dùng hai phương pháp là phương pháp thực nghiệm và phương pháp tính toán. Đối với phương pháp thực nghiệm, người ta đo các thông số của dòng chảy như vận tốc, mực nước… ngoài thực địa hoặc trên mô hình vật lý thiết lập trong phòng thí nghiệm. Trong phương pháp tính toán, người ta tiến hành giải các phương trình vi phân mô tả chuyển động và cũng tìm được các thông số kia. Tuy nhiên do hạn chế của toán học, người ta chưa có thể giải chính xác các phương trình vi phân này mà phải giải gần đúng và phổ biến nhất là giải bằng phương pháp toán số. Để giải toán cho một dòng chảy cụ thể, người ta cũng phải thiết lập mô hình cho nó và gọi là mô hình toán. Điểm yếu của phương pháp mô hình toán so với phương pháp thực nghiệm là độ chính xác nhưng nó có điểm mạnh là chi phí thấp và thời gian đi tới kết quả thường là ngắn. Cùng với sự tiến bộ của khoa học, các phương pháp tính cũng ngày càng được hoàn thiện để gia tăng độ chính xác và mức độ phức tạp của bài toán. Điều đó đã làm cho phương pháp mô hình toán trở nên ngày càng được sử dụng phổ biến.

Việc mô phỏng dòng chảy qua các đoạn sông cong bằng mô hình toán đã được thực hiện bởi các tác giả khác nhau 1 , 2 nhưng chủ yếu bằng mô hình 2D của bộ phần mềm MIKE. Ở bài báo này nhóm tác giả sẽ trình bày mô hình tính toán dòng chảy qua các đoạn sông cong bằng mô hình 3D được phát triển bởi nhóm tác giả có thể mô tả cấu trúc dòng chảy chi tiết hơn.

Mục tiêu của bài báo là trình bày một nghiên cứu về khảo sát cấu trúc dòng chảy của sông Đồng Nai ở khu vực Cù lao Rùa bằng mô hình toán số ba chiều. Nghiên cứu sẽ chỉ ra được sự gia tăng bất thường của vận tốc và ứng suất ma sát đáy ở những đoạn sông cong tăng so với những đoạn sông khác.

Phương pháp NGHIÊN CỨU VÀ MÔ HÌNH TOÁN

Phương trình cơ bản

Trong nghiên cứu này, một mô Figure 3 chiều (3D) đã được thiết lập. Đây là mô hình chi tiết nhất, cho phép nghiên cứu dòng chảy trong không gian 3 chiều. Mô hình được thiết lập bằng phần mềm F28 3 đã được tích hợp thêm module 3D. Trong mô hình này, dòng chảy được giải từ phương trình dòng chảy ba chiều với giả thiết thủy tĩnh. Trong hệ tọa độ biến đổi “sigma”, phương trình được viết:

Trong đó: – mực nước; – hai thành phần vận tốc của dòng chảy trên phương ngang; ω – thành phần vận tốc trên phương thẳng đứng “sigma”. và ; D – độ sâu; – toán tử vi phân trên mặt phẳng “sigma”; F – vector ngoại lực; AH và AV – độ nhớt rối.

Độ nhớt rối ngang, A _H , được tính theo Smagorinsky 4 :

Hệ số C nằm trong khoảng 0.01 - 0.5. Độ nhớt rối theo phương đứng, A _V , được tính theo mô hình Prandtl-Kolmogorov (1942):

Trong đó - hằng số mô hình, xác định trong quá trình hiệu chỉnh mô hình; L – chiều dài xáo trộn; k – động năng rối. Các thông số này được lấy theo Davies và Gerritsen 5 như sau:

Với 0.09 – hằng số mô hình; = 0.4 – hằng số Karman; - vận tốc ma sát trên mặt biển; - dạng sửa chữa của vận tốc ma sát đáy , trong đó ; - giá trị trung bình của các vận tốc ma sát tính từ vận tốc tại các mắt lưới trên đường thủy trực nếu biên dạng vận tốc là logarit; và - độ sâu tương đối. trong tính toán này được lấy bằng 0.1.

Trên mặt thoáng ( = 0)

Tại đáy ( = -1):

Trong đó C _D là hệ số ma sát đáy và được tính :

Với Δ – khoảng cách từ mắt lưới đầu tiên tới đáy; và z ₀ – thông số nhám. Các phương trình (1) và (2) được giải bằng phương pháp thể tích hữu hạn theo sơ đồ được giới thiệu trong tài liệu 6

Mô hình Cù Lao Rùa

Lưới tính

Đoạn sông Đồng Nai khu vực Cù lao Rùa được chia thành 7261 phần tử tứ giác ( Figure 1 ). Kích thước của các phần tử khoảng 12 – 27m. Số phần tử theo chiều ngang sông trên nhánh chính là 16 còn trên nhánh phụ là 8. Theo chiều sâu mô hình được chia thành 10 lớp. Lưới tính này là khá mịn, đủ để có thể mô phỏng một cách chi tiết cấu trúc dòng chảy khu vực nghiên cứu.

Địa hình đáy khu vực cù lao Rùa được tham khảo từ các số liệu khảo sát bởi Viện Thủy công, thực hiện năm 2017 1 . Địa hình đáy của mô hình cũng được giới thiệu trên Figure 2 .

Figure 1 . Lưới tính mô hình 3D sông Đồng Nai đoạn Cù lao Rùa tích hợp trong mô hình 1D2D hệ thống sông Đồng Nai

Điều kiện biên

Để giải quyết vấn đề điều kiện biên, mô hình khu vực Cù lao Rùa đã được tích hợp vào mô hình hệ thống sông Đồng Nai (xem Figure 1 ) kế thừa từ nghiên cứu đề xuất lựa chọn chiến lược quản lý ngập lụt thích hợp trên cơ sở các dự án đã, đang và dự kiến triển khai tại Tp.HCM của Lê Song Giang 7 . Mô hình hệ thống sông Đồng Nai là mô hình tích hợp 1D2D trong đó các sông rạch được làm mô hình 1D còn biển và các vùng đất ngập triều được làm mô hình 2D. Mô hình 1D đoạn sông tại Cù lao Rùa của mô hình hệ thống sông Đồng Nai được loại bỏ và thay bằng mô hình 3D. Hai đầu của mô hình 3D được nối vào vị trí cắt bỏ bằng liên kết siêu nút. Như vậy sau khi tích hợp thêm mô hình Cù lao Rùa thì mô hình cuối cùng là mô hình tích hợp 1D2D3D.

Các biên thượng lưu của mô hình gồm Trị An, Phước Hòa, Dầu Tiếng và Vàm Cỏ Đông, Vàm Cỏ Tây được áp đặt lưu lượng theo số liệu quan trắc. Trên biên biển từ Vũng Tàu tới Gò Công, mực nước được tính theo quy luật tương quan theo số liệu đo tại trạm Vũng Tàu.

Figure 2 . Địa hình đoạn sông Đồng Nai khu vực Cù lao Rùa (đơn vị thang màu: mét)

Các thông số mô hình

Đối với mô hình thủy lực 3D cù lao Rùa, các thông số mô hình như sau: C = 0.1; c _m = 0.09; C ' _m = 0.1 và z ₀ = 0.1-1,0mm. Các thông số này được xác định trong quá trình hiệu chỉnh mô hình.

Bước thời gian tính là Dt = 0.6s. Mặc dù bước thời gian tính khá nhỏ để đảm bảo điều kiện ổn định như tốc độ tính toán khá tốt. Trên máy PC Core i7 3.2GHz, một giờ máy tính tính được khoảng 5 giờ mô phỏng.

Hiệu chỉnh mô hình

Mô hình 1D2D hệ thống sông Đồng Nai đã được hiệu chỉnh trong nghiên cứu 7 và sẽ không được hiệu chỉnh lại. Trong nghiên cứu này chỉ thực hiện hiệu chỉnh phần mô hình 3D. Số liệu đo bởi Viện Thủy Công 1 trong thời gian từ ngày 5 – 8/8/2017 được dùng để hiệu chỉnh. Thời gian này là mùa lũ. Mặc dù trên sông không có lũ nhưng nước ở khu vực này hầu như chỉ chảy theo một chiều từ thượng lưu ra biển. Dòng chảy ngược khá yếu. Vị trí hai trạm đo được giới thiệu trên Figure 3 .

Các Figure 4 và Figure 5 giới thiệu kết quả tính mực nước và lưu lượng ở hai vị trí khảo sát. Kết quả tính khá phù hợp với số liệu thực đo.

Figure 3 . Vị trí các trạm đo mực nước và vận tốc ở Cù lao Rùa

Figure 4 . Mực nước tại cù lao Rùa tính toán trong thời gian từ 5 – 8/8/2017 (đường liền – tính; symbol – đo)

Số liệu đo của Viện Thủy Công 1 cũng cung cấp phân bố vận tốc trên mặt cắt ngang tại các trạm đo ở một số thời điểm. Figure 6 giới thiệu phân bố vận tốc tính tính toán trên mặt cắt ngang tại Trạm 1 vào thời điểm triều xuống mạnh nhất. Kết quả tính cũng khá phù hợp với số liệu đo. Với các kết quả hiệu chỉnh này, mô hình đã đạt được độ tin cậy và có thể sử dụng để tính toán.

Figure 5 . Lưu lượng tại cù lao Rùa tính toán trong thời gian từ 5 – 8/8/2017 (đường liền – tính; symbol – đo)

Figure 6 . Trường vận tốc trên mặt cắt ngang tại trạm 1 lúc 18h ngày 6/8/2019

kết quả tính toán và bàn luận

Trường vận tốc.

Figure 7 giới thiệu trường vận tốc trên mặt thoáng khi triều xuống mạnh nhất xảy ra vào thời điểm 22 giờ ngày 7/8/2017. Đây cũng là hình ảnh đặc trưng dòng chảy tại Cù lao Rùa. Trên nhánh chính và nhánh phụ có nhiều đoạn cong nên ở nhiều đoạn xuất hiện các dòng xoáy thứ cấp và các khu vực tập trung vận tốc. Hai khu vực điển hình trong số đó được chỉ ra trên Figure 7 .

Khu vực 1 là chỗ hẹp nhất của Cù lao (thường được gọi là cổ rùa). Theo 1 , khu vực này bị sạt lở liên tục và cổ rùa ngày càng bị hẹp lại. Figure 8 giới thiệu phân bố vận tốc trên các tầng mặt và tầng đáy của khu vực. Hình ảnh cho thấy trên nhánh chính tại lân cận cổ rùa hướng của vector vận tốc trên các tầng không khác nhau nhiều, tuy nhiên độ lớn của vận tốc phía bờ lõm, phía sát với cổ rùa lớn hơn rõ rệt (khoảng 25%) so với vận tốc ở phía bờ đối diện. Trên nhánh phụ, do bán kính cong nhỏ, dòng xoáy thứ cấp có cường độ rất mạnh. Hướng của vector vận tốc trên tầng đáy lệch hẳn so với hướng của vector vận tốc trên tầng mặt. Hiệu ứng 3 chiều có thể nhìn thấy rõ khi xem phân bố vận tốc trên mặt cắt ngang sông.

Figure 7 . Hình ảnh đặc trưng trường vận tốc trên mặt thoáng tại Cù lao Rùa và 2 khu vực có xoáy thứ cấp phát triển mạnh

Figure 8 . Trường vận tốc tại khu vực cổ rùa

Figure 9 cho thấy tại mặt cắt 1-1, dòng xoáy thứ cấp khá yếu. Tuy nhiên trục động lực của dòng chảy thể hiện bằng thành phần vận tốc dọc trục lệch hoàn toàn về phía bờ phải. Lưu ý là các đường đẳng giá trị thành phần vận tốc pháp tuyến được chuẩn hóa với vận tốc pháp tuyến trung bình trên mặt cắt, V _m .

Figure 9 . Thành phần vận tốc tiếp tuyến và pháp tuyến với mặt cắt 1-1.

Figure 10 giới thiệu hình ảnh dòng chảy trên mặt cắt 2-2. Tại mặt cắt này dòng xoáy thứ cấp khá mạnh. Dòng chảy có vận tốc không lớn lắm nhưng có khuynh hướng khoét vào bờ lõm.

Figure 10 . Thành phần vận tốc tiếp tuyến và pháp tuyến với mặt cắt 2-2

Trên dòng chính, dòng chảy tại khu vục 2 đổi hướng rất gấp ( Figure 11 ). Tại phía trước đoạn cong, phía bờ lồi dòng chảy có vận tốc rất lớn. Sau khúc cong, một xoáy thứ cấp được hình thành với cường độ rất mạnh. Các hình trình bày phân bố vận tốc trên các mặt cắt 3-3 ( Figure 12 ) và 4-4 ( Figure 13 ) cho thấy rõ điều này.

Figure 11 . Trường vận tốc tại khu vực 2

Figure 12 . Thành phần vận tốc tiếp tuyến và pháp tuyến với mặt cắt 3-3

Figure 13 . Thành phần vận tốc tiếp tuyến và pháp tuyến với mặt cắt 4-4

Ứng suất ma sát đáy.

Một thông số quan trọng của dòng chảy là ứng suất ma sát đáy. Thông số này có ý nghĩa quyết định đối với quá trình bồi xói lòng dẫn. Mặc dù để xác định lòng sông tại từng vị trí cụ thể sẽ bồi hay xói sẽ phải cần nhiều thông tin hơn là chỉ ứng suất ma sát đáy nhưng ở mức độ ban đầu thông số này có thể giúp đánh giá sơ bộ khả năng hoạt động của đáy. Figure 14 và Figure 15 giới thiệu kết quả tính toán ứng suất tiếp đáy tại thời điểm triều xuống mạnh nhất và triều lê mạnh nhất vào ngày triều cường (7/8/2017).

Các hình này cho thấy khi triều lên ứng suất ma sát đáy tại Cù lao Rùa khá nhỏ. Ngược lại, khi triều xuống ứng suất ma sát đáy khá lớn, đặc biệt là tại 2 khu vực 1 và khu vực 2 được nói tới ở phần trên. Ứng suất ma sát đáy tại các khu vực này có thể lên tới một vài N/m ² , vượt ứng suất giới hạn xói của các loại vật liệu đáy thông thường.

Figure 14 . Ứng suất ma sát đáy khi dòng triều xuống mạnh nhất

Figure 15 . Ứng suất ma sát đáy khi dòng triều lên mạnh nhất

kẾt luẬn

Bài báo đã trình bày nghiên cứu về cấu trúc dòng chảy ở sông Đồng Nai đoạn chảy qua Cù Lao Rùa bằng mô hình toán số ba chiều. Nghiên cứu cho thấy dòng chảy ở khu vực Cù Lao Rùa có cấu trúc 3 chiều rất rõ rệt. Tại các đoạn cong xuất hiện các xoáy thứ cấp. Kết quả tính cũng tái hiện được điều mà các nghiên cứu trước đây đã ghi nhận nhưng bằng các con số cụ thể là trước khi vào đoạn cong vận tốc phía bờ lồi tăng mạnh và sau đoạn cong thì vận tốc phía bờ lõm tăng mạnh. Kết quả tính toán cũng cho thấy ứng suất ma sát đáy tại các đoạn cong là khá lớn, vượt ứng suất tiếp ngưỡng xói của các loại vật liệu đáy thông thường. Kết quả này góp phần lý giải nguyên nhân sạt lở khu vực.

LỜI CẢM ƠN

Nghiên cứu này được tài trợ bởi Đại học Quốc Gia TpHCM trong đề tài ‘‘Nghiên cứu diễn biến bồi xói và sạt lở lòng dẫn sông Đồng Nai đoạn Cù Lao Rùa” mã số đề tài C2019-20-35.

XUNG ĐỘT LỢI ÍCH

Nhóm tác giả cam đoan rằng không có xung đột lợi ích trong công bố bài báo “Nghiên cứu bằng phương pháp toán số dòng chảy khu vực cù lao Rùa trên sông Đồng Nai”.

ĐÓNG GÓP CỦA TÁC GIẢ

Các tác giả Nguyễn Mộng Giang và Trần Thị Mỹ Hồng thu thập số liệu, dữ liệu và lập mô hình, tính toán. Tác giả Võ Lê Phú phân tích và diễn giải kết quả tính toán. Tác giả Lê Song Giang xây dựng thuật toán và hướng dẫn tính toán.