var("x0 x1 p0 p1 dp0 dp1") var("a0 a1 a2 a3") p(x) = a0 + a1*x + a2*x^2 + a3*x^3 dp(x) = derivative(p(x), x) eq0 = (p(x = x0) == p0) eq1 = (p(x = x1) == p1) eq2 = (dp(x = x0) == dp0) eq3 = (dp(x = x1) == dp1) res = solve([eq0, eq1, eq2, eq3], [a0, a1, a2, a3]) coeff0 = res[0][0].rhs() coeff1 = res[0][1].rhs() coeff2 = res[0][2].rhs() coeff3 = res[0][3].rhs() def cubInterpol(xx0, xx1, y0, y1, dy0, dy1): return ( coeff0(x0 = xx0, x1 = xx1, p0 = y0, p1 = y1, dp0 = dy0, dp1 = dy1) + coeff1(x0 = xx0, x1 = xx1, p0 = y0, p1 = y1, dp0 = dy0, dp1 = dy1)*x + coeff2(x0 = xx0, x1 = xx1, p0 = y0, p1 = y1, dp0 = dy0, dp1 = dy1)*x^2 + coeff3(x0 = xx0, x1 = xx1, p0 = y0, p1 = y1, dp0 = dy0, dp1 = dy1)*x^3 ) def plotCubInterpol( xx0, xx1, y0, y1, dy0, dy1, xmin = -8, xmax = 8, ymin = -6, ymax = 6 ): cubPol = cubInterpol(xx0, xx1, y0, y1, dy0, dy1) return plot( cubPol, x, xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax, aspect_ratio = 1, color = "blue" ) def drawCubicInterpolation( xx0, xx1, y0, y1, dy0, dy1, xmin = -8, xmax = 8, ymin = -6, ymax = 6 ): f = plotCubInterpol( xx0, xx1, y0, y1, dy0, dy1, xmin, xmax, ymin, ymax ); nodes = ( point([xx0, y0], color="red", size=40) + point([xx1, y1], color="red", size=40) ) tangents = ( line([(xx0 - 1, y0 - dy0), (xx0 + 1, y0 + dy0)], color="green") + line([(xx1 - 1, y1 - dy1), (xx1 + 1, y1 + dy1)], color="green") ) return f + nodes + tangents