二次插值法python实现

2019独角兽企业重金招聘Python工程师标准>>>

# -*- coding: utf-8 -*-
'''
二次插值法python实现
f(x)=x^4 - 4x^3 - 6x^2 -16x +4极值
区间[-1,6] e=0.05
'''
import numpy as np
import matplotlib.pyplot as plt

'''
函数表达式
'''
def f(x):
    return 1.0 * (pow(x, 4) - 4 * pow(x, 3) - 6 * pow(x, 2) - 16 * x + 4)

'''
绘制函数图像
'''
def printFunc(f, a, b, x, y):
    t = np.arange(a, b, 0.01)
    s = f(t)
    plt.plot(t, s)
    plt.plot([x], [y], 'ro')
    plt.show()

'''
e为精度
'''
def search(f, x1, x2, x3, e):

    if f(x2) > f(x1) or f(x2) > f(x3):
        print "不满足两头大中间小的性质"
        return
    
    #系数矩阵
    A = [[pow(x1, 2), x1, 1], [pow(x2, 2), x2, 1], [pow(x3, 2), x3, 1]]
    b = [f(x1), f(x2), f(x3)]

    X= np.linalg.solve(A, b)

    a0 = X[0]
    a1 = X[1]
    a2 = X[2]

    x = (-1)*a1/(2 * a0)
    print x
    if abs(x - x2) < e:
        if f(x) < f(x2):
            y = f(x)
            print x
            return (x, y)
        else:
            y = f(x2)
            print x2
            return (x2, y)

    arr = [x1, x2, x3, x]
    arr.sort()

    if f(x2)> f(x):
        index = arr.index(x)
        x2 = x
        x1 = arr[index - 1]
        x3 = arr[index + 1]
    else:
        index = arr.index(x2)
        
        x1 = arr[index - 1]
        x3 = arr[index + 1]
        
    return search(f, x1, x2, x3, e)
        
def regre(f, a, b, e):
    x1 = a
    x3 = b
    x2 = (a + b)/2.0
    p = search(f, x1, x2, x3, e)
    printFunc(f, a, b, p[0], p[1])

a = -1.0
b = 6.0
e = 0.05
regre(f, a, b, e)

    

 

转载于:https://my.oschina.net/wtzheng/blog/745019