divide and conquer

(a  bi)(c  di) =
  (ac-bd)  (bc + ad)i   4 mult, 2 add/sub

x: ac   y: bd   z: (a+b)(c+d)

(a  bi)(c  di) =
  (x-y)  (z-x-y)i       3 mult, 5 add/sub

def bitlength(n): return len(bin(n))-2

def leftright(x,n): # x > 1
  left  = x >> n/2
  right = x - (left << n/2)
  return left, right

def fastmul(x,y): # x,y > 0
  n = bitlength(x)
  if n <= 4:
    return x*y
  xl, xr = leftright(x,n)
  yl, yr = leftright(y,n)
  pl = fastmul(xl,yl)
  pr = fastmul(xr,yr)
  pm = fastmul(xl+xr,yl+yr)
  return (pl << 2*(n/2)) + ((pm - pl - pr) << n/2) + pr

for j in range(400):
  for k in range(400):
    assert( fastmul(j,k) == j*k)

- assume 7 digit inputs
- shift <-  7/2 = 3
- split inputs:  aaaa bbb    cccc ddd
-                 a'   b'     c'   d'
- L <- a' * c'
- R <- b' * d'
- S <- (a' + b') * (c' + d')

- return     L * base^(2*shift)  +
          (S-L-R) * base^shift   +
             R
- correct ?
                        aaaa bbb
                     *  cccc ddd
                        --------
       llllllll  mmmmmmm  rrrrrr
         a'c'   a'd'+b'c'   b'd'

- want L = a'c'                check :)
- want R = b'd'                check :)
- want M = S-L-R = b'c' + a'd'       ???
S-L-R = a'c' + b'c' + a'd' + b'd - a'c' - b'd'
      =        b'c' + a'd'     check :)

return x*y if x <= 10 or y <= 10
(1234,5678)
   (12,56)
      (1,5) 5
      (2,6) 12
      1+2,5+6
      (3,11) 33
      5  33-5-12  12
      500 + 160 + 12
      672
  (34,78)
      (3,7) 21
      (4,8) 32
      3+4,7+8
      (7,15) 105
      21 105-21-32 32
      2100 + 520 + 32
      2652
   12+34,56+78
   (46,134)
       (4,13) 52
       (6,4)  24
       4+6,13+4
       (10,17) 170
       52 170-52-24 24
       5200 + 940 + 24
       6164
    672  6164-672-2652  2652
    6720000 + 284000 + 2652
    7006652

def ms(x): # based on:
#http://stackoverflow.com/questions/18761766/mergesort-python
  result = []
  if len(x) < 2:
    return x
  mid = len(x)/2 ; y = ms(x[:mid]) ; z = ms(x[mid:])
  i,j = 0,0
  while i < len(y) and j < len(z):
    if y[i] > z[j]:
      result.append(z[j]) ; j += 1
    else:
      result.append(y[i]) ; i += 1
  result += y[i:] + z[j:]
  return result

def select(L,k):
  A = []; B = [];  C = []
  v = L[random.randint(0,len(L)-1)]
  for j in L:
    if j<v:
      A.append(j)
    elif j>v:
      C.append(j)
    else:
      B.append(j)
  if k <= len(A):
    return select(A,k)
  if k <= len(A) + len(B):
    return v
  return select(C,k-len(A)-len(B))

# partition, qs:http://hetland.org/coding/python/quicksort.html
def partition(list, start, end):
  pivot = list[end]    # partition around the last value
  bottom = start-1     # initially outside sublist to be partitioned
  top = end            # "         "
  done = 0
  while not done:      # until all elements partitioned...
    while not done:  # until we find an out of place element...
      bottom = bottom+1  # move the bottom up.
      if bottom == top:  # if hit top...
        done = 1       # ... we are done.
        break
      if list[bottom] > pivot:     # bottom out of place?
        list[top] = list[bottom] # put it at top...
        break                    # ... resume search from top
    while not done:        # until we find an out of place element...
      top = top-1        # move top down.
      if top == bottom:  # if hit bottom...
        done = 1       # ... we are done.
        break
      if list[top] < pivot:        # top out of place?
        list[bottom] = list[top] # put it at bottom...
        break                    # ... resume search from bottom
  list[top] = pivot                    # replace pivot
  return top                           # return split index

def quicksort(list, start, end):
  if start < end:        # two or more elements?
    split = partition(list, start, end)
    quicksort(list, start, split-1)
    quicksort(list, split+1, end)
  return