数值计算讨论组

欢迎来到数值计算讨论组。 本论坛旨在增强师生交流的时效性,主要讨论数值计算相关问题,包括初等数值分析、矩阵计算、非线性方程(组)数值解法、微分方程数值解法等方面的内容。


    [数]Newton插值多项式相关两个参考程序

    分享

    Admin
    Admin

    帖子数 : 46
    注册日期 : 10-05-13

    [数]Newton插值多项式相关两个参考程序

    帖子 由 Admin 于 周三 十月 27, 2010 9:25 pm

    课本上给出的那种方法请参考帖子“(Newton插值多项式——cos x)帮忙看下这程序为啥不对好吗?谢谢老师和同学”,不过下面的程序在计算差商或差分时更节省存储空间。

    1 cos(x) (未给出误差分析,请参考插值多项式的误差余项个人完成)
    代码:
    % P37 例3  函数cos(x)的Newton插值逼近
    clear; clc; close all;
    format long;

    X=[0.0  0.2  0.4  0.6  0.8  1.0  1.2];
    Y=[1.0  0.980067  0.921061  0.825336  0.696707    0.540302  0.362358];
    n = length(X);
    A = Y';

    % 计算cos(x)的各阶差商,存入A
    % 在X(1),X(2),..,X(k)处的k-1阶差商存入A(k), k=1,2,...,n
    for i = 1:n-1
        for j = n:-1:i+1
            A(j) = ( A(j) - A(j-1))/(X(j) - X(j-i));
        end
    end
    A

    % 计算Newton插值多项式在x0处的函数值y0
    x0 = 0.3;
    y0 = A(1);
    temp = 1;
    for j = 2:n
        temp = temp * (x0 - X(j-1));
        y0 = y0 + A(j)*temp;
    end

    % 输出
    x0
    y0

    2 Bessel函数
    代码:
    % P44 例4 Bessel函数的Newton向前插值逼近
    clear;clc;
    format long
    x = 2.4:0.1:2.9;
    y = besselj(0,x);
    n = length(x);

    % 计算各阶差分。通过下面的计算,将x0=2.4处的第i阶差分存入y(i+1)
    for i = 1:n-1
        for j = n:-1:i+1
            y(j) = y(j)-y(j-1);
        end
    end

    % 通过Newton向前插值公式计算2.45处的函数值
    f = y(1);
    t = 0.5; % 0<t<1, 步长h = 0.1, t=0.5即现在考虑的点为2.45 + t*h = 2.45
    co = 1;
    for i = 1:n-1 % i 表示差分阶数
        % co 表示组合数(t,i)   
        co = co * (t-i+1)/i;
        % 输出 利用从点2.4及后面i个节点得到插值多项式 在2.45处的函数值
        f = f + co*y(i+1)
    end

    % 2.45处精确函数值
    f_exact = besselj(0,2.45)

      目前的日期/时间是周五 十二月 14, 2018 2:31 am