· 7 years ago · Oct 29, 2018, 02:08 PM
1# File: KnuthMorrisPratt.py
2# Author: Keith Schwarz (htiek@cs.stanford.edu)
3#
4# An implementation of the Knuth-Morris-Pratt (KMP) string-matching algorithm.
5# This algorithm takes as input a pattern string P and target string T, then
6# finds the first occurrence of the string T in the pattern P, doing so in time
7# O(|P| + |T|). The logic behind the algorithm is not particularly complex,
8# though getting it to run in linear time requires a few non-obvious tricks.
9#
10# To motivate KMP, consider the naive algorithm for trying to match a pattern
11# string P against a target T. This would work by considering all possible
12# start positions for the pattern P in the target T, then checking whether a
13# match exists at each of those positions. For example, to match the string
14# ABC against the target string ABABABACCABC, we'd get
15#
16# ABABABACCABC
17# ABX (first two characters match, last does not)
18# X (first character doesn't match)
19# ABX (first two characters match, last does not)
20# X (first character doesn't match)
21# ABX (first two characters match, last does not)
22# X (first character doesn't match)
23# AX (first character matches, second doesn't)
24# X (first character doesn't match)
25# X (first character doesn't match)
26# ABC (match found)
27#
28# This algorithm runs in O(mn) in the worst case, where m = |T| and n = |P|,
29# because it has to do O(n) work to check whether the string matches O(m) times
30# for each spot in the string.
31#
32# However, a lot of this is wasted work. For example, in the above example,
33# consider what happens when we know that the string ABC does not match the
34# first part of the string, ABA. At this point, it would be silly to actually
35# try to match the string at the string starting with the B, since there's no
36# possible way that the string could match there. Instead, it would make more
37# sense to instead start over and try matching ABC at the next A. In fact,
38# more generally, if we can use the information we have about what characters
39# we already matched to determine where we should try to resume the search in
40# the string, we can avoid revisiting characters multiple times when there's no
41# hope that they could ever match.
42#
43# The idea we'll use is to look for "borders" of a string, which are substrings
44# that are both a prefix and suffix of the string. For example, the string
45# "aabcaa" has "aa" as a border, while the string "abc" just has the empty
46# string as a border. Borders are useful in KMP because they encode
47# information about where we might need to pick up the search when a particular
48# match attempt fails. For example, suppose that we want to match ABABC
49# against the string ABABABC. If we start off by trying to match the string,
50# we'll find that they overlap like this:
51#
52# ABABABC
53# ABABx
54#
55# That is, the first four characters match, but the fifth does not. At this
56# point, rather than naively restarting the search at the second character (B),
57# or even restarting it at the third position (A), we can instead note that we
58# can treat the last two characters of the string we matched (AB) as the first
59# two characters of the pattern string ABABC if we just treated it instead as
60# though we had
61#
62# ABABABC
63# ABABx
64# ABABC
65#
66# If we can somehow remember the fact that we already matched the AB at the
67# start of this string, we could just confirm that the three characters after
68# it are ABC and be done. There's no need to confirm that the characters at
69# the front match.
70#
71# In order to make this possible, we'll construct a special data structure
72# called the "fail table." This table stores, for each possible prefix of the
73# string to match, the length of the longest border of that prefix. That way,
74# when we find a mismatch, we know where the next possible start location could
75# be found. In particular, once we have a mismatch, if there's any border of
76# the prefix of the pattern that we matched so far, then we can treat the end
77# of that matching prefix as the start of a prefix of the word that occurs
78# later in the target.
79#
80# The basic idea behind KMP is, given this table, to execute the following:
81#
82# - Guess that the string starts at the beginning of the target.
83# - Match as much of the string as possible.
84# - If the whole string matched, we're done.
85# - Otherwise, a mismatch was found. Look up the largest border of the
86# string that was matched so far in the failure table. Suppose it has
87# length k.
88# - Update our guess of the start position to be where that border occurs
89# in the portion matched so far, then repeat this process.
90#
91# Notice that once we've matched a character against part of the pattern (or
92# found that it can't possibly match), we never visit that character again.
93# This is responsible for the fast runtime of the algorithm (though I'll give a
94# more formal description later on).
95
96# Function: failTable(pattern)
97# Usage: failTable("This is a string!")
98# -----------------------------------------------------------------------------
99# Given a string, constructs the KMP failure table for that string. The values
100# in the table are defined as
101#
102# table[i] = |LongestProperBoundary(pattern[0:i)])|
103#
104# Where the longest proper boundary of a string is the longest proper substring
105# of that string that is both a prefix and a suffix. For example, given the
106# string "abcabc," the longest proper boundary is abc. Similarly, given the
107# string "apple," the longest proper boundary is the empty string.
108#
109# As a sample output of this function, given the string "ababcac", the table
110# would be
111#
112# a b a b c a c
113# * 0 0 1 2 0 1 0
114#
115# This means, for example, that the longest proper boundary of the prefix "aba"
116# has length 1, while the longest proper boundary of the string as a whole is
117# the empty string. Notice that the first entry is *, which we have chosen
118# because there is no mathematically well-defined proper substring of the empty
119# string. We can put anything we want there, and we'll go with None.
120#
121# To compute the values of this table, we use a dynamic programming algorithm
122# to compute a slightly stronger version of the function. We define the
123# function "Extended Longest Proper Boundary" (xLPB) as follows:
124#
125# xLPB(string, n, char) = The longest proper boundary of string[0:n] + char
126#
127# The idea behind this function is that we want to be able to recycle the
128# values of the longest proper boundary function for smaller prefixes of the
129# string in order to compute the longest proper boundary for longer prefixes.
130# To make this easier, the xLPB function allows us to talk about what would
131# happen if we extended the longest proper boundary of some prefix of the
132# string by a single character. Notice that for any nonzero n, we have that
133#
134# LongestProperBoundary(string[0:n]) = xLPB(string, n - 1, string[n])
135#
136# That is, we simply tear off the last character and use it as the final
137# argument to xLPB. Given this xLPB function, we can compute its values
138# recursively using the following logic. As a base case, xLPB(string, 0, char)
139# is the longest proper boundary of string[0:0] + char = char. But this has
140# only one proper boundary, the empty string, and so its value must be zero.
141#
142# Now suppose that for all n' < n we have the value of xLPB(string, n', char)
143# for any character char. Suppose we want to go and compute
144# xLPB(string, n, char). Let's think about what this would mean. Given that
145# n is not zero, we can think of this problem as trying to find the longest
146# proper boundary of this string:
147#
148# +------------+---+------------+------------+---+
149# | LPB | ? | ... | LPB | c |
150# +------------+---+------------+------------+---+
151#
152# ^ ^ ^
153# +----------------------+-------------------+ |
154# | |
155# String of length n New character
156#
157# The idea is that we have the original string of length n, followed by our new
158# character char (which we'll abbreviate c). In this diagram, I've marked the
159# LPB of the string of length n. Notice that right after the LPB at the prefix
160# of the string, we have some character whose value is unknown (since n != 0
161# and the LPB can't be the whole string). If this value is equal to c, then
162# the LPB of the whole string can be formed by simply extending the LPB of the
163# first n characters. There can't be a longer proper boundary, since otherwise
164# we could show that by taking that longer boundary and dropping off the
165# character c, we'd end up with a longer proper boundary for the first n
166# characters of the string, contradicting that we chose the longest proper
167# boundary.
168#
169# By our above argument, remember that the length of the longest proper
170# boundary of the first n characters of the string is given by
171#
172# xLPB(string, n - 1, string[n - 1])
173#
174# Thus we have the first part of our recurrence, which is defined as
175#
176# xLPB(string, n, char) =
177# if n = 0, then 0.
178# let k = xLPB(string, n - 1, string[n - 1])
179# if string[k] == char, return k + 1
180# else, ???
181#
182# Now, suppose that we find that the character after the LPB does not match.
183# If this happens, we can then make the following observation. Below I've
184# reprinted the above diagram:
185#
186# +------------+---+------------+------------+---+
187# | LPB | ? | ... | LPB | c |
188# +------------+---+------------+------------+---+
189#
190# ^ ^ ^
191# +----------------------+-------------------+ |
192# | |
193# String of length n New character
194#
195# Notice that any LPB of this new string must be a prefix of the LPB of the
196# first n characters and a suffix of the LPB followed by the character c.
197# Since by definition the LPB of the first n characters must be a prefix of
198# those n characters, we have the following elegant conclusion to our
199# recurrence:
200#
201# xLPB(string, n, char) =
202# if n = 0, then 0.
203# let k = xLPB(string, n - 1, string[n - 1])
204# if string[k] == char, return k + 1
205# else, xLPB(string, k, char)
206#
207# The reason for this is that xLPB(string, k, char) asks for the longest
208# proper boundary of the LPB of the string formed from the first n characters
209# of the string followed by the character c, which is exactly what we described
210# above.
211#
212# As written, filling in the table of LPB values would take O(n^2) time, where
213# n is the length of the string. However, using dynamic programming and an
214# amortized analysis, we can show that this function can be made to run in
215# O(n) time. In particular, suppose that for all n' < n, we know the value of
216# LPB(string[0:n]). Then in the above formulation of xLPB, the first
217# recursive call is known, and the only recursive call we may actually need to
218# make is the second.
219#
220# However, this doesn't seem to say anything about the runtime of the second
221# recursive call, which seems as though it might cause the evaluation of this
222# function to run in time O(n). This is correct, but in an *amortized* sense
223# the whole table can still be computed in O(n) time overall. To see this,
224# let's define a potential function Phi(k) that associates a potential at each
225# point of the computation of the table. In particular, define Phi(k) as
226#
227# Phi(0) = 0
228# Phi(k + 1) = result[k - 1]
229#
230# Here, result is the resulting table of LPB values. Because of this, we can
231# remark that result[k] < k, since the longest proper border of a string can't
232# be any longer than that string.
233#
234# Let's now show that this potential function gives an amortized O(1) cost for
235# each table entry computation, and thus an O(n) overall runtime for the table-
236# building algorithm. To see this, consider what happens when the logic to
237# compute the next value runs. The runtime for this step is bounded by the
238# number of recursive calls made to a subproblem. However, each subproblem is
239# then of size given by the LPB of a slightly smaller problem. This subproblem
240# must then have size at most the size of that smaller subproblem. In other
241# words, we can say that each recursive call drops the maximum possible value
242# of the LPB for the current prefix by at least one. Consequently, if k
243# recursive calls are made, the LPB of the current prefix is at least k smaller
244# than the LPB of the previous prefix, and so
245#
246# D Phi = -k
247#
248# And so the amortized cost of computing the next term is 1 + k - k = O(1).
249
250def failTable(pattern):
251 # Create the resulting table, which for length zero is None.
252 result = [None]
253
254 # Iterate across the rest of the characters, filling in the values for the
255 # rest of the table.
256 for i in range(0, len(pattern)):
257 # Keep track of the size of the subproblem we're dealing with, which
258 # starts off using the first i characters of the string.
259 j = i
260
261 while True:
262 # If j hits zero, the recursion says that the resulting value is
263 # zero since we're looking for the LPB of a single-character
264 # string.
265 if j == 0:
266 result.append(0)
267 break
268
269 # Otherwise, if the character one step after the LPB matches the
270 # next character in the sequence, then we can extend the LPB by one
271 # character to get an LPB for the whole sequence.
272 if pattern[result[j]] == pattern[i]:
273 result.append(result[j] + 1)
274 break
275
276 # Finally, if neither of these hold, then we need to reduce the
277 # subproblem to the LPB of the LPB.
278 j = result[j]
279
280 return result
281
282# Function: kmpMatch(needle, haystack)
283# Usage: print kmpMatch("0101", "0011001011") # Prints 5
284# -----------------------------------------------------------------------------
285# Uses the KMP algorithm to find an occurrence of the specified needle string
286# in the haystack string. To do this, we compute the failure table, which
287# is done above. Next, we iterate across the string, keeping track of a
288# candidate start point and length matched so far. Whenever a match occurs, we
289# update the length of the match we've made. On a failure, we update these
290# values by trying to preserve the maximum proper border of the string we were
291# able to manage by that point.
292def kmpMatch(needle, haystack):
293 # Compute the failure table for the needle we're looking up.
294 fail = failTable(needle)
295
296 # Keep track of the start index and next match position, both of which
297 # start at zero since our candidate match is at the beginning and is trying
298 # to match the first character.
299 index = 0
300 match = 0
301
302 # Loop until we fall off the string or match.
303 while index + match < len(haystack):
304 print index, match
305
306 # If the current character matches the expected character, then bump up
307 # the match index.
308 if haystack[index + match] == needle[match]:
309 match = match + 1
310
311 # If we completely matched everything, we're done.
312 if match == len(needle):
313 return index
314
315 # Otherwise, we need to look at the fail table to determine what to do
316 # next.
317 else:
318 # If we couldn't match the first character, then just advance the
319 # start index. We need to try again.
320 if match == 0:
321 index = index + 1
322
323 # Otherwise, see how much we need to skip forward before we have
324 # another feasible match.
325 else:
326 index = index + match - fail[match]
327 match = fail[match]
328
329 # If we made it here, then no match was found.
330 return None