Introduction to Programming, Aug-Nov 2008 Lecture 22, Monday 19 Nov 2008 Memoization using arrays ------------------------ Recall that we proposed to modularize the memo table by defining a datatype Table a b to store the values of functions of type a->b supporting the following functions. emptytable :: returns an empty memo table memofind :: tells whether a value exists in the memo table memolookup :: looks up a value in the memo table memoupdate :: adds a new entry to the memo table Here is an implementation of a memo using Arrays. The function "emptytable" initializes all entries in the array to "initval". To check whether a value exists in the table, "memofind" compares the existing value to "initval". For this, we have to store "initval" permanently in the table, so the definition of Table a b includes both "Array a b" for the memo and a separate value of type b, the initial value. Futher, since "memofind" checks equality for values im the array, there is a dependence Eq b. "memolookup" just returns the value at position i in the array using the operator "!" while "memoupdate" updates the value at position i using the "//" operator. module Arraymemo(Table,emptytable,memofind,memolookup,memoupdate) where import Array data (Ix a,Eq b) => Table a b = T (Array a b) b emptytable :: (Ix a,Eq b) => (a,a) -> b -> Table a b emptytable bounds initval = T (array bounds [ (i,initval) | i <- range bounds]) initval memofind :: (Ix a,Eq b) => (Table a b) -> a-> Bool memofind (T table initval) index | table!index == initval = False | otherwise = True memolookup :: (Ix a,Eq b) => (Table a b) -> a -> b memolookup (T table initval) index = table!index memoupdate :: (Ix a,Eq b) => (Table a b) -> (a,b) -> (Table a b) memoupdate (T table initval) (index,value) = T (table//[(index,value)]) initval Now, let us see how to memoize the recursive program to compute the number of paths in a rectangular grid. Recall that we have a rectangular array of intersections with m rows and n columns, with each intersection connected by a one way road going right and another one way road going down, so from intersection (i,j) one can go right to intersection (i,j+1) or down to intersection (i+1,j). Some intersections are blocked and we want to compute the number of paths from (1,1) to (m,n). We assume the input is 0-1 array grid with indices (1,1) to (m,n) such that grid(i,j) = 1 if the intersection is open 0 if the intersection is blocked Let gridpaths(i,j) denote the number of paths from (i,j) to (m,n). In general, we have the recursive definition gridpaths(i,j) = gridpaths(i,j+1) + gridpaths(i+1,j) The base case is gridpaths(m,n) = 1 If grid(i,j) = 0, we have gridpaths(i,j) = 0 since no paths can pass through intersection (i,j) For the bottom row, we have gridpaths(m,j) = gridpaths(m,j+1) and for the rightmost column we have gridpaths(i,n) = gridpaths(i+1,n) Here is a recursive definition of gridpaths in Haskell. The first argument is the input grid and the second argument is the location of an intersection, so hgridpaths grid (i,j) computes the function gridpaths(i,j) described above. hgridpaths :: (Array (Int,Int) Int) -> (Int,Int) -> Int hgridpaths grid (i,j) = | grid!(i,j) == 0 = 0 | i == m && j == n = 1 | i == m = hgridpaths grid (i,j+1) | j == n = hgridpaths grid (i+1,j) | otherwise = (hgridpaths grid (i,j+1)) + (hgridpaths grid (i+1,j)) where ((lr,lc),(ur,uc)) = bounds grid m = ur - lr + 1 n = uc - lc + 1 To memoize this, we write the following. Recall that the memoized function takes the current memo as an extra argument and returns a pair of (value, newmemo). module Gridpaths where import Arraymemo hgridpaths :: (Array (Int,Int) Int) -> (Int,Int) -> Int hgridpaths grid (i,j) = fst (memogridpaths grid (i,j) emptymemo) where ((lr,lc),(ur,uc)) = bounds grid m = ur - lr + 1 n = uc - lc + 1 emptymemo :: Table (Int,Int) Int emptymemo = emptytable ((1,m),(1,n)) (-1) memogridpaths :: (Array (Int,Int) Int) -> (Int,Int) -> (Table (Int,Int) Int) -> (Int,Table (Int,Int) Int) memogridpaths grid (i,j) memo | memofind memo (i,j) = (memolookup memo (i,j),memo) | grid!(i,j) == 0 = (0, memoupdate memo ((i,j),0)) | i == m && j == n = (1, memoupdate memo ((i,j),1)) | i == m = (lastrowvalue, lastrowmemo) | j == n = (lastcolvalue, lastcolmemo) | otherwise = (valc,memoc) where (lastrowvalue,memo1) = memogridpaths grid (i,j+1) memo lastrowmemo = memoupdate memo1 ((i,j),lastrowvalue) (lastcolvalue,memo2) = memogridpaths grid (i+1,j) memo lastcolmemo = memoupdate memo1 ((i,j),lastcolvalue) (vala,memoa) = memogridpaths grid (i+1,j) memo (valb,memob) = memogridpaths grid (i,j+1) memoa valc = vala + valb memoc = memoupdate memob ((i,j),valc) Longest common subsequence -------------------------- The unix utility "diff" compares two text files and reports the difference between them in an optimal way. To do this, it tries to optimally match the two files line by line. For instance, if file1 has 4 lines x1,x2,x3,x4 and file2 has 5 lines y1,y2,y3,y4,y5 such that x1 = y2, x3 = y4 and x4 = y5, diff will match up x1,x3,x4 with y2,y4,y5 and report that file2 has an extra line y1 before this sequence, and between x1=y2 and x3=y4, file1 has an extra line x2 while file2 has an extra line y3. In this scenario, the matching sequences x1,x3,x4 and y2,y4,y5 are subsequences of the two original sequences of lines. What diff does is to compute the longest common subsequence of lines between two files. Abstractly, we can think of the two sequences as consisting of a set of letters, such as V = T G C A T G G A T C W = G T T T G C A C Now we can drop some positions and obtain a subsequence. For instance "T G A C" is a subsequence of both words. "G T G A C" is a longer common subsequence as is "T G C A C". These are (probably!) the longest common subsequences in this example, which shows that the longest common subsequence is not unique. Our aim is to find the longest common subsequence (lcs), but for the moment we concentrate on finding the *length* of the longest common subsequence (llcs) instead. Consider the example above. When we match the first G in W to a G in V, we can match it to either the G at position 2 in V or at position 6 or 7 in V. Suppose we have a solution in which we match the first G in to the G at position 7 in V. Then, all further matches from W will be to the right of position 7 in V. If we take this matching and shift the match for the first G back to the G at position 2, we still have a matching. And, by moving this match back, perhaps we might enable longer matches between W and V. Thus, in general it is best to match a letter as far to the left as possible. Let LLCS(i,j) denote the length of the longest common subsequence starting at position i in V and at position j in W. We have two cases to consider. 1. Suppose V[i] == W[j]. If there is any matches that uses V[i] or W[j], we can shift it so that V[i] matches W[j] and get a new matching that is at least as good. So, we may as well match V[i] to W[j] and continue from (i+1,j+1). This suggests that If V[i] == W[j], LLCS(i,j) = 1 + LLCS(i+1,j+1) 2. On the other hand, suppose V[i] /= W[j] We clearly cannot use both V[i] and W[j] in the final matching --- if we use V[i], it matches with W[k] for some k > j and nothing beyond V[i] can then match W[j]. A symmetric argument shows that if we use W[j] in the final matching, V[i] can play no role. However, a priori, we have no information about which of V[i] or W[j] to drop from the solution, so we try both and take the maximum. If V[i] /= W[j], LLCS(i,j) = max(LLCS(i+1,j), <-- drop V[i] LLCS(i,j+1)) <-- drop W[j] 3. Base cases If V is of length m and W is of length n, we have LLCS(m,n) = 1, if V[m] == W[n] 0, otherwise If we are at the last letter of V, we have LLCS(m,j) = 1, if V[m] == W[j] LLCS(m,j+1), otherwise Symmetrically, at the last letter of W, we have LLCS(i,n) = 1, if V[i] == W[n] LLCS(i+1,n), otherwise We can derive a dynamic programming solution to this problem by filling up the values LCSS[i,j] in an M x N array, starting with the last row and column as base cases.