Search code examples
pythonnumpyscipynumerical-integrationpolynomials

How do you integrate over a curve in 2D space that is parameterized by a polynomial?


I have a simple curve in 2D space (x, y) that is parameterized by a polynomial as a function of t and I would like to know the length of that curve. How do I accomplish that? I looked into scipy.integrate and numpy.polyint but I didn't manage to find a solution. It seems that both of them are only able to integrate over a 1D polynomial. Here is an example for a curve:

import numpy as np
from scipy import integrate


x0, y0 = 0.0, 0.0
vx, vy = 0.1, 0.1
ax, ay = -0.0001, 0
coeff = np.array([[ax, ay], [vx, vy], [x0, y0]])

pos = lambda t: np.polyval(coeff, t)

Solution

  • The arc length is a single-variable polynomial on the curve parameter. You need to define the expression for the differential of the arc length and then you'll be able to integrate over it, as explained in the link in the comments. As you can see there, it can be simply expressed as the Euclidean norm of the vector (dx/dt, dy/dt). Here is therefore how you can implement it:

    import numpy as np
    import scipy
    
    x0, y0 = 0.0, 0.0
    vx, vy = 0.1, 0.1
    ax, ay = -0.0001, 0
    coeff = np.array([[ax, ay], [vx, vy], [x0, y0]])
    
    # Position expression is not really necessary
    pos = lambda t: np.polyval(coeff, t)
    
    # Derivative of the arc length
    def ds(t):
        # Coefficients of polynomial derivative
        coeff_d = coeff[:-1] * np.arange(len(coeff) - 1, 0, -1)[:, np.newaxis]
        # Norm of position derivatives
        return np.linalg.norm(np.polyval(coeff_d, np.expand_dims(t, -1)), axis=-1)
    
    # Integrate across parameter interval
    t_start, t_end = 0, 1
    arc_length, err = scipy.integrate.quad(ds, t_start, t_end)
    print(arc_length)
    # 0.1413506691471052
    

    Of course, you could try to work out the analytical expression of the integral of ds and then you wouldn't need any integration method.