Python学习手札04——Newton插值计算问题

Newton插值计算

说明：该博的目标是求单点插值

计算步骤

1、按照需求定插值函数的阶数；
2、根据阶数选取阶数+1个已知数作为计算需要的数；
3、按照如下公式计算差商，最好列出差商表，方便计算；

4、然后将差商带入如下多项式，得出插值多项式；

5、最后带值计算出插值结果

代码示例

# -*- coding: utf-8 -*-
"""
Created on Fri Oct 11 12:52:34 2019

@author: ZHIFU
"""

import xlrd
def Newton(number, n):
    # number为x的值，n 的值为阶数+1
    # cols_1 = []
    # cols_2 = []
    workbook = xlrd.open_workbook(r"C:\Users\ZHIFU\Desktop\shuzhifenxi.xlsx")
    print(workbook)
    sheet_name =workbook.sheet_names()[0] #  获取sheet表格的名称
    sheet = workbook.sheet_by_index(0) #  根据sheet索引或者名称获取sheet内容
    print(sheet.name, sheet.nrows, sheet.ncols) #  打印名称及行列
    sheet_rows = sheet.nrows
    print("已知参数的个数：%d"% sheet_rows)
    cols_1 = sheet.col_values(0) #  获取第1列数据
    cols_2 = sheet.col_values(1) #  获取第2列数据
    print("已知参数自变量的值：%s"% cols_1)
    print("已知参数应变量的值：%s"% cols_2)
    print("从excle中选取的表格：%s"% sheet_name)
    freelist_1 = [] #  number与已知自变量差值的绝对值
    for i in range(0, sheet_rows):
        freelist_1.append(abs(cols_1[i] - number))
        
    # print("计算值与自变量的差的绝对值：%s"% freelist_1)
    freelist_1.sort()   
    freelist_1_len = len(freelist_1)
    print("freelist_1中元素的个数：%d"% freelist_1_len)
    freelist_index_1 = []  #  寻找最接近计算值的自变量在原列表（cols_1）中的索引值
    for i in range(0, n):
        for j in range(0, sheet_rows):
            if abs(cols_1[j] - number) == freelist_1[i]:
                freelist_index_1.append(j)
                freelist_index_2 = []
                for k in freelist_index_1:
                    if not k in freelist_index_2:
                        freelist_index_2.append(k)
        i += 1
    print("最接近计算值的自变量的索引：%s"% freelist_index_2)
    freelist_index_2.sort()
    freelist_x_value = [] 
    freelist_y_value = [] 
    for i in range(0, n):
        j = freelist_index_2[i]
        freelist_x_value.append(cols_1[j])
        freelist_y_value.append(cols_2[j])
    count = len(freelist_x_value)
    print("最接近计算值的自变量：%s"% freelist_x_value)
    print("最接近计算值的自变量对应的应变量：%s"% freelist_y_value)
    freelist_f = [] # 每一阶差商的第一个差商按阶增加的顺序添加到列表中
    fist_x_y = freelist_y_value[0] # 第一个x对应的f(x)
    m = 0
    while m <= count:
        f1_list = []
        if len(freelist_y_value) == 1:
            break
        for j in range(0, count-m-1):
            for k in range(j+1, j+2):
                for i in range(j+m+1, j+m+2):
                    f1= ((freelist_y_value[k] - freelist_y_value[j])
                        /(freelist_x_value[i] - freelist_x_value[j]))
                    f1_list.append(f1)
        freelist_y_value = f1_list # 使每次得到的差商进入循环
        # print(freelist_y_value)
        freelist_f.append(freelist_y_value[0])
        m += 1
    freelist_f.insert(0, fist_x_y) # 在列表中插入0阶差商
    print("freelist_f的值为%s"% freelist_f) # 列表中除第一个元素为第一个x对应的f(x)外，其余均为差商
    r = 1
    r_list = [r,] # 0阶差商的系数为1
    for i in range(0, len(freelist_f)-1):
        r *= number - freelist_x_value[i]
        r_list.append(r)
    print("r_list对应的值为%s"% r_list) # 列表中除1以外，其余均为x - xi的乘积
    re = 0
    for i in range(0, len(freelist_f)):
        re += freelist_f[i] * r_list[i]
    print("%.2f 的 %d 阶插值结果为 %.4f "% (number, n-1, re))
Newton(0.285, 3)