A solution to the numerical approximation problem from MSU's olympiad

msu_numerical.py

from sympy import symbols, expand, Matrix, solve_linear_system, Rational, N

a, b, c = symbols('a b c')
h = symbols('h')
g1, g2, g3 = symbols('g1 g2 g3')

f = lambda x: a + b * x + c * x**2

approx = g1 * f(-h/6) + g2 * f(0) + g3 * f(5 * h / 6)

print(expand(approx))

#                                                2            2
#                      b⋅g₁⋅h   5⋅b⋅g₃⋅h   c⋅g₁⋅h    25⋅c⋅g₃⋅h
# a⋅g₁ + a⋅g₂ + a⋅g₃ - ────── + ──────── + ─────── + ──────────
#                        6         6          36         36

system = Matrix((
	(        1, 1,              1, 0),
	(     -h/6, 0,        5*h / 6, 0),
	(h**2 / 36, 0, 25 * h**2 / 36, 2)
))


sols = solve_linear_system(system, g1, g2, g3)

y = lambda x: x**3

ans = sols[g1] * y(-h/6) + sols[g2] * y(0) + sols[g3] * y(5 * h / 6)
ans = ans.subs(h, Rational(6, 100))
print(N(ans))
# 0.0800000000000000