如何解决Python中的函数逼近任务?

问题描述 投票:1回答:2

考虑线[1,15]上的复数数学函数:f(x)= sin(x / 5)* exp(x / 10)+ 5 * exp(-x / 2)

enter image description here

度数n的多项式(w_0 + w_1 x + w_2 x ^ 2 + ... + w_n x ^ n)由它通过的任何n + 1个不同点唯一定义。这意味着其系数w_0,... w_n可以从以下线性方程组确定:

enter image description here

其中x_1,...,x_n,x_ {n + 1}是多项式通过的点,并且通过f(x_1),...,f(x_n),f(x_ {n + 1}) - 它必须在这些点上取得的价值。

我正在尝试为第三度的多项式形成线性方程组(即,指定系数矩阵A和自由向量b),其必须与点1,4,10处的函数f一致,并且15.使用scipy.linalg.solve函数解决此系统。

A = numpy.array([[1,1。,1.,1。],[1.,4.,8.,64。],[1.,10.,100.,1000。],[ 1.,15.,225。,3375。]])

V = numpy.array([3.25,1.74,2.50,0.63])

numpy.linalg.solve(A, V)

我得到了错误的答案,这就是enter image description here

所以问题是:矩阵是否正确?

python algorithm matlab numpy math
2个回答
3
投票

不,你的矩阵不正确。

最大的错误是A的第二个子矩阵。第三个条目应该是4**2,这是16,但你有8.不太重要,你的常数数组V只有两个小数位,但你真的应该有更多的精度。线性方程组有时对提供的值非常敏感,因此尽可能精确。此外,最后三个条目的舍入很糟糕:你向下舍入但你应该四舍五入。如果你真的想要两个小数位(我不推荐),那么值应该是

V = numpy.array([3.25, 1.75, 2.51, 0.64])

但更好的是

V = numpy.array([3.252216865271419, 1.7468459495903677,
                 2.5054164070002463, 0.6352214195786656])

随着AV的变化,我得到了结果

array([ 4.36264154, -1.29552587,  0.19333685, -0.00823565])

我得到这两个同情图,第一个显示你的原始函数,第二个使用近似的三次多项式。

enter image description here

enter image description here

他们看起来很接近我!当我在1,4,10和15处计算函数值时,最大绝对误差为15,即-4.57042132584462e-6。这比我预期的要大一些,但可能已经足够了。


0
投票

它来自数据科学课程吗? :)这是我做的几乎通用的解决方案:

%matplotlib inline

import numpy as np;
import math;
import matplotlib.pyplot as plt;

def f(x):
    return np.sin(x / 5) * np.exp(x / 10) + 5 * np.exp(-x / 2)

# approximate at the given points (feel free to experiment: change/add/remove)
points = np.array([1, 4, 10, 15])
n = points.size

# fill A-matrix, each row is 1 or xi^0, xi^1, xi^2, xi^3 .. xi^n
A = np.zeros((n, n))
for index in range(0, n):
    A[index] = np.power(np.full(n, points[index]), np.arange(0, n, 1))

# fill b-matrix, i.e. function value at the given points
b = f(points)

# solve to get approximation polynomial coefficents
solve = np.linalg.solve(A,b)

# define the polynome approximation of the function
def polinom(x): 
    # Yi = solve * Xi where Xi = x^i
    tiles = np.tile(x, (n, 1))
    tiles[0] = np.ones(x.size)
    for index in range(1, n):
        tiles[index] = tiles[index]**index
    return solve.dot(tiles)

# plot the graphs of original function and its approximation
x = np.linspace(1, 15, 100)
plt.plot(x, f(x))
plt.plot(x, polinom(x))

# print out the coefficients of polynome approximating our function
print(solve)
© www.soinside.com 2019 - 2024. All rights reserved.