满足积分和两点的多项式
Polynomial which satisfies integral and two points
Consider two points (x_0, f_0) and (x_1, f_1)
let p(x) be the degree two polynomial for which
p(x_0) = f_0
p(x_1) = f_1
and the integral of p(x) from -1 to 1 is equal to 0
Write a function which accepts two arguments
1. a length 2 NumPy vector 'x' of floating point values, with 'x[i]' containing the value of x_i,
2. a length 2 NumPy vector 'f' of floating point values, with 'f[i]' containing the value of f_i,
and which returns
a length 3 NumPy vector of floating point values containing the power series coefficients, in order from the highest order term to the constant term, for p(x)
我不知道从哪里开始。我最初的想法是有一个微分方程 P(1)=P(-1) 初始值 p(x_0) = f_0 和 p(x_1) = f_1,但我也遇到了实施问题。
使用sympy、Python的符号数学库,问题可以表述如下:
from sympy import symbols Eq, solve, integrate
def give_coeff(x, f):
a, b, c, X = symbols('a, b, c, X')
F = a * X * X + b * X + c # we have a second order polynomial
sol = solve([Eq(integrate(F, (X, -1, 1)), 0), # the integral should be zero (2/3*a + 2*c)
Eq(F.subs(X, x[0]), f[0]), # filling in x[0] should give f[0]
Eq(F.subs(X, x[1]), f[1])], # filling in x[1] should give f[1]
(a, b, c)) # solve for a, b and c
return sol[a].evalf(), sol[b].evalf(), sol[c].evalf()
import numpy as np
coeff = give_coeff(np.array([1, 2]), np.array([3, 4]))
print(coeff)
代码甚至可以展开为任意次数的多项式:
from sympy import Eq, solve, symbols, integrate
def give_coeff(x, f):
assert len(x) == len(f), "x and f need to have the same length"
degree = len(x)
X = symbols('X')
a = [symbols(f'a_{i}') for i in range(degree + 1)]
F = 0
for ai in a[::-1]:
F = F * X + ai
sol = solve([Eq(integrate(F, (X, -1, 1)), 0)] +
[Eq(F.subs(X, xi), fi) for xi, fi in zip(x, f)],
(*a,))
# print(sol)
# print(F.subs(sol).expand())
return [sol[ai].evalf() for ai in a[::-1]]
import numpy as np
coeff = give_coeff(np.array([1, 2]), np.array([3, 4]))
print(coeff)
print(give_coeff(np.array([1, 2, 3, 4, 5]), np.array([3, 4, 6, 9, 1])))
PS:只用numpy求解二次方程,np.linalg.solve
可以用3个方程求解3个未知数的线性方程组。方程式需要 "hand calculated",这更容易出错,并且更详细地扩展到更高的程度。
import numpy as np
def np_give_coeff(x, f):
# general equation: F = a*X**2 + b*X + c
# 3 equations:
# integral (F, (X, -1, 1)) == 0 or (2/3*a + 2*c) == 0
# a*x[0]**2 + b*x[0] + c == f[0]
# a*x[1]**2 + b*x[1] + c == f[1]
A = np.array([[2/3, 0, 2],
[x[0]**2, x[0], 1],
[x[1]**2, x[1], 1]])
B = np.array([0, f[0], f[1]])
return np.linalg.solve(A, B)
coeff = np_give_coeff(np.array([1, 2]), np.array([3, 4]))
print(coeff)
你可以一般地解决这个问题,利用
并将其添加为约束条件。然后你有 3 个方程式 3 个未知数 (a, b, c)。
还有其他有趣的技巧,这是一个巧妙的问题。尝试根据 a(x-b)(x-c)
编写公式,然后得到 3bc + 1 = 0
。任何以点 (x0,y0),(x1,x1)
开头的解决方案对于 (k*x0,k*y0),(k*x1,k*y1)
都有类似的解决方案。
Consider two points (x_0, f_0) and (x_1, f_1)
let p(x) be the degree two polynomial for which
p(x_0) = f_0
p(x_1) = f_1
and the integral of p(x) from -1 to 1 is equal to 0
Write a function which accepts two arguments
1. a length 2 NumPy vector 'x' of floating point values, with 'x[i]' containing the value of x_i,
2. a length 2 NumPy vector 'f' of floating point values, with 'f[i]' containing the value of f_i,
and which returns
a length 3 NumPy vector of floating point values containing the power series coefficients, in order from the highest order term to the constant term, for p(x)
我不知道从哪里开始。我最初的想法是有一个微分方程 P(1)=P(-1) 初始值 p(x_0) = f_0 和 p(x_1) = f_1,但我也遇到了实施问题。
使用sympy、Python的符号数学库,问题可以表述如下:
from sympy import symbols Eq, solve, integrate
def give_coeff(x, f):
a, b, c, X = symbols('a, b, c, X')
F = a * X * X + b * X + c # we have a second order polynomial
sol = solve([Eq(integrate(F, (X, -1, 1)), 0), # the integral should be zero (2/3*a + 2*c)
Eq(F.subs(X, x[0]), f[0]), # filling in x[0] should give f[0]
Eq(F.subs(X, x[1]), f[1])], # filling in x[1] should give f[1]
(a, b, c)) # solve for a, b and c
return sol[a].evalf(), sol[b].evalf(), sol[c].evalf()
import numpy as np
coeff = give_coeff(np.array([1, 2]), np.array([3, 4]))
print(coeff)
代码甚至可以展开为任意次数的多项式:
from sympy import Eq, solve, symbols, integrate
def give_coeff(x, f):
assert len(x) == len(f), "x and f need to have the same length"
degree = len(x)
X = symbols('X')
a = [symbols(f'a_{i}') for i in range(degree + 1)]
F = 0
for ai in a[::-1]:
F = F * X + ai
sol = solve([Eq(integrate(F, (X, -1, 1)), 0)] +
[Eq(F.subs(X, xi), fi) for xi, fi in zip(x, f)],
(*a,))
# print(sol)
# print(F.subs(sol).expand())
return [sol[ai].evalf() for ai in a[::-1]]
import numpy as np
coeff = give_coeff(np.array([1, 2]), np.array([3, 4]))
print(coeff)
print(give_coeff(np.array([1, 2, 3, 4, 5]), np.array([3, 4, 6, 9, 1])))
PS:只用numpy求解二次方程,np.linalg.solve
可以用3个方程求解3个未知数的线性方程组。方程式需要 "hand calculated",这更容易出错,并且更详细地扩展到更高的程度。
import numpy as np
def np_give_coeff(x, f):
# general equation: F = a*X**2 + b*X + c
# 3 equations:
# integral (F, (X, -1, 1)) == 0 or (2/3*a + 2*c) == 0
# a*x[0]**2 + b*x[0] + c == f[0]
# a*x[1]**2 + b*x[1] + c == f[1]
A = np.array([[2/3, 0, 2],
[x[0]**2, x[0], 1],
[x[1]**2, x[1], 1]])
B = np.array([0, f[0], f[1]])
return np.linalg.solve(A, B)
coeff = np_give_coeff(np.array([1, 2]), np.array([3, 4]))
print(coeff)
你可以一般地解决这个问题,利用
并将其添加为约束条件。然后你有 3 个方程式 3 个未知数 (a, b, c)。
还有其他有趣的技巧,这是一个巧妙的问题。尝试根据 a(x-b)(x-c)
编写公式,然后得到 3bc + 1 = 0
。任何以点 (x0,y0),(x1,x1)
开头的解决方案对于 (k*x0,k*y0),(k*x1,k*y1)
都有类似的解决方案。