Add the divided_difference_coefficients() function.

diff --git a/mjo/interpolation.py b/mjo/interpolation.py
index 0cafc283e7b445f0ec8af9487f3964195a291cc5..d5011ea47e154048b79396f4141fb907597bc614 100644 (file)
--- a/mjo/interpolation.py
+++ b/mjo/interpolation.py
@@ -73,6 +73,26 @@ def lagrange_polynomial(x, xs, ys):
     return sigma
 
 
+
+def divided_difference_coefficients(xs):
+    """
+    Assuming some function `f`, compute the coefficients of the
+    divided difference f[xs[0], ..., xs[n]].
+
+    TESTS:
+
+        sage: divided_difference_coefficients([0])
+        [1]
+        sage: divided_difference_coefficients([0, pi])
+        [-1/pi, 1/pi]
+        sage: divided_difference_coefficients([0, pi, 2*pi])
+        [1/2/pi^2, -1/pi^2, 1/2/pi^2]
+
+    """
+    coeffs = [ product([ (QQ(1) / (xj - xi)) for xi in xs if xi != xj ])
+               for xj in xs ]
+    return coeffs
+
 def divided_difference(f, xs):
     """
     Return the Newton divided difference of `f` at the points