Sebastian Cor
Sebastian Cor

Reputation: 9

Polynomial division algorithm

I'm using sage and was trying to implement univariate polynomial division with the pseudocode given by Wikipedia.

But I think it is stuck looping, for example if I ask div(x^2-1,x-1) it doesn't give the immediate answer. It should return (0,x+1) but it does nothing.

The code:

def div(p,q):   
    if q==0:
        return("NaN")
    elif q!=0:
        l=0
        r=p
        while r!=0 and q.degree()<=r.degree():
            t=(r.leading_coefficient())/(q.leading_coefficient())
            l=l+t
            r=r-(t*q)
        return(l,r)

Edit: I was reading the pseudocode wrong and missed that I was not reducing the degree of my polynomial so clearly it just was doing nothing. I 'fixed it' but now it is giving me some new error but I think it's some coercion error.

Any help is appreciated!

The new code:

def div(p,q):
    if q==0:
        return("NaN")
    elif q!=0:
        l=0
        r=p
        while r!=0 and q.degree()<=r.degree():
            t=r.leading_coefficient()/q.leading_coefficient()
            m=x^r.degree()/x^q.degree()
            l=l+t*m
            r=r-(t*m*q)
            print(l,r) #To see when the code fails
        return(l,r)   

Edit 2: Upon inspection of "Polynomials in sage" it says that if you divide two polynomials then it coerces it to an element of the fraction field which would give me an error in the r.degree() line. Anyone knows a workaround to this?

Upvotes: 0

Views: 1403

Answers (2)

Samuel Leli&#232;vre
Samuel Leli&#232;vre

Reputation: 3453

In Sage

  • p / q gives a result in the fraction field
  • p // q gives the quotient in the ring (dropping any remainder)
  • p % q gives the remainder
  • p.quo_rem(q) gives the pair (quotient, remainder) at once

Upvotes: 1

Sebastian Cor
Sebastian Cor

Reputation: 9

Fixed by casting m into the base ring, m=R(m) I guess this is some sloppy patch I would like to know if there is some smarter way to make this algorithm. The fixed code:

x=var('x')
R.<x>=PolynomialRing(QQ)

def div(p,q):
    if q==0:
        return("NaN")
    elif q!=0:
        l=0
        r=p
        while r!=0 and q.degree()<=r.degree():
            t=r.leading_coefficient()/q.leading_coefficient()
            m=x^r.degree()/x^q.degree()
            m=R(m)
            l=l+t*m
            r=r-(t*m*q)
            print(l,r)
        return(l,r)

Upvotes: 0

Related Questions