Độ sai số mô hình phần tử hữu hạn lan truyền mặn trong nước dưới đất trong không gian một chiều sử dụng phần tử tuyến tính của một số bài toán chuẩn

KHOA HỌC CÔNG NGHỆ<br /> <br /> ĐỘ SAI SỐ MÔ HÌNH PHẦN TỬ HỮU HẠN LAN TRUYỀN MẶN<br /> TRONG NƯỚC DƯỚI ĐẤT TRONG KHÔNG GIAN MỘT CHIỀU SỬ DỤNG<br /> PHẦN TỬ TUYẾN TÍNH CỦA MỘT SỐ BÀI TOÁN CHUẨN<br /> <br /> Nguyễn Văn Hoàng<br /> Viện Hàn lâm Khoa học và Công nghệ Việt Nam<br /> Nguyễn Thành Công<br /> Viện Khoa học Thủy lợi Việt Nam<br /> <br /> Tóm tắt: Bài báo trình bày kết quả đánh giá sai số giữa mô hình giải tích và mô hình (MH) số phần<br /> tử hữu hạn (PTHH) sử dụng hàm dáng tuyến tính một số bài toán chuẩn một chiều và các đánh giá<br /> sai số. Kết quả nổi bật cho thấy: 1) Với biên có nồng độ chất ô nhiễm-muối không đổi, khu vực càng<br /> gần biên độ sai số càng nhỏ (sai số tương đối không quá 0.05%-0,06%) 2) Với ranh giới giữa NDĐ<br /> bị ô nhiễm và không bị ô nhiễm, sai số của MH PTHH cao trong khoảng thời gian đầu của quá trình<br /> lan truyền, và giảm dần theo thời gian; 3) Sai số tuyệt đối và sai số tương đối của nồng độ theo MH<br /> PTHH không đồng bộ (sai số tuyệt đối lớn nhưng sai số tương đối lại nhỏ và ngược lại) có thể có<br /> mặt ở những miền có giá trị nồng độ lớn và thấp. Từ kết quả nghiên cứu đã kiến nghị về việc lựa<br /> chọn biên mô hình hợp lý để giảm thiểu sai số trong mô phỏng đồng thời việc đánh giá độ sai số kết<br /> quả MH PTHH cần lưu ý việc sử dụng sai số tuyệt đối hay sai số tương đối (phụ thuộc vào yêu cầu<br /> thực tế của bài toán lan truyền chất rất độc hại nên sử dụng giá trị sai số tương đối, và những chất<br /> không độc hại như độ mặn thì nên dùng giá trị tuyệt đối...).<br /> Từ khóa: Nước dưới đất, lan truyền ô nhiễm, nhiễm mặn, mô hình, phần tử hữu hạn.<br /> <br /> Summary: The paper presents results obtained by analytical and finite element modeling (FEM) using<br /> linear shape function of some one-dimensional problems and assessment of the errors. It has been found<br /> that: (1) For the the boundary of constant concentration of pollutantor saline, the closer to the boundary,<br /> the smaller the error is: the relative error does not exceed 0.06% and decreases very small value, less than<br /> 0,05%; (2) For the boundary between the polluted and unpolluted groundwater, the error values are<br /> greater at the initial period of speeding process, and gradually decrease as time passes; (3) The absolute<br /> error and relative error are not of the same trend (i.e. absolute error is great while relative error is less and<br /> vice versa) may be present at the areas having either high or low concentration. From the study results<br /> some proposals for selection of proper boudaries to reduce error some recommedations of the evaluation<br /> FEM methods has been made: either absolute or relative FEM concentration errors should be well<br /> associated with the requirement of the actual problems (for heavily toxic substances, it should be used the<br /> relative error, and for the non-toxic substances such as salinity absolute error should be used).<br /> Keywords: Groundwater, pollutanttransport, salinity transport, finite element model, linear element.<br /> *<br /> 1. MỞ ĐẦU trong nước, các nghiên cứu nhiễm mặn còn<br /> Nghiên cứu lan truyền mặn trong nước dưới tương đối mới mẻ, hạn chế và thông thường<br /> đất là một lĩnh vực rất quan trọng đối với nước trong các báo cáo đánh giá tài nguyên nước<br /> ta do đường bờ biển dài và nguy cơ xâm nhập dưới đất chủ yếu nêu lên biên giới phân chia<br /> mặn nước dưới đất từ nước biển rất cao. Ở vùng có nước mặn với độ tổng khoáng hóa<br /> bằng 1g/l và đề xuất rằng cần lưu ý đến khả<br /> năng nhiễm mặn vào các công trình khai thác<br /> Ngày nhận bài: 24/02/2016<br /> Ngày thông qua phản biện: 25/3/2016<br /> tại các khu vực này, nhưcác tác giả Đặng Hữu<br /> Ngày duyệt đăng: 20/4/2016 Ơn (1996, 1997), Nguyễn Trường Giang và<br /> <br /> TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016 1<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> nnk (1998), N guyễn Văn Đản và nnk (1998), trợ bởi Quỹ phát triển khoa học và công nghệ<br /> Ngô Ngọc Cát và Đoàn Văn Cánh (1998), N gô quốc gia-NAFOSTED thực hiện xây dựng<br /> Ngọc Cát và nnk (1998), Hồ Vương Bính và phần mềm mô hình phần tử hữu hạn chuyển<br /> nnk (1996). M ột số tác giả có sử dụng mô hình động nước dưới đất và lan truyền các chất ô<br /> số trong đánh giá nhiễm mặn như Phạm Quý nhiễm cũng như nhiễm mặn. Bài viết trình bày<br /> Nhân (2000), Ngô Đức Chân và Trương Thanh kết quả nghiên đánh giá sự đúng đắn và xác<br /> Cường (2005), Bùi Trần Vượng (2008), Phan định độ chính xác của phần mềm chương trình<br /> Chu Nam (2010)... Nghiên cứu sử dụng được xây dựng thông qua việc so sánh các kết<br /> phương pháp phần tử hữu hạn trong mô phỏng quả của mô hình số và kết quả lời giải tích<br /> chuyển động của nước dưới đất, lan truyền các chính xác của một số bài toán chuẩn lan truyền<br /> chất ô nhiễm và nhiễm mặn có tác giả Nguyễn các chất ô nhiễm và nhiễm mặn.<br /> Văn Hoàng qua một số công trình đã công bố 2. MỘT SỐ BÀI TOÁN CHUẨN MỘT<br /> như mô hình chuyển động nước dưới đất và CHIỀU LAN TRUYỀN CHẤT Ô NHIỄM<br /> xâm nhập mặn khu vực Nam Định (Vũ Đình NƯỚC DƯỚI ĐẤT<br /> Hùng chủ nhiệm, 2000) trong đề tài "Nghiên<br /> cứu nguyên nhân làm suy thoái chất lượng Bức tranh lan truyền các chất ô nhiễm trong<br /> nước dưới đất vùng duyên hải Nam Định, Hải nước dưới đất của một tầng chứa nước nhất<br /> Phòng và giải pháp khắc phục", tổng quan định nào đó ngoài các yếu tố quyết định là<br /> nghiên cứu nhiễm mặn đất-nước ven biển và trường vận tốc dòng nước dưới đất và các<br /> một kết quả nghiên cứu ban đầu về giải pháp thông số lan truyền vật chất của tầng chứa<br /> bổ cập nhân tạo nước dưới đất và chống xâm nước, còn bị quyết định bởi kiểu điều kiện<br /> nhập mặn bằng đê ngầm (N guyễn Văn Hoàng biên và giá trị biên. Có hai kiểu điều kiện biên<br /> và Nguyễn Thành Công, 2001), Groundwater phổ biến đối với các nguồn ô nhiễm nước dưới<br /> Artificial Recharge and Salinization đất là nguồn ô nhiễm tồn tại một khoảng thời<br /> Prevention as a Drought-Fighting M easure in gian nhất định nào đó và tồn tại mãi mãi. M ột<br /> Central Coastal Areas of Vietnam (Shaft số bài toán chuẩn lan truyền các chất ô nhiễm<br /> M eister and Nguyen Van Hoang, 2001), và nhiễm mặn với các điều kiện biên đặc trưng<br /> Ground Water Abstraction Potential and Salt được mô tả như sau:<br /> Water Intrusion Issue in Da Nang Coastal 2.1. Nguồn ô nhiễm ở dạng dải trong một<br /> Area (Phạm Quý Nhân, Nguyễn Văn Hoàng, thời gian rất ngắn<br /> 2005), đánh giá tiềm năng nước dưới đất và Đây là trường hợp tầng chứa nước dưới đất có<br /> khả năng nhiễm mặn trong quá trình khai thác chiều dày và chiều rộng không thay đổi, và<br /> trên đảo Vĩnh Thực-Quảng Ninh (N guyễn Văn cho rằng nước dưới đất không chuyển động.<br /> Hoàng và N gô Ngọc Cát, 2005), mô hình phần Ở một ranh giới của tầng tồn tại một dải mà<br /> tử hữu hạn đánh giá xâm nhập mặn công trình nước dưới đất tại đây có nồng độ chất ô<br /> khai thác nước dưới đất ven biển tỉnh Thái nhiễm hoặc muối trong nước là C0. Có thể có<br /> Bình do nước biển dâng (N guyễn Văn Hoàng hai trường hợp là: 1) biên ô nhiễm hoặc mặn<br /> và nnk, 2011)... tiếp giáp với đất đá không thấm nước, và 2)<br /> Tuy nhiên, nghiên cứ đánh giá sai số tính toán biên ô nhiễm hoặc mặn tiếp giáp với nguồn<br /> nồng độ các chất ô nhiễm hoặc muối theo kết nước mặt nhạt. H ai trường hợp này được biểu<br /> quả mô hình số chưa được phân tích đánh giá diễn dưới dạng một chiều và có lời giải giả<br /> cụ thể. Trong chương trình nghiên cứu cơ bản tích tương ứng như sau:<br /> định hướng ứng dụng, đề tài nghiên cứu mang - Tại x=0 tiếp giáp đất đá không thấm, tức là ∂C/∂x=0,<br /> mã số ĐT.NCCB-ĐHUD.2012-G/04 được tài chiều dài đới ô nhiễm hoặc bị mặn là h (hình 1).<br /> <br /> 2 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> Trong trường hợp này thì nồng độ tại một C   x h   x  h <br /> điểm bất kỳ nào đó tính từ gốc tọa độ x vào C  0 erf    erf    (1)<br /> 2   4Dt   4Dt  <br /> thời điểm t xác định theo công thức:<br /> - Tại x=0 tiếp giáp với khối nước mặt có nồng<br /> độ chất ô nhiễm hoặc muối C=0 (hình2)<br /> §Ç u kª nh<br /> <br /> <br /> <br /> <br /> ti Õp ví i h å<br /> n − íc nh ¹ t<br /> § Ç u k ª nh<br /> b Þt kÝ n<br /> <br /> <br /> <br /> <br /> ban ®Çu C 0 ban ®Çu C=0<br /> ban ®Çu C 0 ban ®Çu C=0<br /> <br /> 0 h x<br /> 0 h x<br /> <br /> <br /> Hình 1. Nguồn diện bị chặn phía biên Hình 2. Nguồn diện tiếp xúc nước nhạt phía biên<br /> <br /> <br /> Trong trường hợp này thì nồng độ tại một C ( x, t ) 1  x  ( V / R )t <br />  erfc ;<br /> điểm bất kỳ nào đó tính từ gốc tọa độ x vào C0 2  2 Dt / R <br /> thời điểm t xác định theo công thức: (4)<br /> x  ( V / R )t 2  y 2<br />  <br />  ; erfc(  )  e dy<br /> C   x   x h   xh  2 Dt / R<br /> C 0 2erf  erf erf  (2)<br /> 2   4Dt  4Dt  4Dt III. MÔ HÌNH PHÂN TỬ HỮU HẠN LAN<br /> 2.2. Biên ô nhiễm mãi mãi: TRUYỀN CHẤT Ô NHIỄM HOẶC<br /> NHIỄM MẶN NƯỚC DƯỚI ĐẤT MỘT<br /> Đây là trường hợp: tầng chứa nước dưới đất CHIỀU<br /> có chiều dày và chiều rộng không thay đổi, và<br /> Phương trình lan truyền chất ô nhiễm hoặc<br /> nước dưới đất chuyển động dọc theo chiều dài<br /> nhiễm mặn nước dưới đất một chiều có dạng<br /> của tầng. Ở một ranh giới của tầng nước dưới như sau (Jacob Bear, 1987):<br /> đất tại đây có nồng độ chất ô nhiễm hoặc<br /> muối trong nư ớc là C 0. Điều kiện ban đầu là  2C C C<br /> D 2 V R (5)<br /> nguồn ô nhiễm có nồng độ C0 tồn tại mãi mãi x x t<br /> từ thời gian ban đầu đến hết thời gian tính ở đây: D: hệ số phân tán thủy động lực học<br /> toán, tức là có điều kiện ban đầu và điều kiện theo hướng x (L2/T), C: nồng độ vật chất<br /> biên như sau: trong nước (M /L3), V : vận tốc dòng chảy nước<br /> dưới đất theo hướng x (M /T), R: hệ số chậm<br /> - Điều kiện ban đầu: t=0; C(x)=0; x=0<br /> trễ; t: thời gian (T).<br /> - Điều kiện biên: C(0,t)= C0 với t = 0 ;<br /> Phương trình (5) chỉ có lời giải duy nhất khi có<br /> C(,t) = 0 đầy đủ các điều kiện ban đầu và điều kiện biên<br /> được mô tả như sau:<br /> Điều kiện ban đầu là phân bố nồng độ của vật<br /> chất đang xem xét vào thời điểm ban đầu tùy ý<br /> t=t0 tại mọi vị trí trong miền tính<br /> toán: C  Co ( x ) (6)<br /> Hình 3. Biên là ô nhiễm hoặc nhiễm mặn<br /> Các điều kiện biên có thể là một hoặc đồng<br /> Theo Jacob Bear (1987), sau một khoảng thời<br /> thời các dạng sau:<br /> gian khi mà dòng chảy từ x=0 đã xâm nhập<br /> vào tương đối sâu ta có lời giải sau:  Biên có nồng độ đã biết: C=Cc trênc (7)<br /> <br /> <br /> TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016 3<br />  KHOA HỌC CÔNG NGHỆ<br /> <br />  Biên có gradient nồng độ pháp tuyến với cho từng nút. Khi đó có thể thiết lập công thức<br /> C gần đúng:<br /> đường biên đã biết:  J c trênJc (8)<br /> c M 1<br /> <br /> <br />  Biên có dòng vật chất pháp tuyến với biên<br />  C m (t ) N m ( x); 0  x  L (10)<br /> m 1<br /> <br /> C<br /> đã biết: D  V0(Cq  C)  qc trênqc (9) C m (t ) là giá trị gần đúng nồng độ vào thời<br /> c<br /> điểm t tại nút thứ m, và Nm(x) là hàm dáng.<br /> ở đây: V0, Cq tương ứng là vận tốc dòng chất Các điều kiện biên tại x=0 và x=L có thể cho<br /> lỏng và nồng độ vật chất của chất lỏng này qua bằng chính các giá trị tại các nút ở biên.<br /> biên.<br /> Như vậy với tất cả C1 , C 2 ,..., CM 1 là ẩn số<br /> Như vậy để có thể giải được bài toán lan chưa biết trong bước này. Sử dụng công thức<br /> truyền vật chất (kể cả lan truyền mặn không số dư trọng số:<br /> tính đến ảnh hưởng của tỷ trọng nước mặn,<br /> điều này hoàn toàn phù hợp với nồng độ muối<br /> xL<br />  Φ Φ  2Φ <br /> không cao) cần phải biết trường vận tốc.<br /> x W  R t  V x  D x2 dx  0 ;<br /> o<br /> <br /> Trường vận tốc được xác định bằng mô hình   1, 2,3..., M  1 (11)<br /> chuyển động nước dưới đất thông qua trường<br /> mực nước. Trong W là hàm trọng số và với thành phần<br /> có chứa số dư trên biên R =0<br /> Phương pháp phần tử hữu hạn giải phương<br /> trình (5) như sau (Zienkiewicz O. C. and Áp dụng định lý Green biến đổi được:<br /> Morgan K., 1983). Lấy chiều dài miền mô hình xL<br />  2 x<br /> W  L<br /> <br /> xL<br /> <br /> là L (x=0 đến L). M iền 0 L chia ra M phần tử x W x2 dx  x x x dx W  x x (12)<br /> và M+1 nút; sử dụng hàm dáng toàn cục Nm 0 0 0<br /> <br /> <br /> xL<br /> xL<br />  2 x<br /> W <br /> L<br />   <br /> Áp dụng cho:  WD dx = D  dx  DW  <br /> x0 x 2<br /> x 0<br /> x x  x  x 0<br /> <br /> <br /> <br /> Thay vào (11) ta có:<br /> xL<br /> xL<br />   W   <br />  (WR<br /> x0<br /> t<br />  WV<br /> x<br /> D <br /> x x<br /> )dx  DW<br />  x  x<br /> 0 (13)<br /> 0<br /> <br /> <br /> <br /> Thay (13) vào (10) ta có:<br /> M 1 x L<br /> dNm dW dNm x L<br /> dCm   d<br /> x L<br /> <br /> <br />    {VW dx  D dx dx }Cm (t)dx   RWNm dt dx  DW dx  =0 (14)<br /> m 1  x 0 x 0<br />  x 0<br /> <br /> <br /> <br /> Sử dụng phương pháp Galerkin (Wℓ =N  ). Viết (14) dưới dạng ma trận:<br /> M 1  xL<br /> dC  dN m dN dN m  <br /> KC +E F (15) Trong đó: K =    VN D  dx (16)<br /> dt m1  x0  dx dx dx  <br /> <br /> <br /> M 1 xL     L<br /> x<br /> <br /> E =    RN N m dx (17) F = DN  (18)<br /> m1 <br /> x 0   x  x0<br /> <br /> (1    M  1 )<br /> <br /> <br /> 4 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> Sử dụng hàm dáng tuyến tính, với kí hiệu   x  xi đối với phần tử thứ e ta sẽ có:<br />  he  <br /> N j  N ej  e<br /> ; N i  N ie  e<br /> ; h e  x j  xi (19)<br /> h h<br /> Các hàm thử nghiệm toàn cục khác 0 đối với phần tử e chỉ là những hàm Ni và Nj. Như vậy<br /> N  trong phần tử e nếu không bằng i hoặc j, tức là  không thuộc phần tử e.<br /> E E E<br /> Bởi vì: K m =  K em ; C m   CªmCm ; F    F (20)<br /> e 1 e 1 e1<br /> <br /> <br /> Có thể viết tổng quát cho phần tử e như sau:<br /> he<br />  dN j dN dN j  V D<br /> K  K   VNi<br /> e<br /> ii<br /> e<br /> ij D i dx    e (21)<br /> 0<br /> dx dx dx  2 h<br /> he<br />  dN dN dN i <br /> K  K   VN j i  D j<br /> e<br /> ji<br /> e<br /> jj<br /> dx   V  De (22)<br /> 0<br /> dx dx dx  2 h<br /> he he<br /> he he<br /> C  C   Ni N j dx  ; C ije  C eji  N i N j dx  ;<br /> e<br /> ii<br /> e<br /> jj (23)<br /> 0<br /> 3 0<br /> 6<br /> <br />      <br /> Fi1   DNi  ; F jM  DN j  ; F ie  0 , e  1 & F je  0 , e  M (24)<br />  x  0  x 1<br /> <br /> <br /> Khi có các thành phần của các ma trận phần tử Yếu tố thời gian trong phương trình (15) được<br /> Ke và Fe được hoàn toàn xác định bằng cách biến đổi như sau (Zienkiewicz O. C. and<br /> này, các ma trận K và F sẽ được hình thành Morgan K., 1983):<br /> bằng phép cộng trực tiếp tất cả các phần tử.<br /> <br />    (1   ) K   t  F t  F t  t ;<br />  C  t  t  C 1<br /> 0  1 (25)<br />  K  <br />  t   t  2<br /> <br /> Để đảm bảo độ sai số của mô hình số cần tuân 4.1. Nguồn ô nhiễm ở dạng dải trong một<br /> thủ điều kiện sau (Peter S. Huyakorn & thời gian rất ngắn, nước dưới đất không<br /> George F. Pinder, 1987): chuyển động<br />  x x 2 Dxx Trường hợp là tại x=0 tiếp giáp với khối nước<br /> Số Peclet: Pe   2  x  (28); nhạt, đều có các thông số sau:<br /> Dxx x<br />  t - D=0,02m2/ngày; Miền mô hình dài 50m;<br /> SốCourant: Cr  x  1  x   x t (26) h=10m; Điều kiện ban đầu là C= C0=1g/l trên<br /> x<br /> đoạn x=0m10m, C=0 tại x=10m50m; Điều<br /> 4. ĐÁNH GIÁ S AI SỐ KẾT QUẢ MÔ HÌNH kiện biên: C=0 tại x=0m; Bước lưới trong mô<br /> PTHH SỬ DỤNG PHẦN TỬ TUYẾN TÍNH hình PTHH là 0,02m; Bước thời gian trong mô<br /> SỐ MỘT SỐ BÀI TOÁN CHUẨN LAN hình PTHH là 0,05 ngày.<br /> TRUYỀN CÁC CHẤT Ô NHIỄM<br /> <br /> <br /> TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016 5<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> <br /> <br /> <br /> Hình 4. Nồng độ theo mô hình PTHH: x=0 Hình 5. Nồng độ theo mô hình PTHH: x=0<br /> tiếp giáp đất đá không thấm tiếp giáp nước nhạt<br /> <br /> Kết quả thể hiện trên hai hình 4 và 5 cho thấy Nồng độ tính toán theo phương pháp giải tích<br /> rằng phân bố nồng độ trong hai trường hợp sẽ không thể hiện trong bài viết, mà chỉ thể<br /> điều kiện biên tại x=0: 1) tiếp giáp với khối hiện nồng độ tính toán theo phương pháp<br /> nước nhạt, và 2) tiếp giáp với đất đá không PTHH. Theo hình 5 và 6 cho thấy phân bố<br /> thấm ở khoảng cách x lớn hơn 5m là rất tương nồng độ ở dải x=6-14m có tính đối xứng qua<br /> tự nhau. Vì vậy sẽ chỉ đánh giá độ sai số đối x=10 nên độ sai số ở hai dải x=6-10m và x=10-<br /> với trường hợp điều kiện biên tại x=0 tiếp giáp 14m là tương đương nhau (đúng như thể hiện<br /> với khối nước nhạt vì qua kết quả cũng thể trên hình 6) nên chỉ trình bày kết quả dưới<br /> hiện độ sai số đối với trường hợp điều kiện dạng bảng ở hai dải giá trị x=0-5m (bảng 1) và<br /> biên tại x=0 tiếp giáp với đất đá không thấm. x=8-14m (bảng 2).<br /> <br /> Bảng 1. Nồng độ (g/l) theo MH PTHH trường hợp x=0 tiếp giáp nước nhạt ở dải x=0-5m<br /> x (m) 10ngày 20ngày 30ngày 40ngày 50ngày 60ngày 70ngày 80ngày 90ngày 100ngày<br /> 0,00 0,00000 0,00000 0,00000 0,00000 0,00000 0,00000 0,00000 0,00000 0,00000 0,00000<br /> 0,50 0,57079 0,42391 0,35200 0,30744 0,27641 0,25321 0,23501 0,22022 0,20792 0,19747<br /> 1,00 0,88616 0,73652 0,63879 0,57092 0,52064 0,48155 0,45007 0,42399 0,40195 0,38302<br /> 1,50 0,98231 0,90652 0,82919 0,76444 0,71130 0,66726 0,63017 0,59843 0,57095 0,54686<br /> 2,00 0,99844 0,97468 0,93217 0,88625 0,84283 0,80346 0,76819 0,73661 0,70830 0,68279<br /> 2,50 0,99992 0,99482 0,97755 0,95195 0,92299 0,89354 0,86499 0,83788 0,81247 0,78874<br /> 3,00 1,00000 0,99921 0,99384 0,98232 0,96616 0,94728 0,92710 0,90654 0,88616 0,86627<br /> 3,50 1,00000 0,99991 0,99861 0,99436 0,98670 0,97617 0,96355 0,94955 0,93470 0,91940<br /> 4,00 1,00000 0,99999 0,99974 0,99844 0,99532 0,99015 0,98305 0,97432 0,96428 0,95324<br /> 4,50 1,00000 1,00000 0,99996 0,99962 0,99849 0,99615 0,99237 0,98712 0,98049 0,97267<br /> 5,00 1,00000 1,00000 0,99999 0,99989 0,99940 0,99815 0,99583 0,99229 0,98749 0,98150<br /> <br /> <br /> <br /> 6 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> Bảng 2.Nồng độ (g/l) theo MH PTHH trường hợp x=0 tiếp giáp nước nhạt ở dải x=8-14m<br /> x (m) 10ngày 20ngày 30ngày 40ngày 50ngày 60ngày 70ngày 80ngày 90ngày 100ngày<br /> 8,00 0,99926 0,98771 0,96679 0,94404 0,92246 0,90286 0,88527 0,86950 0,85534 0,84256<br /> 8,50 0,99154 0,95437 0,91603 0,88379 0,85726 0,83524 0,81668 0,80078 0,78701 0,77492<br /> 9,00 0,94490 0,87065 0,82180 0,78777 0,76251 0,74286 0,72702 0,71390 0,70280 0,69326<br /> 9,50 0,79001 0,71576 0,67927 0,65663 0,64085 0,62904 0,61978 0,61225 0,60598 0,60066<br /> 10,00 0,50630 0,50445 0,50363 0,50314 0,50281 0,50257 0,50237 0,50222 0,50210 0,50199<br /> 10,50 0,21922 0,29186 0,32728 0,34919 0,36444 0,37583 0,38477 0,39203 0,39807 0,40319<br /> 11,00 0,05872 0,13412 0,18300 0,21684 0,24187 0,26131 0,27695 0,28991 0,30085 0,31024<br /> 11,50 0,00922 0,04783 0,08683 0,11934 0,14595 0,16798 0,18651 0,20235 0,21606 0,22806<br /> 12,00 0,00082 0,01302 0,03460 0,05778 0,07962 0,09939 0,11707 0,13288 0,14705 0,15979<br /> 12,50 0,00004 0,00268 0,01149 0,02447 0,03909 0,05393 0,06828 0,08189 0,09463 0,10650<br /> 13,00 0,00000 0,00041 0,00316 0,00903 0,01722 0,02675 0,03692 0,04725 0,05747 0,06740<br /> 13,50 0,00000 0,00005 0,00072 0,00289 0,00678 0,01211 0,01846 0,02549 0,03288 0,04045<br /> 14,00 0,00000 0,00000 0,00013 0,00080 0,00238 0,00499 0,00853 0,01283 0,01771 0,02299<br /> <br /> <br /> <br /> <br /> Hình 6. Sai số tuyệt đối giữa mô hình giải Hình 7. Sai số tương đối nồng độ theoMH<br /> tích và mô hình PTHH: x=0-5m PTHH: x=0-5m<br /> <br /> <br /> <br /> <br /> Hình 8. Sai số tuyệt đối nồng độ theo MH Hình 9. Sai số tương đối nồng độ theo MH<br /> PTHH: x=0-16m PTHH: x=0-18m<br /> <br /> <br /> TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016 7<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> Kết quả cho thấy sai số tuyệt đối ở khoảng C=0 tại x=>0m1000m; Điều kiện biên: C=C0<br /> x=1-5m là dải phân tán chất ô nhiễm bị chi tại x=0 và C=0 tại x=1000m; Bước lưới trong<br /> phối bởi điều kiện biên có nồng độ chất ô mô hình PTHH là ∆x=0,5m; Bước thời gian<br /> nhiễm bằng 0 có giá trị sai số tuyệt đối trong trong mô hình PTHH là ∆t=0,1ngày;<br /> thời gian 1-10 ngày rất bé, dưới 0,0002 (khi Kích thước phần tử và bước thời gian được lựa<br /> thời gian khoảng 10 ngày) tới không quá chọn đáp ứng yêu cầu về độ chính xác cần<br /> 0,00002 (khi thời gian tới 100 ngày) (hình 6), thiết đối với các bài toán thực tế nêu trong các<br /> tương ứng với độ sai số tương đối (giá trị tuyệt công thức (28) đến (30) nêu trên.<br /> đối của tỷ số giữa sai số tuyệt đối và nồng độ<br /> tính theo phương pháp giải tích) là từ dưới Trong tính toán nồng độ C các chất ô nhiễm<br /> 0,06% đến 0,005% (hình 7). Tuy nhiên, ở dải thông thường nồng độ được so sánh với nồng<br /> x=6-14m giá trị sai số tuyệt đối lớn nhất tại vị độ C0 tại biên hoặc nguồn ô nhiễm và sử dụng<br /> trí x=10m (là ranh giới ban đầu nước ô nhiễm tỷ số C/C0 gọi là nồng độ tương đối. Trong bài<br /> và không ô nhiễm) đạt tới khoảng trên 0,06g/l toán chuẩn này sử dụng nồng độ tương đối<br /> (thời gian khoảng 10 ngày) và xuống khoảng trong so sánh độ sai số.<br /> 0,002 (thời gian khoảng 100 ngày); càng ra xa Kết quả phân bố nồng độ theo mô hình PTHH<br /> điểm x=10m này giá trị sai số càng giảm (hình cho thấy ở các bước thời gian đầu độ sai số<br /> 8), tương ứng giá trị sai số tương đối lớn nhất tương đối lớn thể hiện qua hình dáng của các<br /> là 5,2% (thời gian khoảng 10 ngày) xuống đường cong nồng độ không được mượt (hình<br /> 2,3% (thời gian khoảng 100 ngày) (hình 9). 10), đặc biệt đối với thời gian trước 8 ngày. Tuy<br /> Trong khoảng thời gian ban đầu, độ sai số lớn nhiên sau đó các đường cong phân bố nồng độ<br /> hơn so với thời gian muộn hơn do độ chênh đã mượt hơn, độ sai số đã nhỏ đi. Vì vậy, đã thể<br /> lệch nồng độ rất lớn ở ranh giới ô nhiễm (có hiện đồ thị phân bố nồng độ tương đố C/C0cho<br /> nồng độ là 1g/l) và không ô nhiễm (có nồng độ các thời gian trước 10 ngày (hình 10) và 91-100<br /> là 0) do giá trị phân tán vật lý nhỏ hơn rất nhiều ngày (hình 11), và số liệu đại diện kết quả tính<br /> sai số trong mô hình số (Huyakorn Pinder, theo MH PTHH thể hiện trong bảng 3 (1-10<br /> 1987). Theo thời gian, sự chênh lệch nồng độ ngày) và bảng 4 (91-100 ngày).<br /> trên giảm dần do quá trình phân tán vật lý nên Sai số tuyệt đối đối với hai khoảng thời gian 1-<br /> sai số của mô hình số giảm đi. Từ đây rút ra 10 ngày và 91-100 ngày được thể hiện tương<br /> mộtkinh nghiệm là trên thực tế các bài toán lan ứng trên hình 12 và 13. Độ sai số nồng độ tương<br /> truyền các chất ô nhiễm và nhiễm mặn, dải biến đối C/C0tuyệt đối đối với thời gian 1-10 ngày<br /> đổi nồng độ ở khu vực lân cận ranh giới cần tương đối lớn, cực đại đến khoảng 0,42 (thời<br /> được xác định, và đưa các giá trị này vào mô điểm 1 ngày) và xuống tới khoảng 0,1 (thời điểm<br /> hình số sẽ loại bỏ được điều kiện biên ranh giới 10 ngày) (hình 12). Vì vậy, sai số tương đối đối<br /> ô nhiễm đột ngột gây sai số lớn trong MH với khoảng thời gian 1-10 ngày là rất lớn và sẽ<br /> PTHH, và có độ sai số sẽ nhỏ đi rất nhiều. không cần thiết được trình bày. Trong khi đó,<br /> 4.2. Biên tiếp xúc ô nhiễm mãi mãi và nước khi thời gian lan truyền đã kéo dài tới 91-100<br /> dưới đất có chuyển động: ngày thì sai số giảm đi rất đáng kể, cực đại chỉ<br /> còn khoảng 0,0053 (thời điểm 91 ngày) xuống<br /> Trường hợp biên tiếp giáp với nguồn ô nhiễm còn 0,005 (thời điểm 100 ngày), tương ứng sai<br /> có nồng độ không đổi mãi mãi có các thông số số tương đối là khoảng 1,5% (thời điểm 91<br /> như sau: ngày) xuống khoảng 1,2% (thời điểm 100 ngày)<br /> - D=1m2/ngày, V=1m/ngày; Miền mô hình dài ở dải nồng độ tương đối C/C0 lớn hơn 0,8<br /> 1000m; Điều kiện ban đầu là C=C0 tại x=0m, (x=115-120m), tăng lên bằng 4%-6% ở dải nồng<br /> <br /> 8 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> độ tương đối C/C0nhỏ hơn 0,02 (x>135m). Một thời gian ban đầu độ sai số lớn hơn so với thời<br /> điều đáng được lưu ý là sai số tuyệt đối sau khi gian muộn hơn, ngoài nguyên nhân do độ<br /> giảm tới rất nhỏ (dưới 0,0001) (ở x=~95m đối chênh lệch nồng độ rất lớn ở ranh giới ô nhiễm<br /> với t=91 ngày và 105m đối với t=100 ngày) đã do giá trị phân tán vật lý nhỏ hơn rất nhiều sai<br /> có xu thế tăng lên (tới 0,0013 ở x=~105m ở thời số trong mô hình số về phân tán, còn do nguyên<br /> điểm t=91 ngày, tới và 0,0011 ở x=~115m ở thời nhân do nước dưới đất chuyển động, mà chiều<br /> điểm t=100 ngày), nhưng do giá trị nồng độ ở dài dòng nước chuyển động trong bước thời<br /> đây thấp nên sai số tương đối rất lớn tới 4%-6%. gian không bằng tích của một số nguyên với<br /> Tương tự như trường hợp trên, trong khoảng bước thời gian (Huyakorn Pinder, 1987).<br /> <br /> <br /> <br /> <br /> Hình 10. Nồng độ tương đố C/C0theo mô hình Hình 11. Nồng độ tương đố C/C0 theo mô hình<br /> PTHH: thời gian 1-10 ngày PTHH: thời gian 91-100 ngày<br /> Bảng 3. Nồng độ tương đối C/C0 theo mô hình PTHH phần tử tuyến tính: 1-10ngày<br /> x (m) 1ngày 2ngày 3ngày 4ngày 5ngày 6ngày 7ngày 8ngày 9ngày 10ngày<br /> 0 1,000 1,000 1,000 1,000 1,000 1,000 1,000 1,000 1,000 1,000<br /> 1 0,9163 0,8375 0,9679 0,9501 0,9859 0,9812 0,9934 0,9921 0,9968 0,9965<br /> 2 0,3265 0,7835 0,7906 0,9197 0,9204 0,9653 0,9664 0,9839 0,9849 0,9922<br /> 3 0,1164 0,4521 0,6973 0,7720 0,8748 0,9000 0,9433 0,9540 0,9730 0,9780<br /> 4 0,0415 0,2227 0,4881 0,6538 0,7554 0,8401 0,8823 0,9232 0,9423 0,9620<br /> 5 0,0148 0,1013 0,2919 0,4949 0,6328 0,7379 0,8146 0,8656 0,9057 0,9308<br /> 6 0,0053 0,0439 0,1577 0,3322 0,4956 0,6210 0,7215 0,7952 0,8498 0,8906<br /> 7 0,0019 0,0184 0,0795 0,2023 0,3561 0,4962 0,6126 0,7072 0,7796 0,8353<br /> 8 0,0007 0,0076 0,0381 0,1143 0,2356 0,3717 0,4972 0,6057 0,6952 0,7663<br /> 9 0,0002 0,0031 0,0176 0,0610 0,1451 0,2602 0,3831 0,4982 0,5998 0,6851<br /> 10 0,0001 0,0012 0,0079 0,0311 0,0843 0,1710 0,2792 0,3921 0,4991 0,5947<br /> 11 0,0000 0,0005 0,0035 0,0153 0,0466 0,1062 0,1925 0,2943 0,3995 0,4997<br /> 12 0,0000 0,0002 0,0015 0,0073 0,0247 0,0628 0,1261 0,2106 0,3068 0,4056<br /> 13 0,0000 0,0001 0,0006 0,0034 0,0127 0,0356 0,0788 0,1438 0,2259 0,3174<br /> 14 0,0000 0,0000 0,0003 0,0015 0,0063 0,0194 0,0472 0,0940 0,1595 0,2391<br /> 15 0,0000 0,0000 0,0001 0,0007 0,0031 0,0103 0,0272 0,0590 0,1082 0,1735<br /> 16 0,0000 0,0000 0,0000 0,0003 0,0015 0,0053 0,0152 0,0357 0,0706 0,1213<br /> 17 0,0000 0,0000 0,0000 0,0001 0,0007 0,0027 0,0083 0,0209 0,0445 0,0818<br /> 18 0,0000 0,0000 0,0000 0,0001 0,0003 0,0013 0,0044 0,0119 0,0271 0,0534<br /> 19 0,0000 0,0000 0,0000 0,0000 0,0001 0,0006 0,0023 0,0066 0,0161 0,0337<br /> 20 0,0000 0,0000 0,0000 0,0000 0,0001 0,0003 0,0011 0,0035 0,0093 0,0207<br /> <br /> <br /> TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016 9<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> Bảng 4 4. Nồng độ tương đố C/C0theo mô hình PTHH phần tử tuyến tính: 91-100ngày<br /> x (m) 91ngày 92ngày 93ngày 95ngày 95ngày 96ngày 97ngày 98ngày 99ngày 100ngày<br /> 50 0,999 0,999 0,999 1,000 1,000 1,000 1,000 1,000 1,000 1,000<br /> 55 0,997 0,997 0,998 0,998 0,999 0,999 0,999 0,999 0,999 0,999<br /> 60 0,989 0,991 0,992 0,993 0,994 0,995 0,996 0,997 0,997 0,998<br /> 65 0,973 0,977 0,980 0,983 0,985 0,987 0,989 0,991 0,992 0,993<br /> 70 0,940 0,948 0,954 0,960 0,965 0,970 0,974 0,977 0,980 0,983<br /> 75 0,882 0,895 0,907 0,917 0,927 0,935 0,943 0,950 0,956 0,961<br /> 80 0,793 0,812 0,830 0,846 0,862 0,876 0,889 0,901 0,912 0,921<br /> 85 0,672 0,697 0,721 0,744 0,766 0,786 0,806 0,823 0,840 0,856<br /> 90 0,530 0,559 0,587 0,615 0,642 0,667 0,692 0,716 0,739 0,760<br /> 95 0,383 0,412 0,442 0,471 0,500 0,529 0,557 0,585 0,612 0,638<br /> 100 0,252 0,278 0,304 0,331 0,358 0,386 0,415 0,443 0,472 0,500<br /> 105 0,150 0,169 0,189 0,211 0,234 0,258 0,283 0,309 0,335 0,362<br /> 110 0,080 0,092 0,106 0,122 0,138 0,156 0,175 0,196 0,217 0,240<br /> 115 0,038 0,045 0,053 0,063 0,073 0,085 0,098 0,112 0,128 0,144<br /> 120 0,016 0,019 0,024 0,029 0,035 0,042 0,049 0,058 0,068 0,079<br /> 125 0,006 0,007 0,009 0,012 0,015 0,018 0,022 0,027 0,032 0,038<br /> 130 0,002 0,003 0,003 0,004 0,006 0,007 0,009 0,011 0,014 0,017<br /> 135 0,001 0,001 0,001 0,001 0,002 0,002 0,003 0,004 0,005 0,007<br /> 140 0,000 0,000 0,000 0,000 0,001 0,001 0,001 0,001 0,002 0,002<br /> <br /> <br /> <br /> <br /> Hình 12. Sai số tuyệt đối nồng độ tương đố Hình 13. Sai số tuyệt đối nồng độ tương<br /> C/C0theoMH PTHH thời điểm 1-10ngày đốiC/C0 theo MH PTHH thời điểm 91-100ngày<br /> <br /> 4. NHẬN XÉT, KẾT LUẬN-KIẾN NGHỊ nhiễm và không bị ô nhiễm, sai số của mô<br /> - Đối với diện gần biên có nồng độ chất ô hình PTHH cao trong khoảng thời gian đầu<br /> nhiễm-muối không đổi, khu vực càng gần biên của quá trình lan truyền, và giảm dần theo thời<br /> độ sai số càng nhỏ; Sai số tương đối không gian: trong thời gian ban đầu (10 ngày trở về<br /> quá 0,06% và giảm xuống giá trị rất nhỏ dưới trước) sai số tương đối tới trên 5% (lớn nhất),<br /> 0,05%; nhưng giảm xuống tới khoảng 1,5% (lớn nhất)<br /> từ ngày thứ 100.<br /> - Đối với ranh giới giữa nước dưới đất bị ô<br /> <br /> 10 TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016<br />  KHOA HỌC CÔNG NGHỆ<br /> <br /> bài toán thực tế là một dải biến đổi nồng độ đều đều<br /> (không đột ngột từ rất cao xuống rất thấp). Điều này<br /> sẽ làm giảm thiểu độ sai số của kết quả MH PTHH;<br /> - Cũng tại khu vực được cho là ranh giới giữa<br /> khu vực ô nhiễm và khu vực không bị ô nhiễm<br /> của tầng chứa nước, nên xác định sự thay đổi<br /> nồng độ trên thực tế vào thời điểm ban đầu bắt<br /> đầu mô hình, mà trên hầu hết các bài toán thực<br /> tế là một dải biến đổi nồng độ đều đều;<br /> - Trong quá trình hiệu chỉnh mô hình, nên lựa<br /> Hình 14. Sai số tương đối nồng độ tương đốiC/C0 chọn các vị trí có nồng độ chất ô nhiễm không<br /> theo MH PTHH thời điểm 91-100ngày thay đổi quá lớn để sử dụng trong hiệu chỉnh<br /> các thông số của mô hình, đồng thời khoảng<br /> - Sai số tuyệt đối và sai số tương đối của nồng<br /> thời gian giữa các thời điểm dùng hiệu chỉnh<br /> độ theo MH PTHH không đồng bộ (tức là sai<br /> đủ dài có sự thay đổi nồng độ tương đối lớn;<br /> số tuyệt đối lớn nhưng sai số tương đối lại nhỏ<br /> và ngược lại) có thể có mặt ở những miền có - Trong các bài toán thực tế, việc đánh giá độ sai<br /> giá trị nồng độ lớn và thấp. số kết quả MH PTHH cần được lưu ý tới việc sử<br /> dụng sai số tuyệt đối hay sai số tương đối tùy<br /> Từ các kết quả trình bày ở trên, một số kiến<br /> thuộc vào yêu cầu thực tế của bài toán (đối với<br /> nghị kết luận đối với MH PTHH lan truyền các<br /> những chất rất độc hại nên sử dụng giá trị sai số<br /> chất ô nhiễm và nhiễm mặn như sau:<br /> tương đối với giá trị nồng độ cho phép làm cơ<br /> - Tại khu vực gần biên ô nhiễm và nhiễm mặn, nên sở, và đối với những chất không độc hại như độ<br /> xác định sự thay đổi nồng độ trên thực tế vào thời mặn thì dùng giá trị tuyệt đối...).<br /> điểm ban đầu bắt đầu mô hình, mà trên hầu hết các<br /> <br /> TÀI LIỆU THAM KHẢO:<br /> <br /> [1] Wen-Hsiung Li, 1972. Differential equations of hydraulic transients, dispersion, and<br /> groundwater flow. Mathematical Methods in Water Resources. Prentice-Hall, Inc.,<br /> Englewood Cliffs, New Jersey.<br /> [2] Zienkiewicz O. C. and Morgan K., 1983. Finite elements and approximation. John Wiley & Sons.<br /> [3] Peter S. Huyakorn & George F. Pinder, 1987. Computational methods in subsurface flow.<br /> Academic Press, Inc. Harcourt Brace Jovanivich, Publishers.<br /> [4] Huyakorn Pinder, 1987, Computational method in subsurface flow, Academic Press.<br /> [5] Jacob Bear and Arnold Verruijt, 1987, Modeling groundwater flow and pollution, D.<br /> Reidel Publishing Company, Dordrecht, Holand.<br /> [6] Jacob Bear, 1979, Hydraulics of groundwater, McGraw-Hill Publishing Company.<br /> [7] Bear, J., 1972. Dynamics of Fluids in Porous Media, Elsevier, New York.<br /> [8] Nguyễn Văn Hoàng, 2014. Báo cáo chuyên đề thuộc ĐT.NCCB-ĐHƯD.2012-G/04: Xây<br /> dựng mô hình bài toán chuẩn lan truyền các chất ô nhiễm trong nước dưới đất sử dụng phần<br /> tử tuyến tính.<br /> <br /> <br /> TẠP CHÍ KHOA HỌC VÀ CÔNG NGHỆ THỦY LỢI SỐ 32 - 2016 11<br />