ĐẠI HỌC QUỐC GIA HÀ NỘI
TRƯỜNG ĐẠI HỌC KHOA HỌC TỰ NHIÊN
....................................

Phan Thành Bắc

MÔ PHỎNG QUÁ TRÌNH LAN TRUYỀN VẬT CHẤT Ô
NHIỄM DƯỚI TÁC ĐỘNG CỦA CÁC YẾU TỐ ĐỘNG LỰC
TẠI VỊNH CAM RANH BẰNG MÔ HÌNH SỐ
Chuyên ngành: Hải dương học
Mã số: 60.44.97

LUẬN VĂN THẠC SĨ KHOA HỌC

Người hướng dẫn khoa học: PGS. TS Nguyễn Minh Huấn

Hà Nội


LỜI CẢM ƠN

Tác giả xin dành những lời đầu tiên bày tỏ lòng biết ơn sâu sắc nhất tới các
thầy, cô giáo trong khoa Khí tượng, Thủy văn - Hải dương học (trường Đại học
Khoa học tự nhiên Hà Nội) và các nhà khoa học tại viện Hải dương học đã tận tình
giúp đỡ, truyền thụ, trao đổi kiến thức chuyên môn cùng tác giả trong thời gian qua.
Luận văn này được hoàn thành ngoài sự nỗ lực làm việc của bản thân còn có công
rất lớn của thầy giáo PGS.TS Nguyễn Minh Huấn, người đã không ngừng đôn đốc,
động viên và truyền thụ kiến thức. Tác giả xin được gửi lời biết ơn chân thành và
sâu sắc nhất đến thầy.
Tác giả cũng xin được gửi lời cảm ơn tới ThS. Nguyễn Chí Công và tất cả
các cán bộ nghiên cứu phòng Vật Lý Biển nói riêng, Viện Hải Dương học – nơi tác

1.2.2 Phương pháp số..........................................................................................................11
1.3 MÔĐUN ECOLAB

14

1.3.1 Cơ sở lý thuyết............................................................................................................14
1.3.2 Ôxy hòa tan (DO) và nhu cầu ôxy sinh hóa (BOD)...................................................16
1.3.3 Các hợp phần của Nitơ...............................................................................................19
1.3.4 Hợp phần của Photpho...............................................................................................22
CHƯƠNG 2. TỔNG QUAN VÙNG NGHIÊN CỨU...................................................................23
2.1 TỔNG QUAN VỀ ĐIỀU KIỆN TỰ NHIÊN

23

2.1.1 Vị trí địa lí...................................................................................................................23
2.1.2 Đặc điểm gió...............................................................................................................24
2.1.3 Đặc điểm thủy, hải văn...............................................................................................24
2.1.4 Đặc điểm nhiệt - muối................................................................................................26
2.1.5 Đặc điểm dòng chảy...................................................................................................27
2.1.6 Đặc điểm thủy triều và dao động mực nước...............................................................28
2.2 ĐẶC ĐIỂM KINH TẾ - XÃ HỘI

28

2.3 HIỆN TRẠNG MÔI TRƯỜNG VỊNH CAM RANH

29

2.3.1 Các nguồn thải............................................................................................................29
2.3.2 Chất lượng nước vịnh Cam Ranh...............................................................................30


126

TÀI LIỆU THAM KHẢO

127


MỞ ĐẦU
Trong những năm gần đây, khu vực đầm thuỷ triều đang đứng trước nguy
cơ ô nhiễm nguồn nước. Đầm Thủy Triều nằm trong vịnh Cam Ranh, thuộc địa
bàn huyện Cam Lâm và thành phố Cam Ranh, tỉnh Khánh Hòa. Nơi đây phong phú
và đa dạng về số lượng cũng như trữ lượng thủy sản. Trong tương lai, đầm Thủy
Triều còn là mắt xích quan trọng trong việc phát triển du lịch của tỉnh Khánh Hòa
khi vịnh Cam Ranh đã được tỉnh này quy hoạch thành trung tâm du lịch biển tầm cỡ
quốc gia và quốc tế đến năm 2025.
Theo nhận định của người dân nơi đây, trong vòng gần chục năm nay, tôm,
cá và các loại nghêu, ốc... trên đầm thường chết hàng loạt, thậm chí “sống dai” như
loài sá sùng biển (gọi là trùn biển) cũng phải chết trắng đầy đầm, môi trường trong
đầm ngày trở nên ngột ngạt, đục ngàu, nước trong đầm có mùi hôi thối nồng nặc
theo thời gian... đã làm cho hệ sinh thái đầm bị biến dạng, nguồn lợi thủy sản cứ thế
không còn nữa. Do vậy, việc đánh khai thác các nguồn lợi trên đầm đã không còn
hiệu quả, đời sống nhân dân lại khốn khó.
Một trong những nguyên nhân gây ô nhiễm đầm là nhà máy đường Cam
Ranh. Trong quá trình vận hành nhà máy, khối nước thải từ nhà máy sau khi được
xử lí sẽ đổ ra đầm qua các cống xả thải. Các kết quả từ phân tích các mẫu nước tại
vị trí cống xả thải và khu vực xung quanh nhà máy đã ghi nhận được sự vượt
ngưỡng của các thông số môi trường xung quanh khu vực này.
Khi khối nước thải được xả ra đầm, quá trình thuỷ động lực (dòng chảy, gió,
quá trình xáo trộn,…) làm khuếch tán các chất đồng thời mang khối nước thải này

truyền một số vật chất có thể gây ô nhiễm từ các nguồn thải của khu công nghiệp,
nuôi trồng thủy sản và khu dân cư trong 2 mùa: mùa mưa và mùa khô. Trong khuôn
khổ của luận văn, mục tiêu của học viên là có thể tính toán, mô phỏng, đưa ra được
bức tranh về quá trình động lực và quá trình truyền tải các vật chất gồm BOD, DO,
NO3-, PO4+, NH3+. Một kịch bản mô phỏng sự lan truyền các vật chất ô nhiễm với
giả thiết có sự gia tăng cực đại nồng độ các chất gây ô nhiễm từ số liệu thực đo tại
cống xả thải và công suất tính tại thời điểm khảo sát từ các nguồn thải của khu công
nghiệp, nuôi trồng thủy sản và khu dân cư để có thể đánh giá mức độ lan truyền và
ảnh hưởng của các vật chất này tới chất lượng nước các bãi tắm khu vực Cam Ranh.
2


Các kết quả nghiên cứu trong luận văn góp phần bổ sung thêm các thông
tin khoa học về những nghiên cứu, đánh giá vai trò và sự tác động của các từ các
nguồn thải của khu công nghiệp, nuôi trồng thủy sản và khu dân cư tác động ngược
lại đối với các khu vực nuôi trồng thủy sản, du lịch sinh thái và các bãi tắm.

3


CHƯƠNG 1. MÔ HÌNH SỐ TRỊ
1.1 TỔNG QUAN TÌNH HÌNH NGHIÊN CỨU
1.1.1 Tổng quan tình hình nghiên cứu trên thế giới
Sử dụng các mô hình số để tính toán, mô phỏng, đánh giá chất lượng môi
trường nước khu vực gần bờ, khu bãi tắm, khu nuôi trồng thủy sản đã được thực
hiện rất phổ biến trên thế giới. Tùy thuộc vào đối tượng và mục đích nghiên cứu,
việc áp dụng các loại mô hình tính toán cũng khác nhau. Có thể liệt kê một số mô
hình thường được áp dụng để đánh giá chất lượng nước trên thế giới.
Mô hình WASP7 (Water Quality Analysis Simulation Program 7) là mô
hình được xây dựng dựa trên mô hình trước đó (WASP – được xây dựng bởi Di

hóa học). Các quá trình thiếu hụt ôxy gần tới giá trị không do các quá trình ôxy hóa,
trong đó quá trình khử nitơ như là bước tương tác đầu tiên. Tính toán thông lượng
trao đổi ôxy hòa tan và các dinh dưỡng giữa trầm tích và nước.
DELFT 3D của Viện nghiên cứu thuỷ lực Hà Lan cho phép kết hợp giữa
mô hình thuỷ lực 3 chiều với mô hình chất lượng nước. Ưu điểm của mô hình này
là việc kết hợp giữa các module tính toán phức tạp để đưa ra những kết quả tính mô
phỏng cho nhiều chất và nhiều quá trình tham gia.
SMS của Trung tâm nghiên cứu và phát triển kỹ thuật của quân đội Mỹ
xây dựng cho phép kết hợp giữa mô hình thuỷ lực 1, 2 chiều với mô hình chất lượng
nước, trong đó module RMA4 là mô hình số trị vận chuyển các yếu tố chất lượng
nước phân bố đồng nhất theo độ sâu. Nó có thể tính toán sự tập trung của 6 thành
phần bảo toàn hoặc không bảo toàn được tính toán theo lưới 1 chiều hoặc 2 chiều.
ECOHAM (phiên bản 1 và 2) là mô hình số 3D kết hợp giữa module thủy
lực với module sinh thái được phát triển bởi nhóm nghiên cứu của Trường đại học
Hamburg (Đức). Mô hình chủ yếu tính toán dựa trên chu trình của các hợp phần của
Nitơ và Photpho trong đó có tính đến cả thực vật và động vật phù du trong nước
biển.
ECOSMO (ECOSystem MOdel) là mô hình cặp ba chiều thủy động lực –
băng biển – sinh địa hóa. Mô hình được phát triển dựa trên mô hình thủy động lực
HAMSOM (HAMburg shelf Ocean Model) đã được liên kết mô đun động lực nhiệt động lực biển - băng (Schrum và Backhaus, 1999) và môđun sinh học
(Schrum, 2006). Môđun sinh học NPZD dựa trên quá trình chuyển đổi giữa mức
5


đầu tiên và thứ hai trong chuỗi thức ăn và được điểu khiển bởi các thông lượng
Nitơ, Photpho và Silic. Điều quan trọng trong tính toán mô hình này là thống nhất
được giới hạn các chu trình dinh dưỡng vĩ mô và động vật phù du như là mô hình
chuẩn đoán biến đổi cho các tương tác phi tuyến trong hệ sinh thái của các mức thứ
nhất và thứ hai trong chuỗi thức ăn. Thêm vào đó, mô hình còn tính toán sự biến đổi
các mảnh vụn và ôxy để có thể đánh giá được lượng còn lại và các quá trình ôxy

chất bổ sung (Flux) và quỹ nguồn (Budget) chạy trên phần mềm chuyên dụng
CABARET of LOICZ (Mỹ) để đánh giá mức độ tích tụ và khuếch tán vật chất tại
một số điểm thuộc vịnh Hạ Long. Sau đó, phương pháp nghiên cứu này còn được sử
dụng tính toán mức độ dinh dưỡng của hệ đầm phá Tam Giang – Cầu Hai (Thừa
Thiên Huế). Tuy nhiên, phương pháp này chưa tính toán đến quá trình khuếch tán
vật chất trong không gian và chỉ giới hạn tại một số điểm nhất định.
Hoàng Dương Tùng (2004), trong phạm vi luận án tiến sĩ, đã sử dụng phần
mềm DELFT 3D - WAQ đánh giá khả năng chịu tải ô nhiễm của Hồ Tây với mục
đích xây dựng căn cứ khoa học trong việc xây dựng kế hoạch bảo vệ và phát triển
Hồ Tây. Nội dung đã xem xét đến khả năng biến động các yếu tố DO, BOD, COD,
NH4+, NO3-, PO4- theo không gian 2 chiều và thời gian.
Trong khuôn khổ đề tài cấp Bộ Thủy sản, Trần Lưu Khanh và các cộng sự
cũng đã tiến hành nghiên cứu sức chịu tải và khả năng tự làm sạch tại khu vực nuôi
cá lồng bè ở Phất Cờ (Quảng Ninh) và Tùng Gấu (Hải Phòng) dựa trên quá trình
chuyển hóa các hợp chất dinh dưỡng, hữu cơ cũng như chế độ thủy động lực tại
thủy vực nghiên cứu.
Trong một số nghiên cứu thuộc chương trình cấp Nhà nước và cấp Bộ, các
đề tài đã triển khai theo hướng: đánh giá nguồn thải (như ô nhiễm biển do sông tải
ra, thuộc đề tài KT.03.07 - 1996), đánh giá tổn thất môi trường do các hoạt động
kinh tế gây ra với vùng ven biển... Tuy nhiên, những nghiên cứu này chưa thể hiện
được mức độ chi tiết cao trong thủy vực nhỏ và số các biến môi trường còn hạn chế,
đồng thời còn mang tính chất vĩ mô cho khu vực nghiên cứu.
Tại khu vực vịnh Cam Ranh,tuy đã có một số công trình nghiên cứu về
môi trường của các đề tài cấp nhà nước và cấp tỉnh do GS-TS. Mai Trọng Nhuận
(2008), Phạm Văn Thơm (2005,2008) đã đánh giá sơ bộ về vịnh chính hoặc hiện
trạng tại khu vực khảo sát. Gần đây nhất việc nghiên cứu liên quan tới sự truyền tải
các vật chất từ các cửa sông, các quá trình tự làm sạch môi trường do PGS-TS. Bùi
7




8

(1.1)


(1.2)

(1.3)

trong đó t là thời gian; x, y là tọa độ Đề Các; η là mực nước bề mặt; d là độ sâu của
nước tĩnh; h = η + d là độ sâu nước tổng cộng; u, v là các thành phần vận tốc theo
phương x và y; f = 2Ωsinθ là tham số Coriolis (Ω là vận tốc góc của Trái đất, θ là vĩ
độ địa lý); Error: Reference source not found tương ứng là các thành phần ứng suất
theo phương x và y tại mặt và tại đáy; g là gia tốc trọng trường; Error: Reference
source not found là mật độ nước; Error: Reference source not found,Error:
Reference source not found, Error: Reference source not found và Error: Reference
source not found là các thành phần tenxơ ứng suất bức xạ; Error: Reference source
not found là nhớt rối theo phương thẳng đứng; Error: Reference source not found là
áp suất khí quyển; Error: Reference source not found là mật độ quy ước của nước; S
là cường độ lưu lượng cung cấp cho các điểm nguồn và (Error: Reference source
not found) là vận tốc tại đó nước được đổ ra môi trường xung quanh.
Biến số có đường gạch ngang biểu thị giá trị trung bình theo độ sâu. Ví dụ,
Error: Reference source not found và Error: Reference source not found là các
thành phần vận tốc trung bình theo độ sâu được xác định bởi:
Error: Reference source not found
(1.4)
Thành phần ứng suất bên Tij (i,j = x,y) bao gồm cả ma sát nhớt, ma sát rối và
chênh lệch bình lưu. Chúng được xác định bằng sử dụng công thức nhớt rối dựa
trên những biến đổi vận tốc trung bình theo độ sâu

Trong tính toán hai chiều Error: Reference source not found là vận tốc
trung bình theo độ sâu và hệ số ma sát đáy có thể được xác định từ hệ số Chezy, C,
10


hoặc hệ số Manning, M
(1.11)
(1.12)
Ứng suất mặt
Ứng suất bề mặt Error: Reference source not found được xác định thông
qua gió bề mặt. Ứng suất mặt được tính toán dựa trên công thức thực nghiệm:
(1.13)
với Error: Reference source not found là mật độ không khí, cd là hệ số cản gió,
Error: Reference source not found là tốc độ gió ở độ cao 10m trên bề mặt biển. Vận
tốc ma sát liên hệ với ứng suất bề mặt được cho bởi công thức:
(1.14)
Hệ số cản cũng có thể là những giá trị không đổi hoặc phụ thuộc vào tốc
độ gió. Công thức bán thực nghiệm được đề xuất bởi Wu (1980, 1984) để xác định
giá trị của hệ số cản:
(1.15)
trong đó, ca, cb, wa và wb là các hệ số thực nghiệm và w 10 là tốc độ gió tại độ cao
10m trên mực nước biển. Giá trị mặc định của các nhân tố thực nghiệm là
ca=1.255x10-3, cb=2.425x10-3, wa=7m/s và wb =25m/s. Các giá trị này cho kết quả
tương đối tốt khi áp dụng cho vùng biển khơi.
1.2.2 Phương pháp số
a. Rời rạc hóa miền không gian
Miền tính được rời rạc hóa bằng phương pháp phần tử hữu hạn. Theo
phương pháp này, miền tính toán được chia nhỏ thành các phần tử liên tục không
chồng nhau. Trong không gian hai chiều, vùng tính toán có thể được rời rạc hóa
thành từng phần tử dạng đa giác, tứ giác hoặc tam giác.

dựa trên phép cầu phương tâm điểm, khi đó phương trình 1.19 được viết lại,
(1.20)
Ở đây Ui và Si tương ứng là các giá trị trung bình của U và S trên toàn bộ
phần tử thứ i và được đặt tại tâm của phần tử, NS là số cạnh của phần tử, n j véctơ
pháp tuyến ngoài đơn vị tại cạnh thứ j và Error: Reference source not foundГj là
chiều dài của giao diện thứ j.
Trong trường hợp 2D phép xấp xỉ Riemann được sử dụng để tính toán các
thông lượng đối lưu tại mặt phân cách của các phần tử. Sử dụng phép giải Roe để
ước lượng cho các biến phụ thuộc phía bên trái và bên phải của của giao diện. Độ
12


chính xác bậc hai theo không gian đạt được bằng cách sử dụng kỹ thuật tái cấu trúc
gradient tuyến tính. Các giá trị gradient trung bình được ước lượng thông qua phép
giải của Jawahar và Kamath, 2000.
Phương trình truyền tải
Các phương trình truyền tải xuất hiện trong mô hình nhiệt – muối, mô hình
rối và mô hình truyền tải. Tất cả các phương trình này đều có dạng chung. Trong
trường hợp 2D, các phương trình truyền tải có dạng tổng quát như phương trình
(1.16) trong đó

(
(1.21)
a. Tích phân theo thời gian
Các phương trình dạng tổng quát được viết:
(
(1.22)
Trong mô phỏng 2D, có hai phương pháp giải cho tích phân theo thời gian
đối với hệ phương trình nước nông và phương trình truyền tải: Phương pháp bậc
thấp và phương pháp bậc cao. Phương pháp bậc thấp là phương pháp Euler hiện bậc

hflood. Thứ hai, độ sâu tổng cộng của nước tĩnh tại cạnh có độ sâu nhỏ hơn h dry và
mực nước bề mặt tại cạnh khác đều phải lớn hơn giá trị 0.
Một phần tử được gọi là khô nếu độ sâu nước nhỏ hơn độ sâu giới hạn khô
hdry, và không một cạnh nào bị ngập. Phần tử này bị loại ra khỏi miền tính toán.
Một phần tử xem như là ngập một phần nếu nếu độ sâu nước lớn hơn h dry
và nhỏ hơn độ sâu giới hạn ướt, hoặc khi độ sâu nhỏ hơn h dry và một trong số các
cạnh khác là biên ngập nước. Trong trường hợp này thông lượng động lượng bằng
không và chỉ có thông lượng khối lượng được tính.
Một phần tử được gọi là ướt nếu độ sâu nước lớn hơn h wet. Trong trường
hợp này cả hai thành phần thông lượng khối lượng và thông lượng động lượng được
tính.
Độ sâu ướt hwet phải lớn hơn độ sâu khô giới hạn h dry và độ sâu giới hạn
ngập hflood, được xác định theo điều kiện hdry < hflood < hwet.

1.3 MÔĐUN ECOLAB
1.3.1 Cơ sở lý thuyết
Động lực học của bình lưu các biến trạng thái trong ECO Lab có thể được
14


mô tả bằng các phương trình truyền tải của vật chất không bảo toàn, có dạng:
(1.25)
trong đó:
c: Nồng độ của biến trạng thái ECO Lab
u, v: Các thành phần vận tốc dòng chảy
Dx, Dy: Các hệ số khuếch tán theo phương x và y
Sc: Nguồn sinh và nguồn mất
Pc: Các quá trình trong ECO Lab
Phương trình truyền tải có thể được viết lại:
(1.26)

thuộc vào nhiệt độ, sự phân huỷ hoá chất, sự quang hợp của tảo và v.v... Khi nồng
độ DO thấp, các loài sinh vật nước giảm hoạt động hoặc bị chết. Do vậy, DO là một
chỉ số quan trọng để đánh giá sự ô nhiễm nước của các thuỷ vực.
Quá trình cân bằng ôxy được xem xét theo các mức độ phức tạp khác nhau
của cân bằng tùy thuộc vào mục đích của người sử dụng. Có 4 mức độ khác nhau
mô tả cân bằng khối DO, trong phạm vi nghiên cứu của luận văn, chỉ tập trung vào
mức cân bằng bậc 3. Mức cân bằng này giả thiết rằng sự biến đổi của nồng độ ôxy
là tổng hợp của các quá tương tác nước - khí quyển (mặt phân cách), quá trình đạm
hóa, nhu cầu ôxy sinh hóa, quá trình quang hợp, quá trình hô hấp, nhu cầu ôxy trầm
tích (chỉ ở đáy). Các quá trình đó được mô tả bằng phương trình cân bằng sau:

(1.29)

trong đó:

16


Quá trình trao đổi ôxy giữa ôxy hòa tan trong nước và khí quyển
(g/m3/ngày). Quá trình này có tính đến mức bão hòa ôxy trong nước C s phụ thuộc
vào nhiệt độ và độ mặn.
(1.29a)
Giá trị Cs được tính thông qua biểu thức thực nghiệm sau:

Tốc độ tương tác K2 (1/s) phụ thuộc vào tốc độ gió W v, tốc độ dòng chảy
V và độ sâu nước:

Quá trình Nitrat hóa (g/m3/ngày), Y1: hệ số bổ sung cho ôxy. Đây là một
quá trình khác ảnh hưởng tới cân bằng ôxy khi ôxy được sử dụng trong quá trình
Nitrat hóa từ amoniac sang nitrite.

biểu thị lượng các chất thải hữu cơ trong nước có thể bị phân huỷ bằng các vi sinh
vật.
Dạng cân bằng của nhu cầu ôxy sinh hóa (BOD) được mô tả bằng phương
trình:
(1.30)

Giải thích các từ ngữ:
S

Độ muối (ppt)

T

Nhiệt độ (0C)

Wv

Tốc độ gió (m/s)

H

Độ sâu nước (m)
18


V

Vận tốc dòng chảy trung bình theo độ sâu (m/s)

HS_nitr


R1

Tốc độ hô hấp của thực vật ở 20oC (g O2/m2/ngày)

θ1

Hệ số nhiệt độ trong quang hợp

R2

Tốc độ hô hấp của động vật và vi khuẩn (dị dưỡng) (g O2/m2/ngày)

θ2

Hệ số nhiệt độ trong hô hấp dị dưỡng

F1(H)

Hàm hấp thụ ánh sáng

k

Hệ số hấp thụ ánh sáng (1/m)

BOD

Nồng độ BOD thực tế (mg O2/l)

K3

Quá trình này được mô tả bằng công thức
(1.31a)
Quá trình Nitrat hóa chuyển từ ammonium sang dạng nitrat. Phương trình
biểu diễn quá trình này có dạng:
(1.31b)
Quá trình hấp thụ ammonium bởi thực vật. Công thức biễu diễn phương
trình này có dạng
(1.31c)
Lượng ammoniac hấp thụ bởi vi khuẩn. Quá trình này được mô tả bằng
công thức sau:
(1.31d)
Quá trình hô hấp của sinh vật dị dưỡng, được mô tả bằng biểu thức:
(1.31e)
Giải thích các ký hiệu:
UNp

Ammoni hấp thụ bởi thực vật (mg N/mg O2)

UNb

Ammoni hấp thụ bởi vi khuẩn (mg N/mg BOD)

HS_NH3-

Nồng độ bán bão hòa của N được hấp thụ bởi vi khuẩn (mg N/l)

b. Hợp phần Nitrite (NO2-)
Dạng cân bằng của hợp phần nitrite được mô tả bằng phương trình:

20


K5

Tốc độ riêng chuyển đổi của nitrite sang nitrate ở 20oC (1/ngày)

θ5

Hệ số nhiệt độ chuyển đổi từ nitrite sang nitrate
c. Hợp phần Nitrate (NO3-)
Các nhân tố ảnh hưởng tới quá trình cân bằng khối nitrite được cho bởi

phương trình
(1.33)

Quá trình chuyển đổi từ nitrite sang nitrate. Quá trình này được mô tả
thông qua biểu thức:

21