| # module 'poly' -- Polynomials |
| |
| # A polynomial is represented by a list of coefficients, e.g., |
| # [1, 10, 5] represents 1*x**0 + 10*x**1 + 5*x**2 (or 1 + 10x + 5x**2). |
| # There is no way to suppress internal zeros; trailing zeros are |
| # taken out by normalize(). |
| |
| def normalize(p):# Strip unnecessary zero coefficients |
| n= len(p) |
| while p: |
| if p[n-1]:return p[:n] |
| n= n-1 |
| return[] |
| |
| def plus(a, b): |
| if len(a)< len(b): a, b= b, a# make sure a is the longest |
| res= a[:]# make a copy |
| for iin range(len(b)): |
| res[i]= res[i]+ b[i] |
| return normalize(res) |
| |
| def minus(a, b): |
| if len(a)< len(b): a, b= b, a# make sure a is the longest |
| res= a[:]# make a copy |
| for iin range(len(b)): |
| res[i]= res[i]- b[i] |
| return normalize(res) |
| |
| def one(power, coeff):# Representation of coeff * x**power |
| res=[] |
| for iin range(power): res.append(0) |
| return res+[coeff] |
| |
| def times(a, b): |
| res=[] |
| for iin range(len(a)): |
| for jin range(len(b)): |
| res= plus(res, one(i+j, a[i]*b[j])) |
| return res |
| |
| def power(a, n):# Raise polynomial a to the positive integral power n |
| if n==0:return[1] |
| if n==1:return a |
| if n/2*2== n: |
| b= power(a, n/2) |
| return times(b, b) |
| return times(power(a, n-1), a) |
| |
| def der(a):# First derivative |
| res= a[1:] |
| for iin range(len(res)): |
| res[i]= res[i]*(i+1) |
| return res |
| |
| # Computing a primitive function would require rational arithmetic... |
