Tài liệu Nghiên cứu thử nghiệm một số phương pháp nội suy trong xử lý số liệu thực nghiệm

ĐẠI HỌC QUỐC GIA HÀ NỘI TRƯỜNG ĐẠI HỌC KHOA HỌC TỰ NHIÊN --------------------- Nguyễn Thị Thuỳ Dương NGHIÊN CỨU THỬ NGHIỆM MỘT SỐ PHƯƠNG PHÁP NỘI SUY TRONG XỬ LÝ SỐ LIỆU THỰC NGHIỆM LUẬN VĂN THẠC SĨ KHOA HỌC Hà Nội – 1/2015 1 ĐẠI HỌC QUỐC GIA HÀ NỘI TRƯỜNG ĐẠI HỌC KHOA HỌC TỰ NHIÊN --------------------- Nguyễn Thị Thuỳ Dương NGHIÊN CỨU THỬ NGHIỆM MỘT SỐ PHƯƠNG PHÁP NỘI SUY TRONG XỬ LÝ SỐ LIỆU THỰC NGHIỆM Chuyên ngành: Vật lý địa cầu Mã số: 60.44.0111 LUẬN VĂN THẠC SĨ KHOA HỌC NGƯỜI HƯỚNG DẪN KHOA HỌC: TS. Nguyễn Đức Vinh Hà Nội – 1/2015 MỤC LỤC Mở đầu Các ký hiệu và chữ viết tắt i Danh mục các hình vẽ ii Danh mục các bảng biểu iv Chương 1 Các phương pháp nội suy thông dụng 5 1.1 Nội suy chính xác tại các điểm nút 6 1.1.1 Nội suy Lagrange 7 1.1.2 Nội suy Newton 8 1.1.3 Nội suy Gauss 13 1.1.4 Nội suy Sterling 14 1.1.5 Nội suy Bessel 14 1.1.6 Nội suy Spline 15 1.2 Nội suy xấp xỉ tại các nút 16 1.2.1 Phương pháp các điểm lụa chọn 17 1.2.2 Phương pháp trung bình 18 1.2.3 Phương pháp bình phương nhỏ nhất 19 Chương 2 Các phép nội suy trong phần mềm SURFER 24 2.1 Vài nét về phần mềm SURFER 24 2.2 Các phương pháp nội suy trong phần mềm SURFER 27 2.2.1 Phương pháp nghịch đảo khoảng cách (Inverse Distance to a 2.2.2 Power) 28 Phương pháp Shepard (Shepard`s Method) 32 2.2.3 Phương pháp lân cận gần nhất (Nearest Neighbor ) 33 2.2.4 Phương pháp trung bình cửa sổ trượt (Moving Average ) 33 2.2.5 Phương pháp hồi qui đa thức (Polynomial Regression ) 34 2.2.6 Phương pháp đa thức địa phương (Local Polynomial ) 35 2.2.7 Phương pháp độ cong tối thiểu (Minimum Curvature ) 36 2.2.8 Phương pháp hàm xuyên tâm cơ bản (Radial Basic Function) 37 2.2.9 Phương pháp Kriging 39 Chương 3 Thử nghiệm một số phép nội suy 41 3.1 Trường hợp nội suy bằng hồi qui đa thức 41 3.1.1 Trường hợp bài toán một biến 41 3.1.2 Trường hợp bài toán hai biến 45 3.2 Trường hợp nội suy không sử dụng hồi qui 52 Kết luận Tài liệu tham khảo Các ký hiệu và chữ viết tắt  - độ dẫn điện [đơn vị Siemens/m]  - hằng số điện môi  - độ thẩm từ  - độ suy giảm của sóng điện từ [dB/m2] E - cường độ điện trường [Volt/m] D - véc tơ cảm ứng điện i DANH MỤC CÁC HÌNH VẼ Trang Hình 1.1. Mô tả phương pháp hồi qui chính xác tại các nút ..................................... 6 Hình 1.2. Kết quả nội suy Spline bằng phần mềm Maple ..................................... 16 Hình 1.3. Ví dụ nội suy xấp xỉ ............................................................................... 17 Hình 1.4. Các điểm thực nghiệm và đương cong hàm thực nghiệm ....................... 21 Hình 2.1: Cửa sổ chính của phần mềm Surfer ........................................................ 24 Hình 2.2. Các tiểu thực đơn của chức năng Grid .................................................... 25 Hình 2.3. Các tiểu thực đơn của chức năng Map ................................................... 26 Hình 2.4. Ví dụ bản đồ dạng các đường đẳng trị .................................................... 26 Hình 2.5 Ví dụ hình vẽ dạng 3D ............................................................................ 27 Hình 2.6. Sơ đồ mô tả phương pháp nghịch đảo khoảng cách ................................ 28 Hình 2.7. Cửa sổ chính của thực đơn Grid ............................................................. 30 Hình 2.8. Cửa sổ Option của nội suy Inverse Distance to a Power ......................... 30 Hình 2.9. Minh hoạ sử dụng đường Breakline ....................................................... 31 Hình 2.10. Minh hoạ dùng đường Fault ................................................................. 32 Hình 2.11. Minh hoạ việc điều chỉnh thông số Anisotropy .................................... 32 Hình 2.12. Cửa sổ điều chỉnh thông số của phương pháp hàng xóm gần nhất ........ 33 Hình 2.13. Cửa sổ của phương pháp trung bình cửa sổ trượt.................................. 34 Hình 2.14. Cửa sổ chọn thông số của phương pháp đa thức địa phương ................ 36 Hình 2.15. Cửa sổ chọn thông số của phương pháp độ cong tối thiểu .................... 37 Hình 2.16. Cửa sổ chọn thông số của phương pháp Radial Basic Function ............ 38 Hình 2.17. Cửa sổ chọn thông số của phương pháp Kriging .................................. 39 Hình 3.1. Đường cong trước hồi qui ..................................................................... 42 Hình 3.2. Đường cong trước và sau khi hồi qui...................................................... 43 Hình 3.3. Đường cong trước hồi qui ...................................................................... 44 Hình 3.4. Đường cong trước và sau khi hồi qui...................................................... 45 Hình 3.5. Bản đồ đẳng trị theo số liệu trên bảng 3.3 .............................................. 47 Hình 3.6. Bản đồ đẳng trị theo số liệu nội suy từ số liệu trên bảng 3.3 ................... 47 ii Hình 3.7. Bản đồ đẳng trị theo số liệu nội suy từ số liệu trên bảng 3.4 ................... 49 Hình 3.8. Bản đồ đẳng trị theo số liệu nội suy từ số liệu trên bảng 3.5 ................... 51 Hình 3.9a. Bản đồ đảng trị mảng số liệu a ............................................................ 52 Hình 3.9b Bản đồ đảng trị mảng số liệu b ............................................................. 52 Hình 3.10. Dùng phương pháp nghich đảo khoảng cách với 2 số liệu hình 3.9 ...... 53 Hình 3.11. Dùng phương pháp Kriging với 2 số liệu hình 3.9 ................................ 53 Hình 3.12. Dùng phương pháp độ cong tối thiểu với 2 số liệu hình 3.9 .................. 54 Hình 3.13. Dùng phương pháp Shepard với 2 số liệu hình 3.9 ............................... 54 Hình 3.14. Dùng phương pháp lân cận gần nhất (Nearest neighbor) ...................... 55 Hình 3.15. Dùng phương pháp hàm xuyên tâm cơ bản (Radial Basic Function).....55 Hình 3.16. Dùng phương pháp trung bình cửa sổ trượt .......................................... 56 Hình 3.17. Dùng phương pháp đa thức địa phương bậc 2 ...................................... 56 Hình 3.18. Dùng phương pháp đa thức địa phương bậc 3 ...................................... 57 Hình 3.19. Sơ đồ các điểm số liệu ban đầu trên 2 mảng số liệu hình 3.9 ................ 57 Hình 3.20. Dùng phương pháp nghịch đảo khoảng cách và Kriging...................... 58 Hình 3.21. Dùng phương pháp độ cong nhỏ nhất và phương pháp Shepard .......... 59 Hình 3.22.Phương pháp lân cận gần nhất và phương pháp xuyên tâm cơ bản. ....... 59 Hình 3.23. Phương pháp nội suy đa thức bậc 2 và 3............................................... 59 iii DANH MỤC CÁC BẢNG BIỂU Trang Bảng 1.1. ................................................................................................................. 5 Bảng 1.2. ................................................................................................................. 8 Bảng 1.3. ............................................................................................................... 11 Bảng 1.4. ............................................................................................................... 15 Bảng 1.5. ............................................................................................................... 20 Bảng 2.1. Ví dụ tệp số liệu chuẩn bị cho Grid ........................................................ 25 Bảng 3.1. Số liệu mô hình tuyến tính ..................................................................... 42 Bảng 3.2. Số liệu mô hình phi tuyến ...................................................................... 44 Bảng 3.3. Số liệu tính lý thuyết .............................................................................. 46 Bảng 3.4. Số liệu tính lý thuyết được cắt bớt ........................................................ 48 Bảng 3.5. Số liệu đã cài nhiễu................................................................................ 50 Bảng 3.6. Hàm lượng chì khu vực X...................................................................... 58 iv MỞ ĐẦU Xử lý số liệu là công việc không thể tránh khỏi trong công tác khảo sát, thực nghiệm. Một trong những khâu xử lý số liệu là nội suy các giá trị theo mạng lưới cần thiết. Toán học tính toán đã cung cấp cho lĩnh vực xử lý số liệu rất nhiều thuật toán nội suy khác nhau. Việc tìm hiểu để ứng dụng một cách hiệu quả những thuật toán này cũng là một bước quan trọng trong qui trình xử lý số liệu. Tin học và máy tính phát triển cũng góp phần đẩy mạnh việc ứng dụng những thuật toán phức tạp hơn, mạnh mẽ hơn. Các hãng sản xuất phần mềm ngày càng đưa ra các phần mềm hoàn thiện hơn, chuyên nghiệp hơn để phục vụ nhu cầu của các nhà xử lý. Cũng như nhiều quốc gia trên thế giới, Việt Nam hàng năm đầu tư một khoản tiền khổng lồ cho công tác điều tra, khảo sát nói chung. Một số lượng lớn thông tin được thu nhập, những thông tin ấy đòi hỏi được xử lý tốt hơn, nhanh hơn, chính xác và rẻ hơn. Những năm gần đây ở các Viện, các Trường, trong đó có Trường Đại học Khoa học Tự nhiên, nhiều đơn vị của Trường đã đóng góp một phần đáng kể trong nghiên cứu, xây dựng các bộ chương trình phục vụ công tác xử lý số liệu. Lĩnh vực nghiên cứu Trái đất cũng như thăm dò khoáng sản rất quan tâm và phát triển mảng nghiên cứu này. Trong khuôn khổ của luận văn này , do thời gian có hạn và hạn chế về mặt kiến thức, học viên được giao đề tài “Nghiên cứu thử nghiệm một số phương pháp nội suy trong xử lý số liệu thực nghiệm”. Nội dung chính của bài luận văn này là tóm tắt hệ thống các phép nội suy thông dụng, nhất là các phép nội suy được cài đặt trong các phần mềm xử lý số liệu như bộ phần mềm SURFER. Bản luận văn được chia làm ba phần: phần thứ nhất trình bày về các phép nội suy dạng đường thông dụng, phần thứ hai trình bày về các phép 3 nội suy dạng mặt được cài đặt trong phần mềm SURFER. Phần cuối trình bày một số kết quả thử nghiệm. Tôi xin chân thành cảm ơn sự giúp đỡ nhiệt tình của các thầy cô giáo ở bộ môn Vật lý địa cầu. Xin cảm ơn các bạn cùng khóa đã tận tình giúp đỡ để chúng tôi có thể hoàn thành luận văn này. 4 CHƯƠNG 1. CÁC PHƯƠNG PHÁP NỘI SUY THÔNG DỤNG Một hàm f(x) bất kỳ có thể được biểu diễn bằng bảng số liệu, đồ thị hoặc công thức. Nếu hàm được cho dưới dạng một công thức thì việc biểu diễn nó dưới dạng bảng hay đồ thị là khá dễ dàng. Trong công tác đo đạc hay thực nghiệm nói chung, hàm chỉ được cho dưới dạng bảng số liệu. Giả sử hàm y=f(x) được cho dưới dạng bảng sau: Bảng 1.1 Trong bảng trên, xi là đối số, yi là hàm số, i là chỉ số. Bài toán nội suy là việc xác định giá trị hàm y ở những chỗ chưa có giá trị của đối số ở trong bảng trên. Nếu hàm y=f(x) được cho dưới dạng công thức, khi ấy với mọi giá trị đối số x ta đều có thể tính được hàm y. Trong trường hợp hàm cho dưới dạng bảng như của chúng ta, bài toán nội suy có thể đặt ra như sau: Cho một số giá trị hàm tại một số giá trị của đối số, cần xác định hàm trên toàn bộ đoạn có chứa các giá trị đối số đã có giá trị hàm. Nói đơn giản hơn, nội suy là xác định hàm tại vị trí không đo được hoặc chưa đo được dựa vào các giá trị đã đo được xung quanh nó. Trong thực nghiệm nói chung và vật lý địa cầu nói riêng, số liệu đo đạc là những bảng số liệu như bảng 1.1. Bài toán nội suy rất quan trọng trong lĩnh vực mô hình hóa, xử lý số liệu thực nghiệm, chuyển các hàm phức tạp sang dạng đơn giản thuận tiện hơn cho các tính toán … 5 Cho đến ngày nay, tồn tại rất nhiều phương pháp giải bài toán nội suy. Có thể phân loại [1,2,4,7] phương pháp theo độ chính xác tại các nút nội suy, theo dạng hàm lựa chọn, theo thủ thuật toán học v.v… 1.1. Nội suy chính xác tại các điểm nút Nội suy chính xác tại các điểm nút là phép nội suy sao cho giá trị hàm đã có tại các điểm nút không bị thay đổi. Trên hình 1.1, các điểm tô đen là các điểm đã cho của hàm, đường cong (x) được tính toán phải đi qua các điểm này. Hình 1.1. Mô tả phương pháp hồi qui chính xác tại các nút Trong thực tế, người ta thường tìm hàm (x) dưới dạng đa thức bậc n dạng y = a0+ a1x + a1x + a2x2 + ….+ anxn vì những lý do: - Một hàm bất kỳ có đạo hàm bậc n có thể khai chuyển thành chuỗi đa thức - Nội suy hàm đa thức sẽ đưa bài toán về việc giải hệ phương trình tuyến tính Giả sử hàm thực nghiệm cho dưới dạng bảng (như bảng 1.1), các giá trị y0, y1, y2 … yn tương ứng với các đối số x0, x1, x2, … xn (gọi là các nút nội suy). 6 Khi ấy, các ẩn số (các hệ số của đa thức) a0, a1, a2, … an của đa thức nội suy sẽ được xác định thông qua việc giải hệ phương trình đại số (1.1)[7]: (1.1) Về mặt toán học, hệ phương trình (1.1) có nghiệm duy nhất, nghĩa là tồn tại một đa thức bậc n, gọi là đa thức nội suy, đáp ứng điều kiện (xi ) = yi ( với i = 1, 2, 3…n). Điều này cũng có nghĩa là trong mặt phẳng (x,y) tồn tại một đường cong đi qua n+1 điểm và là biểu diễn của đa thức bậc n. Điều đó cũng có nghĩa là nếu số điểm thực nghiệm là 2 thì qua 2 điểm đó có duy nhất một đường thẳng y= a0 + a1x, nếu số điểm là 3 thì qua 3 điểm đó sẽ có parabol y= a0 + a1x + a2x2 . Một điểm thú vị là ta có thể nhận được đa thức nội suy mà không phải giải hệ (1.1). Đã có nhiều công thức có thể tìm các đa thức nội suy mà không phải giải hệ phương trình nói trên (ví dụ công thức Lagrange, Newton, Gauss, Bessel, Sterling …). Các công thức nội suy này được gọi là thông dụng theo ý là chúng được trình bày trong nhiều tài liệu giáo khoa. 1.1.1. Nội suy Lagrange Công thức Lagrange có dạng [1,2,4,7]: 7 (1.2) Ta có thể kiểm định xem (1.2) có phải biểu thức cần tìm không. Tại nút x=x0, tất cả các thành phần từ thứ 2 trở đi đều bằng không. Thành phần đầu tiên sau rút gọn chỉ còn lại y0, như vậy y(x0)=y0. Tại nút x=x1, trừ thành phần thứ 2, các thành phần khác bằng không, ta có y(x1)=y1. Tương tự với mọi nút, ta có thẻ khẳng định (1.1) là đa thức cần tìm thích hợp. Ví dụ 1.1: Cho bảng số liệu (bảng 1.2), nhiệm vụ là tìm đa thức nội suy và tính giá trị hàm tại x=3 và x=5. Bảng 1.2 Đưa các giá trị trong bảng 1.2 vào công thức (1.2) ta có: (1.3) Rút gọn (1.3) ta có y(x)= 1- 2x + 3x2 , tại x=3 và x=5 ta tính được y(3)= 22, y(5)=66. Thế mạnh của nội suy Lagrange là áp dụng được cho cả trường hớp nút cách đều và cách không đều. 1.1.2. Nội suy Newton a. Trường hợp nút cách đều 8 Công thức nội suy Newton cho trường hợp nút cách đều có dạng [1,2,7]: (1.4) Ở đây: xi - các nút nội suy với chỉ số i=0, 1, 2, 3 …(n-1) có bước cách đều, nghĩa là x1= x0 +h, x2=x1 +h, …, xn-1 = xn-2 =xn-1 +h; ci – các hệ số của công thức nội suy, với i = 0, 1, 2, 3 …n Từ công thức (1.4) ta thấy hàm y(x) của công thức nội suy này cũng là một đa thức bậc n , nhưng viết dưới dạng khác. Ta có thể xác định các hệ số ci : Với x=x0 , từ (1.4) ta thấy y(x0)=c0 hay c0 =y0 Với x=x1 ta có y(x1)= y1 = c0 + c1(x1-x0) hay c1=(y1-c0)/(x1-x0)= (y1-y0)/(x1x0). Ta đặt y1-y0 = Δy0 , x1-x0=h sẽ nhận được c1 = Δy0/h Với x=x2 thì y(x2)= y2 = c0 + c1(x2-x0) + c2(x2-x0)(x2-x1) , đưa giá trị c0, c1 vào ta có: Tiếp tục với các hệ số khác ta sẽ nhận được bieur thức tổng quát sau : (1.5) với k=1, 2, 3…n 9 Công thức cho hàm hồi qui Newton sẽ có dạng: (1.6) Công thức (1.6) có thể gọn hơn, đặt (x-x0)/h = t khi ấy: Và công thức (1.6) có trở thành: (1.7) Theo khảo sát, nội suy theo công thức (1.7) cho sai số lớn dần về phần cuối bảng số liệu, vì vậy, khi cần độ chính xác hơn ở đoạn cuối có thể sử dụng công thức [7]: (1.8) Tương tự như biến đổi công thức (1.4), ta sẽ có được các hệ số Ci : 10 Thay vào (1.8) ta có nhận được công thức (1.9): (1.9) Đặt (x- xn)/h = t và thay vào (1.9) ta có: (1.10) Công thức (1.9) và (1.10) được gọi là nội suy Newton lùi [7], công thức (1.6), (1.7) gọi là nội suy Newton tiến. Ta xem xét một ví dụ áp dụng nội suy Newton. Ví dụ 1.2: Số liệu giả sử được trình bày trong bảng 1.3 sau: Bảng 1.3 X 1 2 3 4 5 6 y 2 9 22 41 66 97 Δy Δ2 y Δ3 y 7 13 6 19 6 0 25 6 0 31 6 0 11 Bảng số liệu 1.3 thực chất là bảng 1.2 nhưng đầy đủ hơn, bài toán yêu cầu xác định đa thức nội suy theo công thức nội suy Newton tiến. Theo dữ liệu đã cho, ta có y0 = 2, h=1 và các số gia được tính (đã tính trong bảng) Δy0=7, Δ2y0 = 6, Δ3y = 0. Sử dụng công thức (1.6) ta có: Lắp ráp số liệu từ bảng 1.3 ta nhận được: Đơn giản hơn: Như vậy đa thức vừa tìm được này hoàn toàn trùng với đa thức chúng ta nhận được trong ví dụ 1.1 với nội suy Lagrange. b. Trường hợp nút không cách đều: Trường hợp các nút không cách đều (h ≠ const), công thức nội suy Newton được đề xuất như sau [1,2,7]: (1.11) Về mặt hình thức, công thức này giống với công thức trong trường hợp nút cách đều (1.4). Sự khác nhau là cách tính gia số σk y0 (với k =1, 2, 3…n). - Các gia số bậc 1 có dạng: - Các gia số bậc 2 có dạng: 12 - Dạng tổng quát : (1.12) Ta sẽ xem xét việc sử dụng công thức nội suy nút không cách đều (1.11) qua ví dụ dưới đây. Ví dụ 1.3. Giả sử số liệu đã cho trong bảng 1.2 (như ví dụ 1.1). Ta tính các gia số bậc 1 và 2. - Các gia số bậc 1: - Các gia số bậc 2: Các gia số bậc 2 bằng nhau, không có đa thức bậc cao hơn, đa thức cần tìm có bậc là 2. Đưa các gia số và số liệu trong bảng 1.2 vao (1.11) ta nhận được: Kết quả nội suy trùng với kết quả trong các ví dụ phía trên. 1.1.3. Nội suy Gauss Giả sử hàm thực nghiệm y(x) cho dưới dạng bảng có 2n+1 nút, các nút cách đều nhau (h = xi – xi-1 = const) : 13 Hàm nội suy Gauss cần tìm là đa thức với bậc không cao hơn 2n có dạng: (1.13) Ta phải xác định các hệ số ci của hàm (1.13): - Khi x= x0 ta có y(x0) = y0 = c0 Tiếp tục ta nhận được: Đặt (x – x0) /h = t và đưa biến t vào (1.13) ta có: (1.14) Công thức (1.14) được gọi là công thức nội suy Gauss thứ nhất. Cũng như nội suy Newton, còn có công thức nội suy Gauss thứ 2 như sau [7]: 14