5

I am trying to implement the divide-and-conquer algorithm for polynomial multiplication. Here is the pseudocode given in the lecture notes: enter image description here

where A, B are lists of coefficients of each polynomial, n is the size of the problem (degree - 1) and a_l, b_l are indices of the coefficients of interest.

Here is my attempt at implementing it using Python3:

def poly_mult_dc_naive(A, B, n, a, b):
  n = int(n)
  a = int(a)
  b = int(b)
  C = [None] * int(2*n - 1)

  if n == 1:
    C[0] = A[a] * B[b]
    return C[0]

  C[0:n-1] = poly_mult_dc_naive(A, B, n//2, a, b)
  C[n:2*n-1] = poly_mult_dc_naive(A, B, n//2, a + (n // 2), b + (n // 2))

  W = poly_mult_dc_naive(A, B, n/2, a, b + (n // 2))
  V = poly_mult_dc_naive(A, B, n/2, a + n/2, b)
  C[n // 2:n + (n // 2) - 1] += W + V
  return C

However I'm getting strange results. For example let A = [1,2,3,4] B = [4,3,2,1] I get:

[4, None, 8, 3, 6, 12, None, 16, 9, 12, 2, None, 4, 1, 2, None, 8, 3, 4, None, None]

Correct answer is [4, 11, 20, 30, 20, 11, 4]

Could someone please point out where I've gone wrong and how it could be done?

Mr. T
  • 8,325
  • 9
  • 23
  • 44
  • Just a tip about legibility, why not assume n,a,b are integers, since you're using `// 2` floor division. Make the caller responsible for passing integer arguments. – smci Sep 16 '18 at 01:24

2 Answers2

1

Quick update: I think I've managed to debug my code my using a numpy array for C instead of a list. Here is the updated version:

import numpy as np

def poly_mult_dc_naive(A, B, n: int, a: int, b: int):
  C = np.zeros(2*n - 1, dtype=int) # here I changed it from list to np array
  if n == 1:
    C[0] = A[a] * B[b]
    return C[0]
  C[0:n-1] = poly_mult_dc_naive(A, B, n//2, a, b)
  C[n:2*n-1] = poly_mult_dc_naive(A, B, n//2, a + (n // 2), b + (n // 2))
  W = poly_mult_dc_naive(A, B, n//2, a, b + (n // 2))
  V = poly_mult_dc_naive(A, B, n/2, a + (n//2), b)
  C[(n // 2) : (3*n // 2) - 1] += W + V
  return C

BONUS QUESTION: Does anyone know of a better way I could keep the arguments n, a, and b as int types?

I just feel like having to write:

n = int(n)
a = int(a)
b = int(b)

might not be the most elegant way.

  • Tell us what was the bug? – smci Sep 16 '18 at 01:25
  • I don't see that you need to coerce n,a,b from float to int, just use `// 2` rounding division everywhere (i.e. in the W,V lines), it keeps ints as int. If you have a counterexample which goes wrong, show us. – smci Sep 16 '18 at 01:27
  • 1
    Also, the line `C[n // 2:n + (n // 2) - 1]` badly needs parenthesizing, it misreads very easily. I'd write `C[(n//2) : (3*n//2)-1]` – smci Sep 16 '18 at 01:31
  • 2
    1. So the bug was that using a list for C (i.e. the result coefficients) leads to issues when assigning values to slices of it. So it would instead append values instead of assiging them. Changing the list to a numpy array seems to have fixed it. – Alaa A. Latif Sep 16 '18 at 01:45
  • 1
    2. You're right! Thank you! I tried that at first but I still got an error, but it turned out I missed doing // instead of / for a couple of lines. 3. I've parenthesized the line as suggested. Thanks for the tip! :) – Alaa A. Latif Sep 16 '18 at 01:47
  • Also for clarity, you don't need to initialize or compute `C[]` in the base-case (n==1), just directly `return A[a] * B[b]` (assuming that's correct). And strictly you don't need to initialize `C[]` ever, it's directly computed. So that eliminates two more lines. – smci Sep 16 '18 at 10:02
0
  • There should be no need to coerce n,a,b from float to int. Just use // 2 integer division everywhere (i.e. in the W,V lines). That will "keep" ints as ints.
  • The line C[n // 2:n + (n // 2) - 1] badly needs parenthesizing, it misreads very easily. I'd write C[(n//2) : (3*n//2)-1]
  • But I seriously recommend you use numpy vectors, not Python lists. Adding, multiplying etc. are vectorized.
smci
  • 26,085
  • 16
  • 96
  • 138