title: How to speed up Levenshtein distance calculation
body: <p>I am trying to run a simulation to test the average <a href="http://en.wikipedia.org/wiki/Levenshtein_distance" rel="noreferrer">Levenshtein distance</a> between random 
  binary strings.</p>
  
  <p>My program is in python but I am using this <a href="https://github.com/miohtama/python-Levenshtein" rel="noreferrer">C extension</a>. The function that is relevant and takes most of the time computes the Levenshtein distance between two strings and is this.</p>
  
  <pre><code>lev_edit_distance(size_t len1, const lev_byte *string1,
                    size_t len2, const lev_byte *string2,
                    int xcost)
  {
    size_t i;
    size_t *row;  /* we only need to keep one row of costs */
    size_t *end;
    size_t half;
  
    /* strip common prefix */
    while (len1 &gt; 0 &amp;&amp; len2 &gt; 0 &amp;&amp; *string1 == *string2) {
      len1--;
      len2--;
      string1++;
      string2++;
    }
  
    /* strip common suffix */
    while (len1 &gt; 0 &amp;&amp; len2 &gt; 0 &amp;&amp; string1[len1-1] == string2[len2-1]) {
      len1--;
      len2--;
    }
  
    /* catch trivial cases */
    if (len1 == 0)
      return len2;
    if (len2 == 0)
      return len1;
  
    /* make the inner cycle (i.e. string2) the longer one */
    if (len1 &gt; len2) {
      size_t nx = len1;
      const lev_byte *sx = string1;
      len1 = len2;
      len2 = nx;
      string1 = string2;
      string2 = sx;
    }
    /* check len1 == 1 separately */
    if (len1 == 1) {
      if (xcost)
        return len2 + 1 - 2*(memchr(string2, *string1, len2) != NULL);
      else
        return len2 - (memchr(string2, *string1, len2) != NULL);
    }
    len1++;
    len2++;
    half = len1 &gt;&gt; 1;
    /* initalize first row */
    row = (size_t*)malloc(len2*sizeof(size_t));
    if (!row)
      return (size_t)(-1);
    end = row + len2 - 1;
    for (i = 0; i &lt; len2 - (xcost ? 0 : half); i++)
      row[i] = i;
  
    /* go through the matrix and compute the costs.  yes, this is an extremely
     * obfuscated version, but also extremely memory-conservative and relatively
     * fast.  */
    if (xcost) {
      for (i = 1; i &lt; len1; i++) {
        size_t *p = row + 1;
        const lev_byte char1 = string1[i - 1];
        const lev_byte *char2p = string2;
        size_t D = i;
        size_t x = i;
        while (p &lt;= end) {
          if (char1 == *(char2p++))
            x = --D;
          else
            x++;
          D = *p;
          D++;
          if (x &gt; D)
            x = D;
          *(p++) = x;
        }
      }
    }
    else {
      /* in this case we don't have to scan two corner triangles (of size len1/2)
       * in the matrix because no best path can go throught them. note this
       * breaks when len1 == len2 == 2 so the memchr() special case above is
       * necessary */
      row[0] = len1 - half - 1;
      for (i = 1; i &lt; len1; i++) {
        size_t *p;
        const lev_byte char1 = string1[i - 1];
        const lev_byte *char2p;
        size_t D, x;
        /* skip the upper triangle */
        if (i &gt;= len1 - half) {
          size_t offset = i - (len1 - half);
          size_t c3;
  
          char2p = string2 + offset;
          p = row + offset;
          c3 = *(p++) + (char1 != *(char2p++));
          x = *p;
          x++;
          D = x;
          if (x &gt; c3)
            x = c3;
          *(p++) = x;
        }
        else {
          p = row + 1;
          char2p = string2;
          D = x = i;
        }
        /* skip the lower triangle */
        if (i &lt;= half + 1)
          end = row + len2 + i - half - 2;
        /* main */
        while (p &lt;= end) {
          size_t c3 = --D + (char1 != *(char2p++));
          x++;
          if (x &gt; c3)
            x = c3;
          D = *p;
          D++;
          if (x &gt; D)
            x = D;
          *(p++) = x;
        }
        /* lower triangle sentinel */
        if (i &lt;= half) {
          size_t c3 = --D + (char1 != *char2p);
          x++;
          if (x &gt; c3)
            x = c3;
          *p = x;
        }
      }
    }
  
    i = *end;
    free(row);
    return i;
  }
  </code></pre>
  
  <p>Can this be sped up?</p>
  
  <p>I will be running the code in 32 bit ubuntu on an AMD FX(tm)-8350 Eight-Core Processor. </p>
  
  <p>Here is the python code that calls it.</p>
  
  <pre><code>from Levenshtein import distance
  import random
  for i in xrange(16):
      sum = 0
      for j in xrange(1000):
          str1 = bin(random.getrandbits(2**i))[2:].zfill(2**i)
          str2 = bin(random.getrandbits(2**i))[2:].zfill(2**i)
          sum += distance(str1,str2)
      print i,sum/(1000*2**i)
  </code></pre>