Bài giảng Applied numerical methods (Ứng dụng phương pháp tính số): Chương 4 - TS. Ngô Văn Thanh

Chia sẻ: Little Little | Ngày: | Loại File: PDF | Số trang:24

Thêm vào BST

Báo xấu

94
lượt xem 14
download

Download Vui lòng tải xuống để xem tài liệu đầy đủ

Chương 4 - Các phương pháp cho hàm rời rạc. Chương này trình bày về đa thức nội suy (Polynomial interpolation) và phương pháp bình phương tối thiểu. Nội dung cụ thể gồm có: Nội suy Lagrangian, nội suy Newton, nội suy bậc 3 (Cubic spline interpolation), làm khớp dữ liệu cho các hàm tuyến tính, làm khớp dữ liệu cho các hàm phi tuyến.

Chủ đề:

Bình luận(0) Đăng nhập để gửi bình luận!

Lưu

Nội dung Text: Bài giảng Applied numerical methods (Ứng dụng phương pháp tính số): Chương 4 - TS. Ngô Văn Thanh

TS. Ngô Văn Thanh, Viện Vật lý. Cao học vật lý – chuyên ngành Vật lý lý thuyết.
Chương 4. Các phương pháp cho hàm rời rạc. 4.1. Đa thức nội suy (Polynomial interpolation) 4.1.1. Nội suy Lagrangian. 4.1.2. Nội suy Newton. 4.1.3. Nội suy bậc 3 (Cubic spline interpolation). 4.2. Phương pháp bình phương tối thiểu. 4.2.1. Làm khớp dữ liệu cho các hàm tuyến tính. 4.2.2. Làm khớp dữ liệu cho các hàm phi tuyến. @2009, Ngô Văn Thanh - Viện Vật Lý
4.1. Đa thức nội suy (Polynomial interpolation). 4.1.1. Nội suy Lagrangian. ¾ Xét N điểm dữ liệu (xi, yi), i =1, N: Đa thức nội suy bậc N - 1 có dạng: ¾ Đa thức nội suy chỉ có thể tính được các giá trị của @2009, Ngô Văn Thanh - Viện Vật Lý
@2009, Ngô Văn Thanh - Viện Vật Lý
¾ Thuật toán Neville : P1, P2, …, PN là đa thức bậc 0 qua mỗi một điểm P12, P23, …, P(N-1)N là đa thức bậc 1 qua từng cặp điểm Đa thức bậc N – 1 : P123…N Ví dụ với N = 4: xây dựng dạng bảng các đa thức nội suy x1: x2: x3: x4: Biểu thức liên hệ giữa hai đa thức bậc thấp hơn gần nhất: @2009, Ngô Văn Thanh - Viện Vật Lý
¾ Thuật toán Neville (tiếp theo) : Định nghĩa : @2009, Ngô Văn Thanh - Viện Vật Lý
¾ Chương trình: ********************************************************* INTEGER,FUNCTION iminloc(arr) REAL, INTENT(IN) :: arr INTEGER, DIMENSION(1) :: imin imin=minloc(arr(:)) iminloc=imin(1) END FUNCTION iminloc ********************************************************* SUBROUTINE polint(xa,ya,x,y,dy) IMPLICIT NONE REAL, DIMENSION(:), INTENT(IN) :: xa,ya REAL, INTENT(IN) :: x REAL, INTENT(OUT) :: y,dy INTEGER :: m,n,ns REAL, DIMENSION(size(xa)) :: c,d,den,ho n=size(xa) @2009, Ngô Văn Thanh - Viện Vật Lý
c=ya !Initialize the tableau of c's and d's. d=ya ho=xa-x ns=iminloc(abs(x-xa)) !Find index ns of closest table entry. y=ya(ns) !This is the initial approximation to y. ns=ns-1 do m=1,n-1 !For each column of the tableau, den(1:n-m)=ho(1:n-m)-ho(1+m:n) if(any(den(1:n-m)==0.0)) write(6,*)'calculation failure' den(1:n-m)=(c(2:n-m+1)-d(1:n-m))/den(1:n-m) d(1:n-m)=ho(1+m:n)*den(1:n-m) c(1:n-m)=ho(1:n-m)*den(1:n-m) if (2*ns < n-m) then dy=c(ns+1) else dy=d(ns) ns=ns-1 endif y=y+dy enddo END SUBROUTINE polint @2009, Ngô Văn Thanh - Viện Vật Lý
4.1.2. Nội suy Newton ¾ Biểu thức Newton ¾ Xác định hệ số c người ta sử dụng biểu thức: ¾ Dạng tính kiểu truy hồi ¾ Dạng tổng quát: @2009, Ngô Văn Thanh - Viện Vật Lý
4.1.3. Cubic spline interpolation ¾ Xét N điểm dữ liệu (xi, yi), i =1, N. ¾ Xét điểm x nằm trong khoảng hai điểm ¾ Ta có: ¾ Nếu ta có các giá trị đạo hàm bậc hai thì ta có thể viết ¾ Với ¾ Các đạo hàm: @2009, Ngô Văn Thanh - Viện Vật Lý
¾ Chương trình: ********************************************************* SUBROUTINE spline(x,y,yp1,ypn,y2) IMPLICIT NONE REAL, DIMENSION(:), INTENT(IN) :: x,y REAL, INTENT(IN) :: yp1,ypn REAL, DIMENSION(:), INTENT(OUT) :: y2 INTEGER :: n REAL, DIMENSION(size(x)) :: a,b,c,r n=size(x) c(1:n-1)=x(2:n)-x(1:n-1) !Set up the tridiagonal equations. r(1:n-1)=6.0*((y(2:n)-y(1:n-1))/c(1:n-1)) r(2:n-1)=r(2:n-1)-r(1:n-2) a(2:n-1)=c(1:n-2) b(2:n-1)=2.0*(c(2:n-1)+a(2:n-1)) b(1)=1.0 b(n)=1.0 @2009, Ngô Văn Thanh - Viện Vật Lý
if (yp1 > 0.99e30) then r(1)=0.0 c(1)=0.0 else r(1)=(3.0/(x(2)-x(1)))*((y(2)-y(1))/(x(2)-x(1))-yp1) c(1)=0.5 end if if (ypn > 0.99e30) then r(n)=0.0 a(n)=0.0 else r(n)=(-3.0/(x(n)-x(n-1))) & *((y(n)-y(n-1))/(x(n)-x(n-1))-ypn) a(n)=0.5 end if call tridag(a(2:n),b(1:n),c(1:n-1),r(1:n),y2(1:n)) END SUBROUTINE spline @2009, Ngô Văn Thanh - Viện Vật Lý
SUBROUTINE tridag(a,b,c,r,u) … bet=b(1) if(bet==0.) write(6,*) 'Error 1 in tridag' allocate(gam(n)) u(1)=r(1)/bet do j=2,n gam(j)=c(j-1)/bet bet=b(j)-a(j-1)*gam(j) if(bet==0.) write(6,*) 'Error 2 in tridag' u(j)=(r(j)-a(j-1)*u(j-1))/bet enddo do j=n-1,1,-1 u(j)=u(j)-gam(j+1)*u(j+1) enddo deallocate(gam) @2009, Ngô Văn Thanh - Viện Vật Lý
FUNCTION splint(xa,ya,y2a,x) IMPLICIT NONE REAL, DIMENSION(:), INTENT(IN) :: xa,ya,y2a REAL, INTENT(IN) :: x REAL :: splint INTEGER :: khi,klo,n REAL :: a,b,h n=size(xa) klo=max(min(locate(xa,x),n-1),1) khi=klo+1 h=xa(khi)-xa(klo) if (h == 0.0) write(6,*) 'The xa must be distinct' a=(xa(khi)-x)/h b=(x-xa(klo))/h splint=a*ya(klo)+b*ya(khi) & +((a**3-a)*y2a(klo)+(b**3-b)*y2a(khi))*(h**2)/6.0 END FUNCTION splint @2009, Ngô Văn Thanh - Viện Vật Lý
4.2. Phương pháp bình phương tối thiểu. ¾ Xét N điểm dữ liệu (xi, yi), i =1, N với M tham số hiệu chỉnh a1, a2, … aM. ¾ Hàm làm khớp có dạng: ¾ Định nghĩa hàm bình phương tối thiểu ¾ Trong đó độ lệch chuẩn ¾ Cực tiểu của hàm bình phương tối thiểu @2009, Ngô Văn Thanh - Viện Vật Lý
4.2.1. Làm khớp dữ liệu cho các hàm tuyến tính. Xét N điểm dữ liệu (xi, yi), i =1, N: Làm khớp dữ liệu cho một hàm tuyến tính bậc nhất có dạng: ¾ Hàm bình phương tối thiểu ¾ Cực tiểu của hàm bình phương tối thiểu ¾ Định nghĩa các tổng @2009, Ngô Văn Thanh - Viện Vật Lý
¾ Nghiệm của hệ là ¾ Sai số: Định nghĩa Viết lại nghiệm a, b và các sai số: @2009, Ngô Văn Thanh - Viện Vật Lý
Xét N điểm dữ liệu (xi, yi), i =1, N: Làm khớp dữ liệu cho một hàm tuyến tính có dạng: Trong đó Xk(x) gọi là hàm cơ sở ¾ Hàm bình phương tối thiểu ¾ Định nghĩa ma trận design A ¾ Các hệ số làm khớp cũng là một ma trận a @2009, Ngô Văn Thanh - Viện Vật Lý
¾ Lấy đạo hàm của hàm bình phương tối thiểu ¾ Viết lại dạng phương trình ma trận Trong đó ¾ Phương trình ma trận có thể viết dưới dạng: @2009, Ngô Văn Thanh - Viện Vật Lý
¾ Các hệ số a ¾ Độ lệch: @2009, Ngô Văn Thanh - Viện Vật Lý