Go deh!

A better candidate selection process

2011-11-22T06:52:00.001Z

I just read an article that had a good go at explaining why Silicon Valley firms are filled with white middle class males that managed to actually do what it set out to do when it explained why it wanted to take a more reasoned approach. They boil it down to everyone having an in-built bias - not necessarily a good or bad thing, it is just how humans operate. The article then comes up with a very useful experiment:

When selecting from resumes for interview, have the names, age, sex, and origin of the candidates blanked out.

The author found that his resulting selections for interview were a lot different to what he had before.

He goes on to discuss how orchestras changed from being all male affairs to the current mix of genders by them adopting a selection process where candidates played unseen, behind a screen, and where selected based on how well they played.

This got me wondering about the current controversy about Oxbridge admissions. Maybe the universities should adopt similar schemes, candidates should be lead to a separate room, unseen, where voice changers would mask ethnicity and accent, The interview would be audio only and taped, with the selection panel in a separate room. After the interview, the candidate should be given a copy of the tape, and another copy kept for independent review. The idea would be to make the process a better meritocracy rather than the obvious selection of "people like me" that goes on at the moment.

Oh, here's the link to the original article.

Should you worry about a 2x speedup?

2011-11-08T05:04:00.000Z

Let's take as context an implementation of a task for Rosetta Code - a site set up to compare how different programming languages are used to implement the same task for over five hundred tasks and over four hundred languages.

My short answer would be it depends! You need to:

Read all of the task description,
Read some of the solutions in other languages,
And maybe skim the tasks talk page

I.e. ensure you know what the task is asking for and then verify that both solutions solve the task as stated. It can be very easy to misinterpret what the task is asking for, for example, if a task asks for a particular algorithm, do both the examples you are comparing use that algorithm?

As well as comparing two implementations for speed, you should also compare for readability. How well the code reads can have a large impact on how easy the code is to maintain. It has been known for task descriptions to be modified; someone tracking that modification may need to work out if and how code needs to be updated. If an example is overly complex and/or unidiomatic then it could cause problems.

Time complexity. If one version of code works better when given 'bigger' data then you need to know more about when that happens - it could be that the cross-over point in terms of speed of execution is never likely to be met. Maybe the size of data needed to reach cross over is unreasonable to expect, or that other mechanisms come into play that mask predicted gains (in other words you might need to verify using that actual bigger data set to account for things like swapping or caching at the OS and hardware level.

How fast does it need to be? Rosetta code doesn't usually mention absolute speed of execution, but if one example takes ten hours and the other takes five then you might want to take that into account. If one example took 0.2 seconds and the other only 0.1 seconds then I guess there is an unwritten expectation that examples "don't take long to run" where long is related to the expectation and patience of the user.

You need to look at the context. In the case of Rosetta code, it may be best to give a solution using a similar algorithm to other examples, or a solution that shows accepted use of the language.

When you make your considered choice, you might want to squirrel away the losing code with notes on why it wasn't used, - On Rosetta Code we sometimes add more than one solution to a task with comments contrasting the two if they both have merit.

It seems to me that talk about optimising for speed, and speed comparisons tends to dominate on the web over other optimisations, (usually with no extra info on the accuracy of the result. Actually there might be more cases of a revised result that showed not even the first digit of the original answer was right, but more than two digits of precision were shown in the answers)!

Woa - a glitch on Wolfram Mathworld?

2011-06-08T01:12:00.000+01:00

I had come across Kaprekar numbers and decided to write a Rosetta Code task around the series. There were issues about the need for extra explanation of why 1 is a member of the series.
I had originally read the wp description and noted that there is a main part to testing if a number is in the series:

In mathematics, a Kaprekar number for a given base is a non-negative integer, the representation of whose square in that base can be split into two parts that add up to the original number again. For instance, 45 is a Kaprekar number, because 45² = 2025 and 20+25 = 45.

That is, splitting the ordered digits of the square of the number into two parts that sum to the original number.

Wikipedia then goes on to state that runs of all zeros are not considered positive integers and so 100*100 = 10,000 and 100 + 000 = 100 is disallowed.

What got me was that 1 is considered a member of the series. The only way I could shoe-horn it into the general rule is to allow "splitting" the digits of the square before the first or after the last digit and so producing a number with no digits in it as a partial sum that is allowed and has value zero?!

I widened my references and took a look at the Wolfram Mathworld entry. Its explanation for the series does not allow for 1 being a member, but shows it as such anyway. It is confused as it describes a series that would exclude numbers 4879 and 5292, i.e. Sloanes A053816; but points to Sloanes A006886!

Sloanes A006886 gives another definition for the series:

Kaprekar numbers: n such that n=q+r and n^2=q*10^m+r, for some m >= 1, q>=0 and 0<=r<10^m, with n != 10^a, a>=1.

If we try n = 1 then n*n = 1 and there is no m, q, and r that satisfies all the conditions.

I will submit a comment to Mathworld and see what they have to say.

P.S. I could have started this blog entry: "I'm not a mathematician but ..."

What http://perl6.org/ gets right

2011-06-01T09:04:00.000+01:00

One thing I really like about the perl6 homepage is that they politely ask for people to be nice to each other up-front. The psychology goes a bit further: they use a female "voice", and associate the niceness with a butterfly image that is repeated on other pages.

Now that is brilliant! Finding a problem and directly addressing it in such an appropriate manner is something I really like.

Good enough solutions can be great!

2011-05-29T21:12:00.001+01:00

I've had to imlement the solutions to two problems on Rosetta code lately and they both involve heuristics that get around the prolems of exhaustive search.

The first problem was in finding a way to shuffle the letters of a string so that a minimum of letters end up in the same place. My first solution was to go through all letter permutations and return the best. That takes for ever on longer strings such as the word 'antidisestablishmentarianism'.
The algorithm I then turned to just swaps letters if they swap into positions that put different characters into corresponding locations. That algorithm is very fast and seems to give good answers.

The second problem was in creating knights tours of chess boards of varying sizes. The given Wansdorff's algorithm is much quicker than a full depth first search, and was straight-forward to code . I later found evidence that for larger board sizes the algorithm fails . (Although the paper does suggest simple improvements that work for boards sizes with thousands of squares).

Coding those two, makes me remember that results are what count normally. You won't achieve much by searching for the theory of relativity if the queen wants Newtons laws of motion so that sailing ships can trade and plunder!

"That Goes Without Saying (or Does It)". It didn't!

2011-03-12T07:29:00.000Z

I spent the time to watch a talk from Larry Wall the Perl, and now Perl 6 author; where he is using the Rosetta Code site for what it was built for - comparing programming language design, by examples.

Unfortunately he was preaching to his acolytes so there is no meat to the talk. Larry says 'X' about a language compared to Perl or Perl 6 and the audience accepts without even asking for clarifications. Argument never enters their minds.

What he did do well is to show examples first from the Root mean square page and indicate how anyone can compare and contrast languages using it. Larry then categorized Spiral matrix and Zigzag matrix as tasks whose examples in Perl where "fun to write" which I liked as I remember getting a lot out of Zigzag when I first found RC, so much so that I think starting the new task "Spiral matrix" was probably the first task I wrote for RC. (I thought it would be good to learn by doing a variant of an existing task).

Larry showed extended versions of the zigzag and spiral matrix solutions that used text based turtle graphics to animate the solutions.

In summary, his language comparisons were droopy; his use of RC - spot on!

New laptop with an i7-640M processor; and a 32 inch "monitor".

2011-01-11T19:00:00.000Z

After a lot of hard work, and the car scrappage schemes of last year, I got rewarded with bonus!?
Nothing like a bwankers bonus, but enough to buy some fripperies after the essentials, and I decided to replace my aging 17" desktop-replacement laptop and get me a large 28" monitor.

After much searching in the post-xmas sales, I came across a curious machine at a very good price. The Advent Sienna 710 is a PCWorld/Dixons/Currys own-branded 15" laptop that has what, at the time, is the highest speed core for a laptop processor - an Intel 17-640M for under £500.

I guess the laptop has such a low price as it will not appeal to a gamer as it relies on Intels graphics accelerator and is so basic it doesn't have the usual scroll area at the side of the pad., which isn't too hot but that makes it ideal for me as I wanted something fast and multi-core, but don't indulge in much 3D work and can live with the pad deficiencies.

Our second TV had packed up mid-2010 and we were making-do with a very old 19" LCD. Instead of getting the 28" monitor I decided to buy a 32" Sony 1080P TV and use that as my monitor.

Well, the Advent laptop and Sony TV work beautifully together. I have a huge 32" screen at a decent resolution and with the sound coming through the HDMI connection, I don't need extra laptop speakers as the Sony sound is very good (better than on our more expensive main TV).

I ran my improved Python multi-processing example on the laptop and was very impressed with the speed, and with all the python instances appearing in the task manager.

from concurrent import futures
from math import floor, sqrt

NUMBERS = [2**n - 1 for n in range(55, 65)]

def lowest_factor(n, _start=3):
    if n % 2 == 0:
        return 2
    search_max = int(floor(sqrt(n))) + 1
    for i in range(_start, search_max, 2):
        if n % i == 0:
            return i
    return n

def prime_factors(n, lowest):
    pf = []
    while n > 1:
        pf.append(lowest)
        n //= lowest
        lowest = lowest_factor(n, max(lowest, 3))
    return pf

def prime_factors_of_number_with_lowest_prime_factor(NUMBERS):
    with futures.ProcessPoolExecutor() as executor:
        low_factor, number = max( (l, f) for l, f in zip(executor.map(lowest_factor, NUMBERS), NUMBERS) )
        all_factors = prime_factors(number, low_factor)
        return number, all_factors


def main():
    print( 'For these numbers:\n  ' + '\n  '.join(str(p) for p in NUMBERS) )
    number, all_factors = prime_factors_of_number_with_lowest_prime_factor(NUMBERS)
    print('    The one with the largest minimum prime factor is %i:' % number)
    print('      All its prime factors in order are: %s' % all_factors)


if __name__ == '__main__':
    main()

END.

Stupid macro languages

2010-12-07T05:29:00.000Z

In the Python web world, the debate on whether templating languages should offer programming language type abilities has reared its head again. from "Stupid templating languages", to "Not So Stupid Template Languages", "In Response to "Stupid Template Languages"", and discussed more on reddit here.

In the world of Electronic Design Automation and integrating complex EDA software tools, the EDA companies,, on the whole, have embraced TCL as their scripting tool, but sometimes I have inherited a complex design flow, where, to add control, and later flexibility within that control, proprietary initialization file formats and their parsers have grown by adding their own macro language type features and grown to need the power of something like the bash shell. Rather than ditching their proprietary lashed-together, one-off languages for something that is tried and tested and used/debugged by many - such as bash/python/Tcl/Lua they persevere by adding more capability to an in-house, one project language.

In the case of web templating the main argument is that given an intelligent web templating language, people are apt to migrate what should be business logic into the templates thereby muddying the separation of the business logic from the presentation layer. In the EDA field the problem is often the maintenance of the proprietary macro language and the proliferation of languages with which one needs to be familiar with in order to master a design environment.

I won't comment more on the web templating debate, but when it comes to initialization file formats and or initialization scripts then you should stick to pre-existing standards such as .ini files or .csv files in the simpler cases, or use pre-existing scripting languages when you need more power such as bash scripts, Python, or Tcl.

Ipod Tablet Opportunity?

2010-10-04T16:06:00.000+01:00

The ipad is close to £600 in the UK. The ipod touch is less than £200. Given that there are tablet computers on sale for less than £200, and inspired by news of a future ipod touch add-on that gives it phone capabilities; how about designing a ten-inch screen, extra battery, and phone dongle that would take a (partially disassembled), ipod touch inside it and turn it into a ten inch ipad-alike?

One problem I see is that you would probably have to provide a service to take the users ipod touch and disassemble/assemble it inside the ten-inch "dock"; but Apple seems to provide a lot of monetary difference to try and make a profit.

Just a thought!

Regression runners: Whot! No Job Arrays?

2010-09-18T08:34:00.004+01:00

I went from running regressions consisting of greater than 50K simulations nightly using Mentors QuestaSim, where I had written the regression framework and results preparation myself - to my first real Cadece vManager controlled regression project.

I should explain to my readers who usually see Python focused blog entries, that this is about the running of large regressions on compute farms for verifying ASICs (although my regression runner used Python extensively).

I did like the more minimal, HTML-based GUI on the tool - I've spent too much time creating work-arounds for flashy GUI's that don't work outside a narrow area, and so know they can be a double-edged sword. The tool needs a lot of configuration behind the scenes, and they have chosen to do this using their own vsif file format which seems to me to be a design fault. If they had gone for XML (gosh is that me advocating the use of XML), or YAML, then the format would be as readable and more easily used in other tools - no tool is an island in a design flow.

The main problem with vManager, is with its use of LSF. LSF is a separate tool used to harness a network of computers into a compute cluster. LSF runs jobs on the cluster, managing the cluster to give maximum job throughput. vManager can turn each test that needs to run on the design into a possible series of jobs for LSF, but it fails to use job arrays! Job arrays allow for many similar jobs to be submitted and controlled as one array of jobs where the indexing of the individual job in the array is used to introduce a variation in the test.

In my regression framework I created arrays of over 100K jobs, that took several weeks to run, and all the time that the job array at the heart of my regression was running, IT could change the parameters of my job array to control its use of resources according to changing priorities; scale back the maximum number of concurrent simulations during the 9-till-5; ramp up during weekends, (or stop then resume for a 'planned' license outage). IT had all the information they need to show the full extent of the regression and LSF will give the throughput of the array from which I can estimate completion times. Another engineer wishing to start a regression has immediate feedback of what else is running rather than a huge list of hundreds of jobs that it is difficult to make sense of. In short, ASIC verification regressions map naturally to Job Arrays.

The other problem I have with vManager is a missing feature: Randomization of the order that tests are run. If you have 100 tests to run on five areas of the design then it is natural to compose the vsif file that describes each test in some sort of orderly progression through the tests to run. When the regression is running and you are monitoring partial results then it is easier to get a sense of the overall result if the individual tests are run in a random order. In my regression framework I generated tests in an order, then randomized the set of tests based on Pythons excellent Mersenne Twister implementation. By monitoring results during previous runs I could judge how much of a regression had to run if we needed results in a specified time.

A bit of background (re?)search showed that Mentor have their own regression running framework, that also doesn't use job arrays; and that the Grid Engine cluster management tool supports job arrays.

McCarthy's amb operator in Python

2010-08-23T18:00:00.007+01:00

From what I can gather after some reading [1], [2], [3]; McCarthy's amb operator, would allow you to give the possible values for variables; then give functions that constrain the values; and finally extract values for those variables that satisfy the constraints just by calling amb without any arguments.

I decided to work backwards from how I thought it should be used to an implementation.

Given the problem of finding Pythagorean triples for integers <=10, i.e define three variables x, y, and z that can have values 1..10 subject to the constraint that x*x +y*y == z*z. What are the possible values of x, y and z that satisfy the constraint? For an amb-type solution adapted a little to python, being amb-like is all in the way that the problem is expressed. I would like to write something like:


print("\nSmall Pythagorean triples problem")
# Set rage of global variables
x = amb(range(1,11))
y = amb(range(1,11))
z = amb(range(1,11))

# Extract solutions to global variables and print
for _dummy in amb( lambda x, y, z: x*x + y*y == z*z ):
print x, y, z

I took a little liberty in the way I had seen other implementations work, as I had seen examples where the call of amb with an argument of a function just set things up and 'registered' the function, then subsequent calls of amb() without any arguments, would set the global variables used in the lambda function to successive values that satisfy the constraint. I reasoned that calls to amb after registering the constraint function should follow the Python iterator protocol ( with an __iter__ and a next method), so all the values, if any, could be accessed in a loop. I liked the idea of automatically setting global variables so decided to give that a go (but I was prepared to jettison this setting global variables if it proved long-winded).

The solution follows. To modify global names, I noted that functions have f.func_globals which is the globals environment in which the function was created. I get this from the lambda expression. I also cheat a little in requiring the parameter names of the constraint function to be the global names that were previously assigned value ranges.


import itertools as _itertools

class Amb(object):
 def __init__(self):
     self._names2values   = {}       # set of values for each global name
     self._func           = None     # Boolean constraint function
     self._valueiterator  = None     # itertools.product of names values
     self._funcargnames   = None     # Constraint parameter names

 def __call__(self, arg=None):
     if hasattr(arg, 'func_globals'):
         ##
         ## Called with a constraint function. 
         ##
         globls = arg.func_globals
         # Names used in constraint
         argv = arg.__code__.co_varnames[:arg.__code__.co_argcount]
         for name in argv:
             if name not in self._names2values:
                 assert name in globls, \
                        "Global name %s not found in function globals" % name
                 self._names2values[name] = globls[name]
         # Gather the range of values of all names used in the constraint
         valuesets = [self._names2values[name] for name in argv]
         self._valueiterator = _itertools.product(*valuesets)
         self._func = arg
         self._funcargnames = argv
         return self
     elif arg is not None:
         ##
         ## Assume called with an iterable set of values
         ##
         arg = frozenset(arg)
         return arg
     else:
         ##
         ## blank call tries to return next solution
         ##
         return self._nextinsearch()

 def _nextinsearch(self):
     arg = self._func
     globls = arg.func_globals
     argv = self._funcargnames
     found = False
     for values in self._valueiterator:
         if arg(*values):
             # Set globals.
             found = True
             for n, v in zip(argv, values):
                 globls[n] = v
             break
     if not found: raise StopIteration
     return values

 def __iter__(self):
     return self

 def next(self):
     return self()

if __name__ == '__main__':
 if True:
     amb = Amb()

     print("\nSmall Pythagorean triples problem:")
     x = amb(range(1,11))
     y = amb(range(1,11))
     z = amb(range(1,11))

     for _dummy in amb( lambda x, y, z: x*x + y*y == z*z ):
         print x, y, z


 if True:
     amb = Amb()

     print("\nRosetta Code Amb problem:")
     w1 = amb(["the", "that", "a"])
     w2 = amb(["frog", "elephant", "thing"])
     w3 = amb(["walked", "treaded", "grows"])
     w4 = amb(["slowly", "quickly"])

     for _dummy in amb( lambda w1, w2, w3, w4: \
                          w1[-1] == w2[0] and \
                          w2[-1] == w3[0] and \
                          w3[-1] == w4[0] ):
         print w1, w2, w3, w4

 if True:
     amb = Amb()

     print("\nAmb problem from "
         "http://www.randomhacks.net/articles/2005/10/11/amb-operator:")
     x = amb([1, 2, 3])
     y = amb([4, 5, 6])

     for _dummy in amb( lambda x, y: x * y != 8 ):
         print x, y

Program output looks like this:


Small Pythagorean triples problem:
3 4 5
4 3 5
6 8 10
8 6 10

Rosetta Code Amb problem:
that thing grows slowly

Amb problem from http://www.randomhacks.net/articles/2005/10/11/amb-operator:
1 4
1 5
1 6
2 5
2 6
3 4
3 5
3 6

So, what you end up with is a declarative style of programming that is good as an exercise, but the solving algorithm is just a brute-force search over all permutations of possible inputs.

Python and EDA" at EuroPython

2010-07-24T08:38:00.003+01:00

I presented on Thursday at EuroPython. from the feedback at the end it seemed to go well.

I've been told that the video will go up on the EuroPython web site by the end of next week which should be "interesting", as it would be the first time that I've watched myself present, (and I don't do many conference presentations).

What I found different was that I am used to presenting to others within the industry who know about RTL, VHDL, digital simulation, coverage, ... To address this Python-savvy audience I re-wrote a presentation on coverage ranking of Verilog/VHDL and cast it in terms of getting coverage of a Python program and ranking the coverage of three Python tests of the Python program. I then walked them through the ranking algorithm for this much simplified case, then tried to show some of the extra work that needs to be done for the real app.

I decided to not show a more polished version of the script. I showed a version where the main functionality could easily be split into more functions for readability; a psyco import does effectively nothing ... The reason was that what I showed was what I did to convince myself that the algorithm was correct and the speed was fast enough. I would have chucked it at this stage if it was slow. (The need was speed of execution as well as an acceptable algorithm for ranking).

The finish was a "winding down" example of just what you can do when you decide to add documentation for a file format in with the program.

I would like to thank my audience for making someone from an industry they may not have known, feel welcome; and the EuroPython organisers, especially John Pinner, for making it happen.

More on binary search

2010-04-25T11:29:00.004+01:00

I wrote this post: Am I one of the 10% of programmers who can write a binary search? as a reply to: Are you one of the 10% of programmers who can write a binary search? by Mike Taylor.

Since then, Mike has followed up with this: Testing is not a substitute for thinking which starts as a summary of the responses to the earlier post then moves on to his own views on testing.

There is mention of some great work by Darius Bacon who tests around twenty solutions gathered, painstakingly, from the comments. Darius' program is great work, but I thought that my testing strategy might be better so modified his code to:

Add my attempt as function: paddy3118
Graft on a version of my test routine to his list of function variants

The delta code was added to the end and looks like the following:


def paddy3118(data, item):
   if not data:
       return -1
   lo, hi = 0, len(data)-1
   mid = (hi + lo) // 2
   #print(lo,mid,hi)
   while item != data[mid] and lo < hi:
       lo, mid, hi = ( (lo, (lo + mid) // 2, mid)
                       if item < data[mid] else
                       (mid+1, (hi + mid + 1) // 2, hi) )
       #print(lo,mid,hi)
   #0/0
   return mid if item == data[mid] else -1
# PASSES test()
#test(paddy3118)

#print sorted(k for k,v in globals().items() if type(v)==type(dan) and 'test' not in k)
tests = (Max, aaron_maxwell, ashish_yadav, ben_gutierrez, clark,
        dan, dave_r, ed_marshall, guilherme_melo, jasper, martin,
        michael_fogleman, paddy3118, patrick_shields, paul,
        rodrigo_b, scott, si, tomas_henriquez,
        travis_wrapped, yc)

if __name__ == '__main__':
   totpass = 0
   for t in tests:
       fails =0
       for dlimit in range(5):
           data = list(range(dlimit))
           for item in data:
               try:
                   assert t(data, item) != -1
               except:
                   fails += 1
                   print("  ERROR: %s FAIL when looking for %i in %r" %
                         (t.func_name, item, data))
                   break
           for item in list(range(dlimit+1)):
               try:
                   assert t(data, item - 0.5) == -1
               except:
                   fails += 1
                   print("  ERROR: %s FAIL as %3.1f is found in %r" %
                         (t.func_name, item - 0.5, data))
                   break
           if fails:
               break
       if not fails:
           totpass += 1
           print('PASS: %s' % t.func_name)
   print("\n%i/%i = %3.1f%% of functions pass" %
         (totpass, len(tests), 100 * totpass / float(len(tests))))
   '''
    PASS: Max
      ERROR: aaron_maxwell FAIL as -0.5 is found in []
    PASS: ashish_yadav
    PASS: ben_gutierrez
      ERROR: clark FAIL when looking for 0 in [0]
      ERROR: dan FAIL as -0.5 is found in []
      ERROR: dave_r FAIL as -0.5 is found in []
      ERROR: ed_marshall FAIL as -0.5 is found in []
      ERROR: guilherme_melo FAIL as -0.5 is found in []
      ERROR: jasper FAIL as -0.5 is found in []
      ERROR: martin FAIL as -0.5 is found in []
    PASS: michael_fogleman
    PASS: paddy3118
      ERROR: patrick_shields FAIL as -0.5 is found in []
    PASS: paul
      ERROR: rodrigo_b FAIL as -0.5 is found in []
      ERROR: scott FAIL as 1.5 is found in [0, 1, 2]
    PASS: si
      ERROR: tomas_henriquez FAIL when looking for 0 in [0, 1]
      ERROR: tomas_henriquez FAIL as 1.5 is found in [0, 1]
      ERROR: travis_wrapped FAIL when looking for 0 in [0]
    PASS: yc

    8/21 = 38.1% of functions pass
    '''

I find that i cannot agree with a lot of Mikes writing on test driven development as he seems to take an extreme position by separating test from development too much.
When I was creating my function, i new that there were likely to be issues with off-by-ones when modifying range indices, and with handling of an empty data list. I chose to leave a 'mental place-holder' on the off-by-one issues as I thought i could easily write tests to exercise those corners and fill them in. The alternative would be to mentally interpret those corner cases to write them down in 'development', then still have to thoroughly test for them anyway.

When writing my tests:

I thought about the division of the range of the indices- and thought that I should use data of both odd and even lengths, and of sufficient length so that a range might be divided more than one.
I thought about all the comparisons - and decided to test for every item that is in the data.
I thought about the significant ways that an item might not be in the data - and decided to test for items between each and every data item, an items below the smallest and above the largest data items.

My paddy3118 function passed the Darius tests. My test routine failed the guilherme_melo test that was passed by Darius.

The only conclusion I can bring, is that you need to do some white box testing, i.e. tests based on knowledge of the code too. Mike Taylor needs to recognise that imposing rigid categories might well lead to extremes. I guess, their is no substitute for experience, but the unexperienced need to be given rules until they can learn enough to break them - Do we teach testing up front? Or not? (Remembering that given half a chance, testing by the novice will be inadequate).

- Paddy.

Am I one of the 10% of programmers who can write a binary search?

2010-04-20T00:58:00.002+01:00

I just read this post where the blogger read a passage in a book that stated that only 10% of programmers could write a binary search from the description. I read up to the description of the binary search algorithm given in the book then thought I would try and implement it without finishing the reading of the rest of his blog entry.

I took around half an hour to come up with the following:

'''
 from: http://reprog.wordpress.com/2010/04/19/are-you-one-of-the-10-percent/
 "Only 10% of programmers can write a binary search"

    Binary search solves the problem [of searching within a pre-sorted
    array] by keeping track of a range within the array in which T
    [i.e. the sought value] must be if it is anywhere in the array.
    Initially, the range is the entire array.  The range is shrunk by
    comparing its middle element to T and discarding half the range.
    The process continues until T is discovered in the array, or until
    the range in which it must lie is known to be empty.  In an
    N-element table, the search uses roughly log(2) N comparisons.

'''

def binsearch(data, item):
 if not data:
     return False
 lo, hi = 0, len(data)-1
 mid = (hi + lo) // 2
 #print(lo,mid,hi)
 while item != data[mid] and lo < hi:
     lo, mid, hi = ( (lo, (lo + mid) // 2, mid)
                     if item < data[mid] else
                     (mid+1, (hi + mid + 1) // 2, hi) )
     #print(lo,mid,hi)
 #0/0
 return item == data[mid]

if __name__ == '__main__':
 for dlimit in range(5):
     data = list(range(dlimit))
     print(data)
     for item in data:
         assert binsearch(data, item)
     for item in list(range(dlimit+1)):
         print(item)
         assert not binsearch(data, item - 0.5)

What happens with an empty list, and the little plus ones were added after I put in the assert statements, and I was satisfied enough to go on and complete the reading of the blog.

Extra Restrictions

The blog goes on to state that the rules of test are that no testing is to be done. A compiler can be used to ensure that their are no mechanical errors however.

Whoa! I was shocked and then secretly pleased. I use Python, and automatically added the tests as part of my development process. In no way did I think that I could do the task without passing some tests. I am still not saying that my implementation is correct, but I have tested for a lot of corner-cases such as:

Lists of length zero, one and more than one.
Items in every position of the list.
Items outside the list: less than, greater than and between all members of the list.

I of course assume that items are comparable, and that data is sorted and subscript-able.

All I can conclude is that if given the test and access to a compiler/interpreter of your choice then one should argue convincingly for testing any submission as you would never deign to submit any code as complete unless tested. If the testers remain unconvinced then it may be time to stick to your guns as you don't gain by lowering your standards, especially in an interview situation when you know that testing is right!

P.S. I have since had a look at the code of the Python bisect module and the algorithm given on Wikipedia, and they both seem more elegant than mine.

The Securities and Exchange Commission is to mandate the use of Python?

2010-04-16T21:14:00.001+01:00

Read on

Smartphone Programming, Compare and Contrast

2010-04-10T14:14:00.002+01:00

Apple: Steve Jobs Has Just Gone Mad
Android: Android Scripting Environment brings scripting languages to Android.
Nokia: Python for S60

Do you still want that Ipad?

Why Python?

2010-03-30T20:22:00.002+01:00

Choose for yourself. Have a look at the examples on the Rosetta Code site. You might like to follow the links to a few pages that I did the task descriptions for as they tend to be meatier tasks, but the examples in different languages should allow you to compare how different languages are used to solve the same problems.

Unfortunately it cannot eliminate the issue of the competence of the individual programmers, but it might provide some useful info.

Have fun :-)

Expressions to Truth-Tables

2010-03-16T20:43:00.006Z

I had a need for the data in the truth tables for a standard cell library, rather like that shown here if you look at say page 35 you will see that for for every Cell it gives the cell name; the cell function as a boolean expression, and the truth-table for the cell.

My cell library was an Infineon one and I had the pdf version of the databook. A quick squirt through pdftotext scrambled the truth table, but did preserve the boolean expression of each cells' function, and the cell name.

I wrote a Python script that used a multi-line regexp to extract each cells info including the boolean equation.

The equation from the Infineon databook uses +,*, and ! for boolean or, and and not operators, but also used lots of parentheses so that operator precedence did not affect the value of the infix expression. I realised that I could substitute for the Python bitwise operators, and form a python expression that could then be evaluated with appropriate values for any input names to create a truth table.

I have left out the parsing and present a function to create a truth-table from a boolean expression:


from itertools import product

def truthtable(expressionstring):
  comp = compile(expressionstring.strip(), '-', 'eval')
  inputs = sorted(comp.co_names)
  table = {}
  for values in product( *([(0,1)] * len(inputs)) ):
      valuedict = dict(zip(inputs, values))
      answer = eval(comp, valuedict)
      table[values] = answer
  return expressionstring.strip(), inputs, table

def pptable(ex, inputs, table):
  print("\nTRUTH-TABLE FOR EXPRESSION: %r\n" % ex)
  print(' ' + ' '.join('%2s' % inp for inp in inputs) + ': OUTPUT')
  fmt = '  %i' * len(inputs) + ':  %i'
  for inp, ans in sorted(table.items()):
      print(fmt % tuple(list(inp) + [ans]))

ex, inputs, table = truthtable('(A & (~C)) | D')
pptable(ex, inputs, table)

Going through function truthtable:

compile() gives access to all the names in the compiled expression via attribute .co_names - this is all the names used in the expression.
product() allows you to quickly generate the 1/0 value combinations for the input names.
valuedict assigns a value to each name used in the expression.
eval evaluates the expression with the given values on each name

The output of the program is:


TRUTH-TABLE FOR EXPRESSION: '(A & (~C)) | D'

  A  C  D: OUTPUT
  0  0  0:  0
  0  0  1:  1
  0  1  0:  0
  0  1  1:  1
  1  0  0:  1
  1  0  1:  1
  1  1  0:  0
  1  1  1:  1

- No separate parser required!

NPR 2002 Puzzle

2010-03-04T22:59:00.006Z

A puzzle I copied from here via reddit:
1 2 3 4 5 6 7 8 9 = 2002
Put any combination of plus, times, and spaces between the digits on the left to make the identify true.
The answer was easy and was quick to compute.

I realised that you could evaluate all expressions where the numbers 1 through nine are separated by the permutations of either a null string, '', '+' or '*'

The answer came from the following:


>>> from itertools import product
>>> eqn = '%s'.join(list('123456789')) + '==2002'
>>> eqn
'1%s2%s3%s4%s5%s6%s7%s8%s9==2002'
>>> for ops in product(*[['', '+', '*']]*8):
    if eval(eqn % ops):
        print(eqn % ops)

        
1*23+45*6*7+89==2002
1*2+34*56+7+89==2002
>>>

So, there are two solutions.

What is a first-class function?

2009-12-03T04:44:00.002Z

There is a debate about it on Wikipedia, and I am curious:

Does a language have first-class functions if it can use its functions to solve particular problems in a certain manner, or
Does a language have first-class functions if functions are just as flexible as other types in the language itself?

It probably means 1. above, as the second definition would omit most static languages, as they might be able, for example, to take textual input and change it into an internal representation of an integer, but not do the same for a function with similar effort. Someone who exclusively thinks in terms of static languages may not have realized a difference.

Hamming Numbers (and Perl: Why I just ignored it)

2009-11-28T08:38:00.001Z

I had read the following post: "Perl: Love it, or hate it, but don’t ignore
it.", by Phillip Smith, in which the author wishes Perl
solutions were credited as much as solutions in other languages such as
Python and Ruby.

I agreed that Perl wasn't mentioned as much as it was, and a day or so
later I came across an example of why....

I monitor whats happening in Python, and came across a new recipe added
to Activestate code called: "Recipe 576961: Technique for cyclical
iteration", by Raymond Hettinger. On first, second, and third
readings, i just couldn't understand the need for the deferred_output
function; the whole recipe remained opaque:

Raymond's
code:

from itertools import tee, chain, islice
from heapq import merge

def hamming_numbers():
   # Generate "5-smooth" numbers, also called "Hamming numbers"
   # or "Regular numbers". See: http://en.wikipedia.org/wiki/Regular_number
   # Finds solutions to 2**i * 3**j * 5**k for some integers i, j, and k.

   def no_repeats(it):
       last = None
       for i in it:
           if i is not last:
                yield i
                last = i

    def deferred_output():
       for i in output:
           yield i

   result, p2, p3, p5 = tee(deferred_output(), 4)
   m2 = (2*x for x in p2) # multiples of 2
   m3 = (3*x for x in p3) # multiples of 3
   m5 = (5*x for x in p5) # multiples of 5
   merged = merge(m2, m3, m5)
   combined = chain([1], merged) # prepend a starting point
   output = no_repeats(combined) # eliminate duplicates

   return result

if __name__ == '__main__':
   print(list(islice(hamming_numbers(), 1000))) # Show first 1000 hamming numbers

I was hooked by the series however and so went looking for more
information.

Perl

I came across a mention to the following page: "Infinite
lists in Perl" by Mark Dominus which is copyright 1997 by the
Perl Journal. Normally I wouldn't mention the article, as it spends
more time creating an implementation of generators in Perl as it does
on addressing the computation of Hamming numbers. It also uses the term
"stream" instead of generators. but maybe this is down to the age of
the article, i.e. 1997?

Haskell

Looking a bit further, I found this Dr Dobbs Article, and this
commented Haskell
solution:

hamming = 1 : map (2*) hamming `merge` map (3*) hamming `merge` map (5*) hamming
    where merge (x:xs) (y:ys)
           | x < y = x : xs `merge` (y:ys)
           | x > y = y : (x:xs) `merge` ys
           | otherwise = x : xs `merge` ys

The Haskell solution lead me to think that maybe Python has an issue
with defining recursive generators, and some experimentation lead me to
believe that yep, maybe that was why Raymond has the deferred_output
function.

My Python

Since I have trouble understanding Raymond s I went away and wrote my
own solution. It works. It can generate very large values of the
sequence quite quickly. But it does need memory to store
multiples:

from itertools import islice
from pprint import pprint as pp


def regularnumbergen(start = 1, multipliers=(2, 3, 5)):
   '''\
    Generates the Regular Number sequence defined here:
        http://en.wikipedia.org/wiki/Regular_number
    The sequence starts with 1 and consists of all multiples by 2, 3, or 5
    of any number in the sequence in non-repeating, increasing  order.
    '''
   multiples = [[] for i in multipliers]
   this = start
   while True:
       yield this
       # get multiples
       for multiple, multiplier in zip(multiples, multipliers):
           multiple.append(this * multiplier)
       # get next
       this = min(m[0] for m in multiples)
       # remove this from multiples
       for multiple in multiples:
           if multiple[0] == this:
               del multiple[0]

if __name__ == '__main__':
   print(list(islice(regularnumbergen(), 10)))

Generating the following large values of the sequence took only seconds:

>>> print(list(islice(regularnumbergen(), 99999, 100001)))
[290142196707511001929482240000000000000, 290237644800000000000000000000000000000]
>>>

Wrapup

Well, I've learnt a bit more about generators, the limitatins in their
Python implementation w.r.t. Haskel, and I've put more meat on why I
tend to not quote much Perl In this case, i couldn't see the wood for
the trees! (But I do quote Perl).

Curious about Go.

2009-11-12T08:20:00.003Z

Verilog is to SytemVerilog as C is to C++ .
Verilog is to __________ as C is to Go ?

My main criticism of SytemVerilog was the size of the language compared to Verilog. Maybe we'll get SystemGo in five years time. (Or maybe Go itself will warp into some monstrosity; my crystal ball needs a new backlight).

More Words From Hex Letters (with extra zero's)

2009-11-04T18:55:00.004Z

STOP PRESS! Someone commented that zeros could be used as ones so an update is at the bottom.

It was pointed out to me that I needn't restrict myself to four character words, and could use the number '1' as letter 'l', and the number '5' as the letter 'S'.

So here is the program and output of words of four or more characters that can be made from hexadecimal characters.


import urllib

words = urllib.urlopen('http://www.puzzlers.org/pub/wordlists/unixdict.txt'
                  ).read().split()
#len(words) == 25104

def findwords(hexchars):
hexchars = set(hexchars)
hexwords = [w for w in words if 4 <= len(w) and all(ch in hexchars for ch in w)]
return hexwords

def printwords(hexwords):
for i in range(0, len(hexwords), 5):
   print ('  ' + ' '.join("%-8s" % w for w in hexwords[i:i+5]))

hexchars = 'abcdef'
#hexwords = findwords(hexchars)
#print "\nHEXWORDS"
#printwords(hexwords)

# Now to use one: 1 as ell: l; and five: 5 as ess: s
hexwords = findwords(hexchars + 'ls')
hexwords = [word.replace('s','5').replace('l','1') for word in hexwords]
print "\nHEXWORDS USING one: 1 AS ell: l; AND five: 5 AS ess: s"
printwords(hexwords)

The list of words found have some gems, such as babble, for when you know your generating data that no one will read; baffle, baseball, blab, bleed, bless, (and blessed if it were in the dictionary), cabal, cascade, cease, debacle, decade, decease, (and should have deceased), false, flee, label, ...

(Someone at work suggested we use 1abe1 or label to prefix a section of writes to a peripheral in an assembler prog which got me going).

The words:


HEXWORDS USING one: 1 AS ell: l; AND five: 5 AS ess: s
aaa5     ababa    aba5e    abba5    abbe
abed     abe1     ab1e     ab5ce55  accede
acce55   added    add1e    ade1e    affab1e
a1ba     a1ec     a1fa1fa  a11e1e   a55e55
babb1e   babe     babe1    bade     baff1e
ba1d     ba1e     ba11     ba11ad   ba11ed
ba15a    ba5a1    ba5e     ba5eba11 ba5e1
ba55     bead     bead1e   beebe    beef
befa11   befe11   be1a     be11     be11a
be11e    be55     be55e1   b1ab     b1ade
b1ed     b1eed    b1e55    caba1    cab1e
cafe     ca1eb    ca1f     ca11     ca11a
ca5cade  ca5e     cea5e    cede     ce1ebe5
ce11     c1ad     c1a55    dabb1e   dacca
dada     dade     da1e     da11a5   dead
deaf     dea1     debac1e  deba5e   decade
deca1    decca    decea5e  deed     deface
de11     de11a    ea5e     ea5e1    ecc1e5
efface   effaceab1e e1ba     e11a     e15e
fab1e    facade   face     fade     fa11
fa15e    feeb1e   feed     fee1     fe11
f1ea     f1ed     f1ee     f1eece   1abe1
1ace     1ad1e    1a5e     1a55     1ead
1eaf     1ea5e    1eed5    1e55     1e55ee
5ab1e    5accade  5add1e   5afe     5a1ad
5a1e     5a11e    5cab     5ca1a    5ca1d
5ca1e    5ea1     5ecede   5eeab1e  5eed
5eedbed  5e1f     5e11     51ab     51ed

Additional Words with an 'o'

Without Michaels comment I would have skipped all words that use zero for an o, so I made the following two changes:

hexwords = findwords(hexchars + 'lso')
hexwords = [word.replace('s','5').replace('l','1').replace('o', '0') for word in hexwords]

And created the following extra words:

  ab0de    acc01ade ad0be    a1c0a    a110cab1e
a10e     a100f    a150     ba1b0a   ba550 
b10b     b10c     b100d    b0bb1e   b0ca  
b0de     b01d     b01dface b01e     b010  
b05e     b055     caca0    c10d     c105e 
c0a1     c0a1e5ce c0bb     c0bb1e   c0b01 
c0ca     c0c0     c0c0a    c0da     c0dd1e
c0de     c0ed     c0ffee   c01a     c01d  
c01e     c01055a1 c001     c05ec    dec0de
d0bb5    d0dd     d0d0     d0ff     d01ce 
d01e     d011     d00d1e   d05e     ec01e 
ee0c     fa110ff  f10c     f10e     f100d 
f0a1     f0ca1    f01d     f00d     f001  
f055     1a05     1a550    10ad     10af  
10be     10b0     10ca1    10ca1e   10eb  
10e55    101a     1011     1005e    1005e1eaf
105ab1e  105e     1055     0a5e5    0be5e 
0b0e     0b5e55   0de55a   0ffa1    0ff10ad
0ff5add1e 01af     00d1e5   0510     5caff01d
5c0ff    5c01d    5eaf00d  510b     510e  
50da     50fa     501ace   501d     501e  
5010

The 24 game

2009-11-01T07:53:00.003Z

I've just created a new task over on Rosetta Code for the 24 game, and a separate task to create a solver for the game.

It's a brute -force solver, but I find the game to be playable, so will get my sprogs to test their arithmetic skills with it later.

I've just realised that the solver could be easily modified to solve "The numbers game" from the UK Channel 4 program Countdown. (Although their could be an issue with division as, when given the digits 3 3 7 7, the solver finds the solution: ( 3 + 3 / 7 ) * 7 which may not be allowed by the countdown rules which state that "only whole numbers may be used").

Template for forking a Python program

2009-10-18T20:48:00.001+01:00

I am doing a large regression run of some 30,000 simulations producing aroung 30 Gigs of logs to extract test parameter and simulation result from each.

I have a program that goes through each log in turn extracting results and this may take tens of minutes to complete.

A future optimisation will be to extract pass/fail data for each regression run in the same job that ran the simulation, but at present, there may be changes to the pass/fail criteria, so it is run as a separate task after the regression simulations are all done.

I normally run the log extraction process on an 8core, 16 thread machine with a fast connection to the file server, so was thinking of ways to parallelise the task.

Enter this article, part of Doug Hellmann's PyMOTW series, about forking a process.

After running it, I decided to expand on it to be more like my situation where you have a list of similar tasks to perform and want to split it between multiple processes via the Posix fork() call. The following is just a framework - I hope to get time to acctually use it and see what the problems with it are.

 1 import os
 2 import sys
 3 import time
 4 
 5 maxprocs = 4 # procs - 1 children and the parent to do the work
 6 jobarray = range(25) # Work to be split between proccesses
 7 
 8 def do1job(job):
 9  'Process one item from the jobarray'
10  return job
11 
12 def picklesave(results):
13  'save partialresults to file'
14  time.sleep(3)
15 
16 def unpickleread(proc):
17  'read partial results previousely saved by process number procc'
18  return [] # dummy
19 
20 def dowork(proc, maxprocs=maxprocs, jobarray=jobarray):
21  ' Split jobarray tasks between maxprocs by proci number'
22  time.sleep(proc) # convenience processing time 
23  results = []
24  for count, job in enumerate(jobarray):
25  if count % maxprocs == proc:
26  results.append(do1job(job))
27  picklesave(results)
28  print ("Process: %i completed jobs: %s" % (proc, results))
29  return results
30 
31 
32 def createworkerchildren(maxprocs=maxprocs):
33  for proc in range(1, maxprocs):
34  print 'PARENT: Forking %s' % proc
35  worker_pid = os.fork()
36  if not worker_pid:
37  'This is a child process'
38  dowork(proc)
39  # Children exit here!
40  sys.exit(proc)
41 
42 # Start and don't wait for child proccesses to do their work
43 createworkerchildren()
44 
45 # Do parent processes share of the work
46 results = []
47 results += dowork(0) # don't have to pickle/unpickle Parent processes share
48 
49 # wait for children
50 for proc in range(1, maxprocs):
51  print 'PARENT: Waiting for any child'
52  done = os.wait()
53  print 'PARENT: got child', done
54 
55 # read and integrate child results
56 for proc in range(1, maxprocs):
57  results += unpickleread(proc)
58 
59 # Calculate and print data summary
60 pass
61

Running the above on cygwin produced:

bash$ python testfork.py
PARENT: Forking 1
PARENT: Forking 2
PARENT: Forking 3
Process: 0 completed jobs: [0, 4, 8, 12, 16, 20, 24]
PARENT: Waiting for any child
Process: 1 completed jobs: [1, 5, 9, 13, 17, 21]
PARENT: got child (5536, 256)
PARENT: Waiting for any child
Process: 2 completed jobs: [2, 6, 10, 14, 18, 22]
PARENT: got child (3724, 512)
PARENT: Waiting for any child
Process: 3 completed jobs: [3, 7, 11, 15, 19, 23]
PARENT: got child (4236, 768)
bash$

Notice how each child gets to do a fair poition of the jobs, (assuming all jobs need the same resources).

Extend functools.partial to nested function calls?

2009-10-11T10:18:00.004+01:00

I am reading the 1977 ACM Turing Award Lecture "Can Programming Be Liberated from the von Neumann Style? A Functional Style and Its Algebra of Programs" by John Backus.

In section 5.1 he is contrasting a von Neumann program for inner product, which he gives as:

c := 0
for i := 1 step 1 until n do
   c := c + c[i] x b[i]

The above converts quit nicely to the python equivalent:

>>> v1, v2 = [1,2,3], [6,5,4]
>>> def ip_von_neumann(a, b):
    c = 0
    for i in range(3):
       c = c  + a[i] * b[i]
    return c

>>> ip_von_neumann(v1, v2)
28

Backus then goes on to contrast the von Neumann style with his functional program:

Def Innerproduct ≡ (Insert +)o(ApplyToAll x)oTranspose

Which I converted to the following functional style Python:

>>> from operator import add, mul
>>> from functools import reduce as insert
>>>
>>> def ip_functional(*x):
    return insert(add, map(mul, *x))

>>> ip_functional(v1, v2)
28

Well, the Python works, but it still has this name x that I would like to eliminate.

I think I can go part of the way by using functools.partial. For example, for the single function call of map(mul, *x) I can do:

>>> from functools import partial
>>> p1 = partial(map, mul)
>>> assert p1(v1, v2) == map(mul, *(v1, v2))
>>> map(mul, *(v1, v2))
[6, 10, 12]
>>>

But I don't know how to go that extra step and remove x from the definition fully, ie end up with just some assignment like:

ip_purefunctional = {A pure function of functions without any argument placeholders}

All help would be appreciated, Thanks.

The Otter Algorithm

2009-10-02T23:02:00.000+01:00

I was reading about a Richard Dawkins program that through random mutation and selection of best fit, would gradually converge a random string into a target string. It is called the Weasel Program,and a rough outline of it is:

Start with a random string of 28 characters
Copy this string 100 times, with a 5% chance per character of that character being replaced with a random character
Compare each new string with the target "METHINKS IT IS LIKE A WEASEL", and give each a score
If any of the new strings has a perfect score, halt
Otherwise, take the highest scoring string, and go to step 2

After coding it up, I found it just would not converge to a final match for me. I suspected I had a mistake somewhere , but just did not like how long it was taking to converge.

I decided to change things.

It seemed to me that if you are further away from the ideal, then you should change more letters; and, after seeing the score for how good the current best string is sometimes decrease, I made a modification to choose the best of the copies and the best string that they are mutations of.

With that modification, convergence was swift.

Oh, I call it the Otter program, as it is more than the Weasel.

The prog

 1 
 2 '''
 3 From: http://en.wikipedia.org/wiki/Weasel_program
 4       With modifications
 5 '''
 6 
 7 from string import ascii_uppercase
 8 from random import choice, random
 9 
10 target = list("METHINKS IT IS LIKE A WEASEL")
11 chance = .05
12 times  = 100
13 
14 charset      = ascii_uppercase + ' '
15 iterations   = 0
16 highscore    = [choice(charset) for _ in range(len(target))]
17 perfectscore = len(target)
18 
19 def score(trial):
20     return sum(t==h for t,h in zip(trial, target))
21 
22 while highscore != target:
23     thischance =  1-((perfectscore - score(highscore)) / perfectscore * (1 - chance))
24     #thischance = chance
25     iterations += 1
26     if iterations % 50 == 0:
27         print (iterations, score(highscore), ''.join(highscore))
28     copies = [ [(ch if random() <= thischance else choice(charset)) for ch in highscore]
29                for _ in range(times) ]  + [highscore]
30     highscore = max(copies, key=score)
31 print (iterations, ''.join(highscore))

Function score at line 19 uses the fact that booleans True and False are also numbers and calculates how many characters are in the right position as its result.

Line 23 modifies the amount of randomness injected, based on how bad the last highscore was: Change more if you are further away from the target - less as you approach.
Name thischance is used in line 28. random returns a float between 0.0 and 1.0. if it is less than thischance then a particular character is not changed.

The function for thischance in line 23 is equal to chance when no characters match, and 1.0 if all characters were to match. This forces less change as you approach the target string.

Here is a little snippet to show how thischance varies with the score (sc):

>>> for sc in range(0,29,7):
        print( "score=%2i, thischance=%5.2f" % (sc,
                1-((perfectscore - sc) / perfectscore * (1 - chance))) )


score= 0, thischance= 0.05
score= 7, thischance= 0.29
score=14, thischance= 0.53
score=21, thischance= 0.76
score=28, thischance= 1.00

In line 29 the addition of highscore into copies means that the max taken in line 30 can never decrease.

P.S. The post is only tangentially linked to Darwinism. I am secure in the knowledge that my children's school teaches it, and I am not really interested in it being mentioned in any comments. Treat it just as a way to go from a random choice of letters and home in on a target string by repeated random alterations, and selection of the "best".

Natural Sorting

2009-09-24T06:40:00.000+01:00

Natural sorting orders items in some way that is more than just the
order dictated by computer character codes. It is called natural
because the new ordering is supposed to mean more to a human consumer
than an ordering on base character codes .

I guess, what is becoming the canonical example is the ordering of
numbered files in a directory listing. A normal sorted directory
listing might show:

style="font-family: monospace;">bash$ for (( i=0;
i<111;i+=5)); do touch f$i.txt; done style="font-family: monospace;">
bash$ ls -1 style="font-family: monospace;">
f0.txt style="font-family: monospace;">
f10.txt style="font-family: monospace;">
f100.txt style="font-family: monospace;">
f105.txt style="font-family: monospace;">
f110.txt style="font-family: monospace;">
f15.txt style="font-family: monospace;">
f20.txt style="font-family: monospace;">
f25.txt style="font-family: monospace;">
f30.txt style="font-family: monospace;">
f35.txt style="font-family: monospace;">
f40.txt style="font-family: monospace;">
f45.txt style="font-family: monospace;">
f5.txt style="font-family: monospace;">
f50.txt style="font-family: monospace;">
f55.txt style="font-family: monospace;">
f60.txt style="font-family: monospace;">
f65.txt style="font-family: monospace;">
f70.txt style="font-family: monospace;">
f75.txt style="font-family: monospace;">
f80.txt style="font-family: monospace;">
f85.txt style="font-family: monospace;">
f90.txt style="font-family: monospace;">
f95.txt

It might be more friendly if the files were ordered so that the numeric
section increased steadily.

In other situations, such as in Movie listings, you sometimes see
common words such as "the" ignored for sorting purposes (or even
transferred to the end of the item so "The Silence of the Lambs" would
be listed as "Silence of the Lambs, The")

Ignoring case, and treating white-space characters as the same - and
multiple strings of whitespace characters as a single whitespace, might
also be useful.

Since my initial tentative steps into working with non-ASCII characters
I tried to consider sorting accented characters as if the accent was
not present; ligatures as separate characters, and made a stab at
tackling what I call 'replacements' - where certain
characters, such as the German href="http://en.wikipedia.org/wiki/Scharfes_S">scharfes S
might be sorted as SS,.

I had a quick stab at creating a natural sort routine and
added several natural sorting abilities. They might not all be useful
at once, but it should be straight forward to pick and choose what you
want n your natural sort..

As always, when working with non-ASCII character sets, I find it
difficult to translate my working multi-byte Python program to HTML,
but there should be enough surrounding information to ensure that any
multi-byte character can be worked-out.

  1 
 color="#804040">  2 # -*- coding: utf-8 -*-
 color="#804040">  3 # Not Python 3.x (Can't compare str and int)
 color="#804040">  4 
 color="#804040">  5 '''
 color="#804040">  6     Natural sorting: Sorting of text that does more than rely on the
 color="#804040">  7     order of individual characters codes to make the finding of
 color="#804040">  8     individual strings easier for a human reader.
 color="#804040">  9     
 color="#804040"> 10     Some 'natural' orderings might include:
 color="#804040"> 11         1) Ignore leading, trailing and multiple adjacent spaces
 color="#804040"> 12         2) All space types equivalent.
 color="#804040"> 13         3) Sorting without regard to case.
 color="#804040"> 14         4) Sorting numeric portions of strings in numeric order:
 color="#804040"> 15             foo9.txt before foo10.txt
 color="#804040"> 16             As well as ... x9y99 before x9y100, before x10y0
 color="#804040"> 17             ... (for any number of groups of integers in a string).
 color="#804040"> 18         5) Title sorts: without regard to a leading, very common, word such
 color="#804040"> 19            as 'The' in "The thirty-nine steps".
 color="#804040"> 20         6) Sort letters without regard to accents.
 color="#804040"> 21         7) Sort ligatures as separate letters.
 color="#804040"> 22         8) Replacements:
 color="#804040"> 23             Sort german scharfes S (ÃŸ) as ss
 color="#804040"> 24             Sort Å¿,  LATIN SMALL LETTER LONG S as s
 color="#804040"> 25             Sort Ê’,  LATIN SMALL LETTER EZH as s
 color="#804040"> 26             ...
 color="#804040"> 27 '''
 color="#804040"> 28 
 color="#804040"> 29 from itertools  color="#a020f0">import groupby
 color="#804040"> 30 from unicodedata  color="#a020f0">import decomposition, name
 color="#804040"> 31 
 color="#804040"> 32 commonleaders = [' color="#ff00ff">the'] # lowercase leading words to ignore
 color="#804040"> 33 replacements = {' color="#ff00ff">ÃŸ': 'ss',   color="#0000ff"># Map single char to replacement string
 color="#804040"> 34                 ' color="#ff00ff">Å¿': 's',
 color="#804040"> 35                 ' color="#ff00ff">Ê’': 's',
 color="#804040"> 36                 }
 color="#804040"> 37 
 color="#804040"> 38 hexdigits = set(' color="#ff00ff">0123456789abcdef')
 color="#804040"> 39 decdigits = set(' color="#ff00ff">0123456789')    color="#0000ff"># Don't use str.isnumeric
 color="#804040"> 40 
 color="#804040"> 41 def  color="#008080">splitchar(c):
 color="#804040"> 42     ' color="#ff00ff"> De-ligature. De-accent a char'
 color="#804040"> 43     de = decomposition(c)
 color="#804040"> 44     if de:
 color="#804040"> 45          color="#0000ff"># Just the words that are also hex numbers
 color="#804040"> 46         de = [d  color="#804040">for d  color="#804040">in de.split()
 color="#804040"> 47                    color="#804040">if all(c.lower()
 color="#804040"> 48                           color="#804040">in hexdigits  color="#804040">for c  color="#804040">in d)]
 color="#804040"> 49         n = name(c, c).upper()
 color="#804040"> 50          color="#0000ff"># (Gosh it's onerous)
 color="#804040"> 51          color="#804040">if len(de)> 1  color="#804040">and ' color="#ff00ff">PRECEDE'  color="#804040">in n:
 color="#804040"> 52              color="#0000ff"># E.g. Å‰  LATIN SMALL LETTER N PRECEDED BY APOSTROPHE
 color="#804040"> 53             de[1], de[0] = de[0], de[1]
 color="#804040"> 54         tmp = [ unichr(int(k, 16))  color="#804040">for k  color="#804040">in de]
 color="#804040"> 55         base, others = tmp[0], tmp[1:]
 color="#804040"> 56          color="#804040">if ' color="#ff00ff">LIGATURE'  color="#804040">in n:
 color="#804040"> 57              color="#0000ff"># Assume two character ligature
 color="#804040"> 58             base += others.pop(0)
 color="#804040"> 59     else:
 color="#804040"> 60         base = c
 color="#804040"> 61     return base
 color="#804040"> 62 
 color="#804040"> 63 
 color="#804040"> 64 def  color="#008080">sortkeygen(s):
 color="#804040"> 65     ''' color="#ff00ff">Generate 'natural' sort key for s
 color="#804040"> 66 
 color="#804040"> 67     Doctests:
 color="#804040"> 68         >>> sortkeygen('  some extra    spaces  ')
 color="#804040"> 69         [u'some extra spaces']
 color="#804040"> 70         >>> sortkeygen('CasE InseNsItIve')
 color="#804040"> 71         [u'case insensitive']
 color="#804040"> 72         >>> sortkeygen('The Wind in the Willows')
 color="#804040"> 73         [u'wind in the willows']
 color="#804040"> 74         >>> sortkeygen(u' color="#6a5acd">\462 ligature')
 color="#804040"> 75         [u'ij ligature']
 color="#804040"> 76         >>> sortkeygen(u' color="#6a5acd">\335\375 upper/lower case Y with acute accent')
 color="#804040"> 77         [u'yy upper/lower case y with acute accent']
 color="#804040"> 78         >>> sortkeygen('foo9.txt')
 color="#804040"> 79         [u'foo', 9, u'.txt']
 color="#804040"> 80         >>> sortkeygen('x9y99')
 color="#804040"> 81         [u'x', 9, u'y', 99]
 color="#804040"> 82     '''
 color="#804040"> 83     # Ignore leading and trailing spaces
 color="#804040"> 84     s = unicode(s).strip()
 color="#804040"> 85     # All space types are equivalent
 color="#804040"> 86     s = ' color="#ff00ff"> '.join(s.split())
 color="#804040"> 87     # case insentsitive
 color="#804040"> 88     s = s.lower()
 color="#804040"> 89     # Title
 color="#804040"> 90     words = s.split()
 color="#804040"> 91     if len(words) > 1  color="#804040">and words[0]  color="#804040">in commonleaders:
 color="#804040"> 92         s = ' color="#ff00ff"> '.join( words[1:])
 color="#804040"> 93     # accent and ligatures
 color="#804040"> 94     s = ''.join(splitchar(c)  color="#804040">for c  color="#804040">in s)
 color="#804040"> 95     # Replacements (single char replaced by one or more)
 color="#804040"> 96     s = ''.join( replacements.get(ch, ch)  color="#804040">for ch  color="#804040">in s )
 color="#804040"> 97     # Numeric sections as numerics
 color="#804040"> 98     s = [ int("".join(g))  color="#804040">if isinteger  color="#804040">else "".join(g)
 color="#804040"> 99            color="#804040">for isinteger,g  color="#804040">in groupby(s,  color="#804040">lambda x: x  color="#804040">in decdigits)]
 color="#804040">100 
 color="#804040">101     return s
 color="#804040">102 
 color="#804040">103 def  color="#008080">naturalsort(items):
 color="#804040">104     ''' color="#ff00ff"> Naturally sort a series of strings
 color="#804040">105 
 color="#804040">106     Doctests:
 color="#804040">107         >>> naturalsort(['The Wind in the Willows','The 40th step more',
 color="#804040">108                          'The 39 steps', 'Wanda'])
 color="#804040">109         ['The 39 steps', 'The 40th step more', 'Wanda', 'The Wind in the Willows']
 color="#804040">110 
 color="#804040">111     '''
 color="#804040">112     return sorted(items, key=sortkeygen)
 color="#804040">113 
 color="#804040">114 
 color="#804040">115 



 color="#804040">116 # Extras...
 color="#804040">117 
 color="#804040">118 '''
 color="#804040">119 Help on built-in function decomposition in module unicodedata:
 color="#804040">120 
 color="#804040">121 decomposition(...)
 color="#804040">122     decomposition(unichr)
 color="#804040">123     
 color="#804040">124     Returns the character decomposition mapping assigned to the Unicode
 color="#804040">125     character unichr as string. An empty string is returned in case no
 color="#804040">126     such mapping is defined.
 color="#804040">127 
 color="#804040">128 Help on built-in function name in module unicodedata:
 color="#804040">129 
 color="#804040">130 name(...)
 color="#804040">131     name(unichr[, default])
 color="#804040">132     Returns the name assigned to the Unicode character unichr as a
 color="#804040">133     string. If no name is defined, default is returned, or, if not
 color="#804040">134     given, ValueError is raised.
 color="#804040">135 '''
 color="#804040">136 
 color="#804040">137 
 color="#804040">138 def  color="#008080">_temp():
 color="#804040">139     # for c in u'\370\157\117\463\462\511\733\337\337':
 color="#804040">140     for c  color="#804040">in u' class="Constant">Ã¸oOÄ³Ä²Å‰ÇÃŸÃÃ½':
 color="#804040">141         de = decomposition(c)
 color="#804040">142          color="#804040">if de:
 color="#804040">143             de = [d  color="#804040">for d  color="#804040">in de.split()
 color="#804040">144                        color="#804040">if all(c.upper()
 color="#804040">145                           color="#804040">in ' color="#ff00ff">0123456789ABCDEF'  color="#804040">for c  color="#804040">in d)]
 color="#804040">146             n = name(c, c).upper()
 color="#804040">147              color="#804040">if len(de)> 1  color="#804040">and ' color="#ff00ff">PRECEDE'  color="#804040">in n:
 color="#804040">148                  color="#0000ff"># E.g. Å‰  LATIN SMALL LETTER N PRECEDED BY APOSTROPHE
 color="#804040">149                 de[1], de[0] = de[0], de[1]
 color="#804040">150             tmp = [ unichr(int(k, 16))  color="#804040">for k  color="#804040">in de]
 color="#804040">151             base, others = tmp[0], tmp[1:]
 color="#804040">152          color="#804040">else:
 color="#804040">153             base, others = c, []
 color="#804040">154          color="#804040">print(" color="#ff00ff">%s %o %s\n color="#ff00ff">  %s  %s" %
 color="#804040">155                   (c,ord(c), name(c,' color="#ff00ff">-'),
 color="#804040">156                    name(base, base),
 color="#804040">157                    ' color="#ff00ff"> '.join(name(o,o)  color="#804040">for o  color="#804040">in others)))
 color="#804040">158 
 color="#804040">159 
 color="#804040">160

Easy Command-line Parallelism

2009-09-14T00:11:00.001+01:00

I am trying out a dual chip Xeon server with 8 cores that shows up as
16 CPU's to Redhat Linux. I have my href="http://paddy3118.blogspot.com/search?q=process">process
runner script, that I use with a list of 100000 simulations
to do overnight and have experimented with running different numbers of
simulations in parallel to get good throughput.

The process runner is good, I still use it, but I also wanted something
for more general command-line use. for example, I have around 750
directories, each containing files relating to one test. I organize the
directories to have regular names, and am used to using shell for loops
to process the dataset.

For example, to gzip all the */transcript files I would do:

foreach f (*/transcript)
    gzip $f
end

(Yep, I know it's cshell, but that is the 'standard' at work.)

This isn't making much use of the processing power at hand, so I went
googling for a lightweight way to run such tasks in parallel.

I found two methods. One, to use make and its -P option to run several
jobs in parallel would need the creation of a make file in each case,
which is onerous. The other, which was new to me, rather than a
forgotten feature, is that xargs has a -P option that does a similar
thing - in combination with other options it will run up to a
given number of jobs in parallel:

/bin/ls */transcript \
  | xargs -n1 -iX -P16 gzip X

The above lists all the files, pipes them to xargs which takes them one
at a time (-n1) to form a job by substituting all occurrences of
character X (-iX) in its arguments after its options. a maximum of 16
(-P16) jobs are then run at any one time.

When I want to simulate all 750 tests and create a coverage listing of
each test I use xargs to call the bash shell and write something like:

/bin/ls -d * \
 | xargs -n1 -iX -P16 bash \
 -c 'cd X && vsim tb -c -coverage -do "coverage save -onexit -code t coverX.ucdb;do ../run.do" 2>&1 >/dev/null && vcover report -file coverX.rpt coverX.ucdb 2>&1 >/dev/null'

Yep, I do create such long monstrosities at times. It's because I am
exploring new avenues at the moment, and things are likely to change,
plus I need to know, and change, the details quite frequently. Although
i do make a note of these long one-liners, (and keep a lot of history
in my shells), I tend to only wrap things in a (bash), script when it
becomes stable.

Note:

Following the links from the Wikipedia page, it seems that the -P
option to xargs is not mandated by href="http://www.opengroup.org/onlinepubs/9699919799/utilities/xargs.html">The
Open Group, and does not appear in href="http://docs.sun.com/app/docs/doc/816-5165/xargs-1?a=view">Suns
version of xargs. I work in a mixed Solaris/Linux environment
and have been using Unix for many years which makes it hard to keep up
with the GNU extensions to what were familiar commands. But such is
progress. Bring it on!

J for a Py Guy

2009-09-02T01:19:00.004+01:00

I just started reading a little bit more about the J programming language, after a lot of J activity on RC. I read the forward to "J for C Programmers" where I picked out an acknowledgement to Ken Iverson whose name rang a bell.

It seems that the J-language is from the APL family tree, but uses the ASCII character set and with later ideas from function-level programming added. It is good for MIMD machines it says, so I wonder how well it works for todays cheaper four to sixteen cores in a box machines, or whether it needs massively MIMD machines to really shine?

It has tickled my fancy. I'll read some more of J for C.

Saw this, and thought of you...

2009-08-27T08:01:00.005+01:00

Here in the UK, it is back-to-school time with many parents buying new school computers for their children. Windows 7 is not released yet, but here is an alternative opinion on the operating system from a competitor:

http://windows7sins.org/#7

The Story of the Regexp and the Primes

2009-08-26T07:47:00.000+01:00

If you are all sitting comfortably, then I shall begin...

Once upon a time, in a land an internet away, A python programmer was
browsing the web when he StumbledUpon a post with a curious regular
expression purporting to test numbers for primality. On further
investigation, the regexp was attributed in more than one place to
"Abigail" and took the form:

style="font-weight: bold; font-family: monospace;">^1?$|^(11+?)\1+$

Wanting to find out how it worked, the programmer was stumped, as all
the posts mentioning the regexp asked you to visit a page with a
supposedly excellent explanation of its inner workings, but the domain
name was no longer present.

The programmer decided to work it out for himself and came up with the
following...

A prime number is a
positive integer that is only divisible by itself and one. Zero and one
are not prime.

Let's say a number N, greater than one, is not prime. Then
there exists two numbers, S and T, that when multiplied together equal
N. We have:

S x T = N

Lets us further assume
that T is greater than, or equal to S. (They can be swapped if
necessary to make this so). Then a small bit of manipulation gives:

(S x (T-1)) + S = N

(You take one less S in
the multiplication, then you need to add the S back).

Lets swap the terms and write this as equation (a):

style="font-weight: bold;">S + (S x (T-1)) = N style="font-weight: bold;"> when S,T,N are integers greater
than one, and N is style="font-style: italic; font-weight: bold;">not style="font-weight: bold;"> prime

Abigail's regexp relies on representing integers as strings of that
number of ones:

Zero would be
    ''

One would be     '1'

two becomes     '11'

and so on...

The regexp gives a match if the string of ones does not represent a
prime.

To match zero or one occurrences of a 1 we use the first half of the
regexp:

style="font-family: monospace;">^1?$

Where ^anchors the following to a match from the beginning of the
string, ? will match zero or one occurrences of the preceding 1, and
the trailing $ matches the end of the string, i.e. it matches
a string that contains only zero or one 1 characters.

Going back to equation (a), greater than one, is the same as two or
more.

So S will be represented by two or more 1's - we don't know how many,
but we do know that it is two or more. S, in a regexp, becomes
something like:

style="font-family: monospace;">11+

1 followed by one or more 1's. S is less than or equal to T,
which translates to using a less greedy form of the zero-or-one
matcher, which is +? giving a representation of S as:

style="font-family: monospace;">11+?

We want to find S followed by T-1 identical copies of S, so lets form a
group out of what matches S:

style="font-family: monospace;">(11+?)

Then by referring to this group, it is the first group in the regexp so
can be referred to as \1, we can search for extra copies of the group
by:

style="font-family: monospace;">\1+

S and its copies must cover the whole of the number so must start at ^
and finish at $:

style="font-weight: bold; font-family: monospace;">^(11+?)\1+$

Together with the earlier regexp fraction that covers the numbers zero
and one, they can be or'd together to cover all non primes as
the following:

style="font-family: monospace;">^1?$|^(11+?)\1+$

...Not knowing when to quit, the programmer opted to display some of
his Python code that defines a primality tester using the regexp:

>>> import re

>>> def isprime(n):

    return not re.match(r' style="font-weight: bold;">^1?$|^(11+?)\1+$', '1' * n)

>>> [i for i in range(40) if isprime(i)]

[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]

And another function that shows S and T from the regexp match:

>>> def is style="font-weight: bold;">notprime(n):

    notprime = re.match(r'^(1?)$|^(11+?)(\2+)$', '1' * n)

    if notprime:

        if n <=1:

            return (len( notprime.groups()[0]), )

        else:

            S, SxT1 = notprime.groups()[1:]

            s = len(S)

            t = 1 + len(SxT1)//s

            return (s,t)

    else:

        return ()

>>> [(i, isnotprime(i)) for i in range(15)]

[(0, (0,)), (1, (1,)), (2, ()), (3, ()), (4, (2, 2)), (5, ()), (6, (2, 3)), (7, ()), (8, (2, 4)), (9, (3, 3)), (10, (2, 5)), (11, ()), (12, (2, 6)), (13, ()), (14, (2, 7))]

... And so to bed.

THE END.

Words From Hex Letters

2009-08-08T09:54:00.000+01:00

Enter the six hexadecimal letters into an href="http://rosettacode.org/wiki/Anagrams#Python">anagram
solver with a large scrabble-like href="http://www.puzzlers.org/pub/wordlists/unixdict.txt">dictionary
and you can generate a list of sixteen 'English' words that are also
hex numbers. Good for sprinkling around your assembler listings to
spice up the usual DEADBEEF. CAFÉABBÉ anyone? (confess over a coffee).

style="text-align: left; width: 189px; height: 648px; margin-left: 40px;"
border="1" cellpadding="2" cellspacing="2">

WORD
MEANING

ABBE
priest/clergyman

ABED
in bed/ asleep

BABE
infant

BADE
bid, appeal

BEAD
droplet

BEEF
meat from cow. complaint.

CAFE
small restaurant

CEDE
surrender, abdicate

DADA
father, dad

DADE
<None found>

DEAD
"pining for the href="http://en.wikipedia.org/wiki/Fjords" title="Fjords"
class="mw-redirect">fjords", stunned,
snuffed-it.

DEAF
"hard of hearing"

DEED
contract. accomplishment

FACE
boldness. front

FADE
disappear gradually.

FEED
food. give food

P.S. Anyone have a meaning for DADE?

Facebook User-Clustering Puzzle

2009-08-02T18:32:00.000+01:00

I came across the following href="http://www.facebook.com/careers/puzzles.php?puzzle_id=8">puzzle
at facebook. The task is to generate clusters of users based on who
they send emails to, as extracted from their large log file.

They give a strict format for the log file entries of:

style="font-family: monospace;"><date> '\t'
<send email address> '\t' <to email address>

Some lines from the log might look like:

Thu Dec 11 18:17:01 PST 2008    Finley.Cox@xx.com    James.Harrison@xx.com
Thu Dec 11 18:21:01 PST 2008    Alfie.Khan@xx.com    Tyler.Blanchet@xx.com
Thu Dec 11 18:23:01 PST 2008    Harry.White@xx.com   Finley.Ali@xx.com

They then state that a cluster exists if everyone in the group has sent
at least one email to everyone else in the group, so a cluster of users
A B and C, would have had to have sent the following emails at some
time:

from: A to B

from: A to C

from: B to A

from: B to C

from: C to A

from: C to B

You should extract the clusters from the log and report only maximal
clusters, i.e. in the example above, don't report A,C as a cluster if
you are also reporting A,B,C.

Their is also a format for how clusters are to be reported:

Only clusters of three or more people are to be reported

Reports are to be printed to stdout.

A cluster to a line.

Emails in the cluster to appear in sorted order; style="font-style: italic;">separated by the
two characters of a single comma followed by a single space.

The lines of clusters to also be sorted on the characters
appearing in the each cluster line.

An output for some log file might look like:

Aaron.Lewis@xx.com, Callum.Bennett@xx.com, Jayden.Mills@xx.com, Lucas.Lewis@xx.com, Matthew.Kennedy@xx.com
Alfie.Clarke@xx.com, Amy.Russell@xx.com, Connor.King@xx.com, Ethan.Thomas@xx.com, Joshua.Mcdonald@xx.com

The log files are said to be large.

The above is the original problem, re-phrased.

The Reverse Problem

Early on in thinking about how to solve the problem I thought I would
need much more than the ten lines or so of sample log lines
that they gave. I then thought that style="font-style: italic;">generating a log
file would be a good task too and set out to solve that
instead!

I need to:

Generate pseudo-random email addresses

Generate clusters of users

Generate a log from the clusters

And in that order too. When I have emails, I can create
clusters of emails of different sizes, then expand the clusters into a
log file structure, then print the log in the correct format.

The code follows, this description.

Line 9 is a debug remnant.

Although have only run the code in Python 3, I think it
should work in Python 2.6 as well.

Lines 10-13 are because intern is moved to the sys module in Python 3.

Lines 15-40: I needed thousands of names so trawled WP for href="http://en.wikipedia.org/wiki/List_of_most_popular_given_names">first
names and href="http://en.wikipedia.org/wiki/List_of_most_common_surnames">family
names to create their product.

Lines 27 and 41: I intern'd the strings because I could -
It's a possible premature style="text-decoration: line-through;">ejacu
optimization, but hey.

Lines 43-45: I will be building the log without times but
in order, then adding times to the ordered log entries. to do that i
need the time to start, and some possible time increments between log
entries that I will randomly choose between.

Line 48; The size of cluster I generate is a random choice
between the sizes mentioned here.. For example, a cluster size of 2, ( style="font-family: monospace;">[2]*20), will
predominate.

Line 51 controls how many style="font-style: italic;">extra times an
email might be sent between two people.

Line 53: namegen
just pairs every first name with every family name. If you want more,
you could use two first names for example, as well as extending the
name lists. (or hyphenated family names). Anyone for Imogen Ethan
Islam-Singh :-)

Line 60: name2email
changes the name from namegen into an email address string

Lines 64-69: clustergen
generates clustercount clusters of names. It randomly chooses how many
names to put in each cluster using clustersizefreq;, then generates
that many names per cluster. nextname feeds new names are required, and
the sorts ensure the clusters are generated in the correct order for
printing to satisfy the original puzzles output criteria.

Lines 71-74: clusterprint will correctly format and print
the clustergen data structure.

Lines 76-95: logggen.
The largest function

Line 81 creates a number of clusters

Line 83 takes the cluster of users and turns each cluster
in turn into the expanded set of correspondances between every member
of the cluster, in order.

Line 86 randomly duplicates some entries of the log using
repeatmsgfraction.

Line 89 adds some authenticity by randomising the order
of the emails

Lines 91-94 add the time information to the log.

Line 95 returns the cluster info as well as the log.

Lines 97-101 function logprint can correctly format and
print the log data. Note the handling of the PST timezone in line 100.

 1 
 color="#804040"> 2 '''
 color="#804040"> 3  Random Data generator for:
 color="#804040"> 4   href="http://www.facebook.com/careers/puzzles.php?puzzle_id=8">http://www.facebook.com/careers/puzzles.php?puzzle_id=8
 color="#804040"> 5 '''
 color="#804040"> 6 
 color="#804040"> 7 
 color="#804040"> 8 import itertools, random, datetime
 color="#804040"> 9 from pprint  color="#a020f0">import pprint  color="#a020f0">as pp
 color="#804040"> 10 try:  color="#0000ff"># Python 3/2 compatibility
 color="#804040"> 11  from sys  color="#a020f0">import intern
 color="#804040"> 12 except:
 color="#804040"> 13  pass
 color="#804040"> 14 
 color="#804040"> 15 firstname = ''' color="#ff00ff"> 
 color="#804040"> 16  Lola Imogen Jasmine Leah Daniel Isla Millie Lilly Benjamin Evie
 color="#804040"> 17  Finlay Ella Jacob Grace Jayden Freya Summer James Adam Henry
 color="#804040"> 18  Poppy Joshua Charlotte Thomas Noah Maisie Tyler Liam George Mason
 color="#804040"> 19  Matthew Joseph Nathan Stan Scarlett Jessica Samuel Leo William Olivia
 color="#804040"> 20  Ruby Isabella Emma Isabelle Charlie Phoebe Alfie Owen Erin Oscar
 color="#804040"> 21  Cameron Katie Harvey Aaron Ethan Lily Finley Riley Jamie Megan
 color="#804040"> 22  Caitlin Amy Connor Alexander Ben Eva Bethany Brooke Layla Harry
 color="#804040"> 23  Chloe Harrison Lauren Ava Oliver Mia Abigail Ellie Holly Jake
 color="#804040"> 24  Hannah Luke Molly Jack Amelia Dylan Alex Lucas Lucy Ryan
 color="#804040"> 25  Daisy Sophie Sophia Archie Callum Lewis Isabel Logan Max Emily
 color="#804040"> 26  '''.split()
 color="#804040"> 27 firstname = [intern(n)  color="#804040">for n  color="#804040">in firstname]
 color="#804040"> 28 
 color="#804040"> 29 familyname = list(set( ''' color="#ff00ff"> 
 color="#804040"> 30  Wilson Russell Smith Phillips Dixon White Gill Brown Statham Oneill
 color="#804040"> 31  Harris Doyle James Ahmad Roberts Mann Moore Thomas Saunders Barker
 color="#804040"> 32  Bright Allen George Hill Jackson Rose Malley Mills Sheppard Mason
 color="#804040"> 33  Wright Holt Pearce Bull Richards Ward Lopez Mathey King Stone
 color="#804040"> 34  Johnson Baker Dupont Green Davis Blanchet Taylor Dean Donnelly Evans
 color="#804040"> 35  Woods Patel Meacher Simon Clarke Young Ellis Palmer Lynch Alexander
 color="#804040"> 36  Sutton Ahmed Adams Kennedy Davies Walker Austin Fernandez Mustafa Rodriguez
 color="#804040"> 37  Bennett Murray Powell Harrison Campbell Cole Khan Mcdonald Edwards Wood
 color="#804040"> 38  Thompson Singh Williams Power Price Watson Hall Jones Cox Chapman
 color="#804040"> 39  Arnold Reid Garcia Morgan Mitchell Ali Lewis Stubbs Scott Islam
 color="#804040"> 40  '''.split() ))
 color="#804040"> 41 familyname = [intern(n)  color="#804040">for n  color="#804040">in familyname]
 color="#804040"> 42 
 color="#804040"> 43 starttime = datetime.datetime(2008, 12, 11, 17, 53, 1, 0)
 color="#804040"> 44 # Choice of delta times between log entries
 color="#804040"> 45 timedeltafreq = [ datetime.timedelta(0, secs)  color="#804040">for secs  color="#804040">in [60]*3 + [120]*10 + [240]*10 ]
 color="#804040"> 46 
 color="#804040"> 47 # Choice of cluster sizes
 color="#804040"> 48 clustersizefreq = [2]*20 + [3]*8 + [4]*3 + [5]*2 + [6]*1 +[7]*1 + [8]*1
 color="#804040"> 49 
 color="#804040"> 50 # Many messages between the same two people?
 color="#804040"> 51 repeatmsgfraction = 2.0  color="#0000ff"># Extra 200%
 color="#804040"> 52 
 color="#804040"> 53 def  color="#008080">namegen():
 color="#804040"> 54  'Generates (first, family) name tuples'
 color="#804040"> 55  names = list(itertools.product(firstname, familyname))
 color="#804040"> 56  random.shuffle(names)
 color="#804040"> 57  for name  color="#804040">in names:
 color="#804040"> 58  yield name
 color="#804040"> 59 
 color="#804040"> 60 def  color="#008080">name2email(name):
 color="#804040"> 61  'format (first, family) name as email address'
 color="#804040"> 62  return ' color="#ff00ff">%s.%s@xx.com' % name
 color="#804040"> 63 
 color="#804040"> 64 def  color="#008080">clustergen(clustercount, clustersizefreq=clustersizefreq,
 color="#804040"> 65  firstname=firstname, familyname=familyname):
 color="#804040"> 66  'Generate clustercount clusters of unique users'
 color="#804040"> 67  nextname = namegen().__next__
 color="#804040"> 68  return sorted( sorted(nextname()  color="#804040">for i  color="#804040">in range(random.choice(clustersizefreq)))
 color="#804040"> 69  for j  color="#804040">in range(clustercount) )
 color="#804040"> 70 
 color="#804040"> 71 def  color="#008080">clusterprint(clusters):
 color="#804040"> 72  'Print already sorted clusters in output format'
 color="#804040"> 73  for cluster  color="#804040">in clusters:
 color="#804040"> 74  print ( ' color="#ff00ff">, '.join(name2email(n)  color="#804040">for n  color="#804040">in cluster) )
 color="#804040"> 75 
 color="#804040"> 76 def  color="#008080">loggen(clustercount, clustersizefreq=clustersizefreq,
 color="#804040"> 77  firstname=firstname, familyname=familyname,
 color="#804040"> 78  starttime=starttime, timedeltafreq=timedeltafreq,
 color="#804040"> 79  repeatmsgfraction=repeatmsgfraction):
 color="#804040"> 80  'Generate a log of clustercount clusters of unique users'
 color="#804040"> 81  clusters = clustergen(clustercount)
 color="#804040"> 82  # clustered emails
 color="#804040"> 83  log = sum([ list(itertools.permutations(cluster, 2))  color="#804040">for cluster  color="#804040">in clusters ],
 color="#804040"> 84  [])
 color="#804040"> 85  # repeats
 color="#804040"> 86  for i  color="#804040">in range(int(repeatmsgfraction*len(log))):
 color="#804040"> 87  log.append(random.choice(log))
 color="#804040"> 88  # dispersed emails
 color="#804040"> 89  random.shuffle(log)
 color="#804040"> 90  # log times for emails
 color="#804040"> 91  logtime = starttime
 color="#804040"> 92  for n, (send, to)  color="#804040">in enumerate(log):
 color="#804040"> 93  log[n] = (logtime, send, to)
 color="#804040"> 94  logtime += random.choice(timedeltafreq)
 color="#804040"> 95  return clusters, log
 color="#804040"> 96 
 color="#804040"> 97 def  color="#008080">logprint(log):
 color="#804040"> 98  'Print log in input format'
 color="#804040"> 99  for logtime, send, to  color="#804040">in log:
 color="#804040">100  print (' color="#ff00ff">%s\t color="#ff00ff">%s\t color="#ff00ff">%s' % (logtime.strftime(" color="#ff00ff">%a %b %d %H:%M:%S PST %Y")
 color="#804040">101  , name2email(send), name2email(to)))
 color="#804040">102

The Program in use

Switching to the command
line shell, I generate a log and its clustering and print the results.
Because of the way the log is generated from the clustering, the
clustering of the log should be as stated.

style="font-family: monospace;">>>> style="font-weight: bold;">clusters, log = loggen(3);
clusterprint(clusters); print(); logprint(log); print(); len(log) style="font-family: monospace;">
Alfie.Khan@xx.com,
Lily.Rose@xx.com, Tyler.Blanchet@xx.com style="font-family: monospace;">
Finley.Ali@xx.com,
Harry.White@xx.com

Finley.Cox@xx.com,
James.Harrison@xx.com

Thu Dec 11 17:53:01
PST 2008
Finley.Cox@xx.com    James.Harrison@xx.com style="font-family: monospace;">
Thu Dec 11 17:55:01
PST 2008
Lily.Rose@xx.com    Alfie.Khan@xx.com style="font-family: monospace;">
Thu Dec 11 17:59:01
PST 2008
Finley.Cox@xx.com    James.Harrison@xx.com style="font-family: monospace;">
Thu Dec 11 18:03:01
PST 2008
Finley.Ali@xx.com    Harry.White@xx.com style="font-family: monospace;">
Thu Dec 11 18:07:01
PST 2008
Tyler.Blanchet@xx.com    Alfie.Khan@xx.com style="font-family: monospace;">
Thu Dec 11 18:09:01
PST 2008
Tyler.Blanchet@xx.com    Alfie.Khan@xx.com style="font-family: monospace;">
Thu Dec 11 18:13:01
PST 2008
Alfie.Khan@xx.com    Lily.Rose@xx.com style="font-family: monospace;">
Thu Dec 11 18:17:01
PST 2008
Finley.Cox@xx.com    James.Harrison@xx.com style="font-family: monospace;">
Thu Dec 11 18:21:01
PST 2008
Alfie.Khan@xx.com    Tyler.Blanchet@xx.com style="font-family: monospace;">
Thu Dec 11 18:23:01
PST 2008
Harry.White@xx.com    Finley.Ali@xx.com style="font-family: monospace;">
Thu Dec 11 18:27:01
PST 2008
Harry.White@xx.com    Finley.Ali@xx.com style="font-family: monospace;">
Thu Dec 11 18:29:01
PST 2008
James.Harrison@xx.com    Finley.Cox@xx.com style="font-family: monospace;">
Thu Dec 11 18:33:01
PST 2008
Alfie.Khan@xx.com    Tyler.Blanchet@xx.com style="font-family: monospace;">
Thu Dec 11 18:37:01
PST 2008
Alfie.Khan@xx.com    Lily.Rose@xx.com style="font-family: monospace;">
Thu Dec 11 18:39:01
PST 2008
James.Harrison@xx.com    Finley.Cox@xx.com style="font-family: monospace;">
Thu Dec 11 18:43:01
PST 2008
Tyler.Blanchet@xx.com    Alfie.Khan@xx.com style="font-family: monospace;">
Thu Dec 11 18:45:01
PST 2008
Alfie.Khan@xx.com    Tyler.Blanchet@xx.com style="font-family: monospace;">
Thu Dec 11 18:47:01
PST 2008
Harry.White@xx.com    Finley.Ali@xx.com style="font-family: monospace;">
Thu Dec 11 18:49:01
PST 2008
Harry.White@xx.com    Finley.Ali@xx.com style="font-family: monospace;">
Thu Dec 11 18:53:01
PST 2008
James.Harrison@xx.com    Finley.Cox@xx.com style="font-family: monospace;">
Thu Dec 11 18:57:01
PST 2008
Finley.Cox@xx.com    James.Harrison@xx.com style="font-family: monospace;">
Thu Dec 11 18:59:01
PST 2008
Harry.White@xx.com    Finley.Ali@xx.com style="font-family: monospace;">
Thu Dec 11 19:03:01
PST 2008
Lily.Rose@xx.com    Alfie.Khan@xx.com style="font-family: monospace;">
Thu Dec 11 19:04:01
PST 2008
Tyler.Blanchet@xx.com    Lily.Rose@xx.com style="font-family: monospace;">
Thu Dec 11 19:06:01
PST 2008
James.Harrison@xx.com    Finley.Cox@xx.com style="font-family: monospace;">
Thu Dec 11 19:10:01
PST 2008
Tyler.Blanchet@xx.com    Alfie.Khan@xx.com style="font-family: monospace;">
Thu Dec 11 19:12:01
PST 2008
Lily.Rose@xx.com    Tyler.Blanchet@xx.com style="font-family: monospace;">
Thu Dec 11 19:16:01
PST 2008
Finley.Cox@xx.com    James.Harrison@xx.com style="font-family: monospace;">
Thu Dec 11 19:20:01
PST 2008
Tyler.Blanchet@xx.com    Lily.Rose@xx.com style="font-family: monospace;">
Thu Dec 11 19:22:01
PST 2008
Alfie.Khan@xx.com    Tyler.Blanchet@xx.com style="font-family: monospace;">

30

To Do

The cluster info for the puzzle should, of course, have the
clusters of size two filtered out

I haave a solution to the original puzzle, but: " style="font-style: italic;" title="Insert schoolboy excuse here">Blue
Martians forcibly wiped it from my brain" :-)

- Paddy.

The case of the disappearing over-bar

2009-07-28T09:16:00.000+01:00

Rosetta-code has a task
which asks you to reverse a Unicode string
correctly. I had glanced at the Python
solution and thought nothing of
it until someone did something similar in the R language and stated
that it may be incorrect for the given pattern which includes an
over-bar (my description) over the f in "as⃝df̅" (See the
orginal article for the true Unicode string).

I cut-n-pasted the Python solution and the test string into Python 3.1
idle and decided to test it:

Python 3.1 (r31:73574, Jun 26 2009, 20:21:35) [MSC v.1500 32 bit (Intel)] on win32
Type "copyright", "credits" or "license()" for more information.
>>> # I cut and paste the following string with an over-bar on the f
>>> x = input()
as⃝df̅
>>> # Just showing x gives the over-bar on the f
>>> x
'asâƒdfÌ…'
>>> # Print it and it is fine.
>>> print(x)
asâƒdfÌ…
>>> # Reverse x though, and the over-bar movesover the quote!
>>> x[::-1]
'Ì…fdâƒsa'
>>> # print the reversed x and it disappears altogether!
>>> print(x[::-1])
Ì…fdâƒsa
>>>

|NowCut and paste the unicode characters from firefox for the input
statement, and make sure that:

>>> ['%x' % ord(char) for char in x]
['61', '73', '20dd', '64', '66', '305']

As the unicode might otherwise be lost in my editing to try and present
it to you above.

A bit of scary reading introduced me to
Unicode character classes. It seems that Unicode characters are
sometimes composed, and that if I know what character gets composed
with another, (always the last non-composed character to the left),
then you should be able to form composed groups of characters and then
reverse them.

The WP article gave me the vocabulary, so a quick search in Pythons
library gave me the unicodedata
module which has the name
function. I am not sure if this will work in every case, but from the
WP article and this experimentation:

>>> [(c,'%s' % unicodedata.name(c,0xfffff)) for c in x]
[('a', 'LATIN SMALL LETTER A'), ('s', 'LATIN SMALL LETTER S'), ('⃝', 'COMBINING ENCLOSING CIRCLE'), ('d', 'LATIN SMALL LETTER D'), ('f', 'LATIN SMALL LETTER F'), ('̅', 'COMBINING OVERLINE')]

I think I can group by the presence of the word COMBINING in the name
of a character and so produced the following reversal function:

'''
  Reverse a Unicode string with proper handling of combining characters
'''

import unicodedata

def ureverse(ustring):
    '''
    Reverse a string including unicode combining characters

    Example:
        >>> ucode = ''.join( chr(int(n, 16))
                             for n in ['61', '73', '20dd', '64', '66', '305'] )
        >>> ucoderev = ureverse(ucode)
        >>> ['%x' % ord(char) for char in ucoderev]
        ['66', '305', '64', '73', '20dd', '61']
        >>> 
    '''
    groupedchars = []
    uchar = list(ustring)
    while uchar:
        if 'COMBINING' in unicodedata.name(uchar[0], ''):
            groupedchars[-1] += uchar.pop(0)
        else:
            groupedchars.append(uchar.pop(0))
    # Grouped reversal
    groupedchars = groupedchars[::-1]

    return ''.join(groupedchars)

if __name__ == '__main__':
    ucode = ''.join( chr(int(n, 16))
                     for n in ['61', '73', '20dd', '64', '66', '305'] )
    ucoderev = ureverse(ucode)
    print (ucode)
    print (ucoderev)

It works for the given example text (Try running it to see the output -
I've given up trying to work out what characters you might see).

Gosh, an attractive Unicode issue! Whatever next.

Type based API's and Duck Typing Shocker!!

2009-06-24T22:18:00.004+01:00

The Voidspace tech blog has href="http://www.voidspace.org.uk/python/articles/duck_typing.shtml">an
article in which the author, Michael Ford, blames duck-typing
for a problem in trying to find the type of an argument.

After a very good explanation of how duck typing in Python
comes about Michael starts his section on the problem by
stating the principal of duck typing as:

"The
principle of duck typing says that you shouldn't care what type of
object you have - just whether or not you can do the required action
with your object."

All well and good but then later on Michael states:

"If we
want our code to treat different types of object differently then
the approach in example two fails. This isn't contrived - this is
exactly the situation we found ourselves in with ConfigObj."

So, type is important in his case. Fine. style="font-weight: bold;">But why lay the blame at duck
typings door ?

His example interface relies on knowing what constitutes a "listy"
object, or a "stringy" object, etc, in Python 2.x, which he notes is a
pain, but I think the cause of the problem w.r.t. duck typing is
choosing an API that relies on an ill-defined idea of what constitutes
close enough to a dict/list/string or whatever, then not having the
language provide an easy solution.

Abstract
base classes are an attempt to help out in such cases that is
new in Python 2.6 and 3.0, but if you want to link functionality with
type then you can't expect to get duck typing too.

Where Duck Typing Does Not fit.

If a function has an expensive or non-reversible operation that depends
on an attribute of an argument, then it may be wise to check
all arguments for compatibility up-front; but even then, whats wrong
with reading the code? Defensive programming can creep up on
you...

Killer Script

2009-06-22T19:55:00.000+01:00

Following on from my href="http://paddy3118.blogspot.com/2009/03/batch-process-runner-in-bash-shell.html"
target="_blank">Batch process runner script, I got
thinking that a nice capability would be to have a maximum runtime
limit for the jobs.

My first thought was to create an href="http://unixhelp.ed.ac.uk/CGI/man-cgi?at" target="_blank">at
command set to execute the reuired time after the start of each job,
one for each job, that would kill the job.

Nah.

Discussions with a friend lead to consideration of using href="http://docs.sun.com/app/docs/doc/816-5165/ulimit-1?l=en&a=view&q=ulimit"
target="_blank">ulimit in the script to limit the
run time, but unfortunately on the Unix box I was testing on, ulimit
could only limit the maximum CPU time, not run time. If the job was
stuck in a non-busy wait then it would not be auto-killed.

I decided to try creating a script to run a time eating processing task
that created a background kill task so the script would be self
killing.

The script

Line 1:    The KILLAFTER environment
variable will kill the script after 3 seconds

Line 5:    Store the process ID for the
'bash -c' script for killing later

Line 8:    This is the slleping background
sub-process that wakes after KILLAFTER seconds and kills the script.

Line 9:    But if the normal script commands
finish, store the sub-process PID so it can be cleanly killed itself.

Line 12: Some dummy command that counts up every second.

Line 13:   And save its exit status so it can be
restored after removing the background kill process.

Line 16:   Terminate the background kill process so
we are not left waiting.

Line 18:   Exit with the saved exit status.

Line 19:   Prints the exit status of the whole 'bash
-c' script

Lines 20-24 show the output of the first run. Notice how although the
for loop in line 12 is set to count to 5, the count gets killed after
printing 3, i.e. the three seconds of $KILLAFTER.

Lines 27 onwards show a second run where KILLAFTER is set to more time
than is used by the for loop. Notice how the script prints the full
count up to 5, and has a return status of zero.

 1 bash-paddy:  style="font-weight: bold;">env  style="font-weight: bold;" color="#008080">KILLAFTER style="font-weight: bold;">= style="font-weight: bold;" color="#ff00ff">3 style="font-weight: bold;"> bash  style="font-weight: bold;" color="#6a5acd">-c style="font-weight: bold;"> '
 color="#804040"> 2 # New bash script shell environment
 color="#804040"> 3 
 color="#804040"> 4 # Its PID
 color="#804040"> 5 export  color="#008080">thisshellpid color="#804040">= color="#a020f0">$$
 color="#804040"> 6 
 color="#804040"> 7 # Sleeping killer will blast this script after given time
 color="#804040"> 8 (sleep  color="#a020f0">$KILLAFTER  color="#804040">;  color="#804040">kill  color="#ff00ff">-9 $thisshellpid color="#804040">) color="#804040">&
 color="#804040"> 9 killerpid= color="#a020f0">$!
 color="#804040">10 
 color="#804040">11 # Any time consuming command
 color="#804040">12 for y  color="#804040">in  color="#ff00ff">1 2  color="#ff00ff">3 4  color="#ff00ff">5;  color="#804040">do sleep  color="#ff00ff">1;  color="#804040">echo color="#ff00ff"> $y color="#804040">;  color="#804040">done
 color="#804040">13 trueexitstatus= color="#a020f0">$?
 color="#804040">14 
 color="#804040">15 # If time consuming command worked OK then...
 color="#804040">16 kill  color="#ff00ff">-9 $killerpid
 color="#804040">17 
 color="#804040">18 exit  color="#a020f0">$trueexitstatus
 color="#804040">19  style="font-weight: bold;">'  style="font-weight: bold;" color="#804040">; style="font-weight: bold;">  style="font-weight: bold;" color="#804040">echo style="font-weight: bold;" color="#ff00ff"> Returned status:  style="font-weight: bold;" color="#a020f0">$?
 color="#804040">20 1
 color="#804040">21 2
 color="#804040">22 3
 color="#804040">23 Killed
 color="#804040">24 Returned status:  color="#ff00ff">137
 color="#804040">25 bash-paddy:
 color="#804040">26 bash-paddy:
 color="#804040">27 bash-paddy:  style="font-weight: bold;">env  style="font-weight: bold;" color="#008080">KILLAFTER style="font-weight: bold;">= style="font-weight: bold;" color="#ff00ff">10 style="font-weight: bold;"> bash  style="font-weight: bold;" color="#6a5acd">-c style="font-weight: bold;"> '
 color="#804040">28 # New bash script shell environment
 color="#804040">29 
 color="#804040">30 # Its PID
 color="#804040">31 export  color="#008080">thisshellpid color="#804040">= color="#a020f0">$$
 color="#804040">32 
 color="#804040">33 # Sleeping killer will blast this script after given time
 color="#804040">34 (sleep  color="#a020f0">$KILLAFTER  color="#804040">;  color="#804040">kill  color="#ff00ff">-9 $thisshellpid color="#804040">) color="#804040">&
 color="#804040">35 killerpid= color="#a020f0">$!
 color="#804040">36 
 color="#804040">37 # Any time consuming command
 color="#804040">38 for y  color="#804040">in  color="#ff00ff">1 2  color="#ff00ff">3 4  color="#ff00ff">5;  color="#804040">do sleep  color="#ff00ff">1;  color="#804040">echo color="#ff00ff"> $y color="#804040">;  color="#804040">done
 color="#804040">39 trueexitstatus= color="#a020f0">$?
 color="#804040">40 
 color="#804040">41 # If time consuming command worked OK then...
 color="#804040">42 kill  color="#ff00ff">-9 $killerpid
 color="#804040">43 
 color="#804040">44 exit  color="#a020f0">$trueexitstatus
 color="#804040">45  style="font-weight: bold;">'  style="font-weight: bold;" color="#804040">; style="font-weight: bold;">  style="font-weight: bold;" color="#804040">echo style="font-weight: bold;" color="#ff00ff"> Returned status:  style="font-weight: bold;" color="#a020f0">$?
 color="#804040">46 1
 color="#804040">47 2
 color="#804040">48 3
 color="#804040">49 4
 color="#804040">50 5
 color="#804040">51 Returned status:  color="#ff00ff">0
 color="#804040">52 bash-paddy:
 color="#804040">53

XKCD Simpler Knapsack Solution

2009-06-15T01:01:00.000+01:00

After reading the feedback on my href="http://paddy3118.blogspot.com/2009/06/xkcd-knapsack-solution.html"
target="_blank">earlier solution to the XKCD knapsack problem,
I decided to re-write it without itertools.product and with integer
maths.

I came up with a solution that used explicit while loops that worked,
and used much less traversals through nested loops, but I then boiled
that down to a recursive solution which retains the ability to easily
modify the number of dishes , whilst restricting the number of times
through the loops to the bare minimum.

In the previous solution I worked out the maximum number of each dish
that cost less than or equal to the amount and looped from zero to that
many times. In the solution below the while loops cut off looping on
any particular dish when adding another would put the cost over the
maximum.

And now the code:

items = ( ('MIXED FRUIT',   2.15),
          (' color="#ff00ff">FRENCH FRIES',  2.75),
          (' color="#ff00ff">SIDE SALAD',    3.35),
          (' color="#ff00ff">HOT WINGS',     3.55),
          (' color="#ff00ff">MOZZ STICKS',   4.20),
          (' color="#ff00ff">SAMPLER PLATE', 5.80),
          (' color="#ff00ff">BARBEQUE',      6.55) )

exactcost = 15.05

dishes, prices = zip(*items)
 color="#0000ff"># All calcs in integer cents.
prices = [int(price*100)  color="#804040">for price  color="#804040">in prices]
ecost = int(exactcost*100)

 color="#0000ff"># counts of each dish
dishcounts = [0]*len(prices)
 color="#0000ff">possibleorders = []
cost = 0
maxnesting = len(dishcounts)

 color="#804040">def  color="#008080">enumeratedish(nest, cost, possibleorders):
     color="#804040">global dishcounts, prices, ecost, maxnesting
     color="#804040">while cost <= ecost:
         color="#804040">if nest+1 < maxnesting:
            enumeratedish(nest+1, cost, possibleorders)
         color="#804040">elif cost == ecost:
            possibleorders.append(dishcounts[:])
             color="#804040">break
        dishcounts[nest] += 1
        cost += prices[nest]
         color="#0000ff">#print (dishcounts)
    cost -= prices[nest] * dishcounts[nest]
    dishcounts[nest] = 0

enumeratedish(0, 0, possibleorders)
 color="#804040">print(' color="#6a5acd">\nALL POSSIBLE CHOICES OF MENU ITEMS THAT COST %.2f, ONE PER LINE' % exactcost)
 color="#804040">for order  color="#804040">in possibleorders:
    order = zip(dishes, order)
     color="#804040">print( ' color="#ff00ff">  ',
           ' color="#ff00ff">, '.join(' color="#ff00ff">%s: %i' % dishcount  color="#804040">for dishcount  color="#804040">in order  color="#804040">if dishcount[1]))

feel free to compare the two solutions.

XKCD Knapsack Solution

2009-06-13T21:34:00.001+01:00

A brute force solution to the following comic from href="http://xkcd.com/287/">XKCD:

title="General solutions get you a 50% tip."
src="http://imgs.xkcd.com/comics/np_complete.png">

Assuming the one sandwich represents the end of the menu, I used an
exhaustive search algorithm. Things to note are the use of zip(*table)
to transpose the items and itertools.product instead of nested for
loops.

'''
from:  href="http://xkcd.com/287/">http://xkcd.com/287/
'''

 color="#a020f0">from itertools  color="#a020f0">import product

items = ( (' color="#ff00ff">MIXED FRUIT',   2.15),
          (' color="#ff00ff">FRENCH FRIES',  2.75),
          (' color="#ff00ff">SIDE SALAD',    3.35),
          (' color="#ff00ff">HOT WINGS',     3.55),
          (' color="#ff00ff">MOZZ STICKS',   4.20),
          (' color="#ff00ff">SAMPLER PLATE', 5.80),
          (' color="#ff00ff">BARBEQUE',      6.55) )

exactcost = 15.05

dishes, prices = zip(*items)
rangeofdishes = [tuple(range(1+int( (exactcost+0.005)/price )))  color="#804040">for price  color="#804040">in prices]

possibleorders = [tuple(zip(dishes, numbers))
                   color="#804040">for numbers  color="#804040">in product(*rangeofdishes)
                   color="#804040">if int(exactcost*100)
                       == int(100* sum(num*price
                                        color="#804040">for num,price  color="#804040">in zip(numbers, prices)))
                  ]
 color="#804040">print(' color="#6a5acd">\nALL POSSIBLE CHOICES OF MENU ITEMS THAT COST %.2f, ONE PER LINE' % exactcost)
 color="#804040">for order  color="#804040">in possibleorders:
     color="#804040">print( ' color="#ff00ff">  ',
           ' color="#ff00ff">, '.join(' color="#ff00ff">%s: %i' % dishcount  color="#804040">for dishcount  color="#804040">in order  color="#804040">if dishcount[1]))

The answers I get are:

style="font-family: monospace;">ALL POSSIBLE CHOICES OF MENU
ITEMS THAT COST 15.05, ONE PER LINE style="font-family: monospace;">

MIXED FRUIT: 1, HOT WINGS: 2, SAMPLER PLATE: 1 style="font-family: monospace;">

MIXED FRUIT: 2, MOZZ STICKS: 1, BARBEQUE: 1


MIXED FRUIT: 7

(I'd go with the hot wings).

Zed and community

2009-05-31T09:48:00.001+01:00

Who is Zed?

Zed it seems knows about his blogging style: asking href="http://www.zedshaw.com/blog/2009-02-27.html">questions
as statements;and href="http://www.zedshaw.com/rants/rails_is_a_ghetto.html">ranting.but
it seems that he has written good code for the Rails community in the
past and has earned some peoples respect. Now he is a new member of the
Python community I would hope that current members realize his style
and try not to emulate it.

Zed may be here for a long stay, or gone tomorrow, but I wouldn't want
the lasting effect to be a rise in the abrasiveness of the Python
community, or an unthinking
abandonment of community practices. Even if Zed does produce
Pythons Holy Grail, we have to remember that that is software.
Community matters, and we have to work with the
Zeds, the Alex's, the Xah's, the Tim's without being seen to sink to
the lowest common denominator.

Pipe Fitting with Python Generators

2009-05-23T04:50:00.001+01:00

I was investigating the Hofstadter-Conway href="http://mathworld.wolfram.com/Hofstadter-Conway10000-DollarSequence.html"
target="_blank">$10,000 sequence which, for
mathematics, has a very href="http://www.nytimes.com/1988/08/30/science/intellectual-duel-brash-challenge-swift-response.html"
target="_blank"> interesting story. To
'follow the money', you have to generate the sequence of values and
calculate maxima and minima of the sequence values so I decided to
create a generator of the sequence values and pass it through filters
that would emit local maxima/minima. I ended up with something like:

high_generator =
maxima(sequence_generator())

highs = [high_generator.next() for i in range(1000)]

Now I just hated the nested call of generators. If it were a Unix
shell, I would have hoped to have written:

sequence_generator |
maxima |head 1000 > highs

Notice how the data flows from source at the left to sink at the right.
I like that.

I spent some time trying to use classes and redefine bitwise or for
class instances, but it was looking messy and my heart wasn't really in
a class based solution.

I left it for a week. And in that week I remembered seeing a function
called pipe being used to pipe things so decided to write my own pipe
function for combining generators in a more natural way.

Pipe

What I wanted was a function called pipe, that when called with several
generators, returned a nested call of all its arguments, with the
leftmost argument being in the center of the nest. For example, a call
of:

pipe(a,b,c,d)

Would return the equivalent of:

d(c(b(a)))

The eventual definition turned out to be quite straight-forward:

 color="#804040">def  color="#008080">pipe(*cmds):
    ''' color="#ff00ff"> pipe(a,b,c,d, ...) -> yield from ...d(c(b(a())))
 color="#ff00ff">    '''
    gen = cmds[0]
     color="#804040">for cmd  color="#804040">in cmds[1:]:
        gen = cmd(gen)
     color="#804040">for x  color="#804040">in gen:
         color="#804040">yield x

Pipe in use

Lets create some simple generators; style="font-family: monospace;">count is just
itertools.count; sqr
yields the square of items on its input iterator

from itertools  color="#a020f0">import count

 color="#804040">def  color="#008080">sqr(x):
 color="#804040">    for val  color="#804040">in x:
 color="#804040">        yield val*val

I can now generate the squares of integers:

>>> p = pipe(count(), sqr)
>>> [p.next()  color="#804040">for i  color="#804040">in range(5)]
[0, 1, 4, 9, 16]

But I need a more complex definition for head as it needs to have a
parameter of the maximum number of items to pass through.

def  color="#008080">head(n):
    ''' color="#ff00ff"> return a generator that passes through at most n items
 color="#ff00ff">    '''
     color="#804040">def  color="#008080">head(x):
        i = n
         color="#804040">for val  color="#804040">in x:
             color="#804040">if i>0:
                i -= 1
                 color="#804040">yield val
             color="#804040">else:
                 color="#804040">raise StopIteration
     color="#804040">return head

And now with head I can more naturally express nested generator
functions as I can stop an infinite generator stream easily for example:

>>> list(pipe(count(), sqr, head(4)))
[0, 1, 4, 9]
>>>  color="#804040">for i  color="#804040">in pipe(count(), sqr, head(10)):
     color="#804040">print i,

0 1 4 9 16 25 36 49 64 81

Generalized function generator and filter

I thought it might be handy to have a function that skipped an input if
the prefilter(input) was false; applied a given function, fn(input) to
the input stream and passed on fn(input) to the output if
postfilter(fn(input)) was true. Think of a stream of data. You can
filter it, apply a function to it, or filter after applying the
function - all via this pipe filter. I called it fnf for function
filter:

def  color="#008080">fnf(fn, prefilter=None, postfilter=None):
    ''' color="#ff00ff"> functionfilter:
 color="#ff00ff">            for s in x:
 color="#ff00ff">                yield fn(s) subject to prefilter(s), postfilter(fn(s)) returning True if present,
 color="#ff00ff">                or skip this s
 color="#ff00ff">    '''
     color="#804040">def  color="#008080">_fnf(x):
         color="#804040">for s  color="#804040">in x:
             color="#804040">if prefilter  color="#804040">is None  color="#804040">or prefilter(s):
                fns = fn(s)
                 color="#804040">if postfilter  color="#804040">is None  color="#804040">or postfilter(fns):
                     color="#804040">yield fns
    _fnf.__doc__ = ' color="#ff00ff">_fnf = fnf(%s, %s, %s)' %(fn, prefilter, postfilter)
     color="#804040">return _fnf

I'll show fnf in action.

>>> # turns the input stream of numbers, x into a stream of tuples x,x/2.0
>>>  color="#804040">def  color="#008080">n_nby2(x):  color="#804040">return x, x/2.0

>>> list(pipe( count(), sqr, head(5),  style="font-weight: bold; text-decoration: underline;">fnf(n_nby2) ))
[(0, 0.0), (1, 0.5), (4, 2.0), (9, 4.5), (16, 8.0)]
>>>  color="#0000ff"># Lets create a prefilter
>>>  color="#804040">def  color="#008080">even(x):  color="#804040">return x%2 == 0

>>> list(pipe( count(), sqr, head(5), fnf(n_nby2, style="font-weight: bold; text-decoration: underline;">even) ))
[(0, 0.0), (4, 2.0), (16, 8.0)]
>>>  color="#0000ff"># And now a postfilter which must work on the result of applying fn
>>>  color="#804040">def  color="#008080">frac(x):  color="#804040">return int(x[1]) != x[1]

>>> list(pipe( count(), sqr, head(5), fnf(n_nby2, style="font-style: italic; font-weight: bold;">None, style="font-weight: bold; text-decoration: underline;">frac) ))
[(1, 0.5), (9, 4.5)]
>>>  color="#0000ff"># Note the presence of  style="font-weight: bold; font-style: italic;">None in the above call as we don't use a prefilter

Thats it for now with the pipe fitting. I do have a tee function for
printing intermediate results, but it's style="text-decoration: line-through;">late
early and I'm tired....

- Paddy.

Extended unpacking in Python 3.0

2009-05-13T22:15:00.000+01:00

After reading about href="http://www.rosettacode.org/wiki/Pattern_Matching"
target="_blank">pattern matching on Rosetta Code,
I wondered, (and still am wondering), if it has anything to do with
list unpacking in Python, which lead me to do some simple explorations
of the target="_blank">extended iterable unpacking that
is new in Python 3.0.

Introducing a star

The feature builds on the pre-existing capability of being able to
unpack a list or tuple into multiple names in Python 2.x, so in Python
2.x you
could do:

>>> import platform
>>> platform.python_version()
' color="#ff00ff">2.6.1'
>>> lst = [1,2,3]
>>> a,b,c = lst
>>> print a,b,c
1 2 3
>>> tpl = (1,2,3)
>>> a,b,c = tpl
>>> print a,b,c
1 2 3
>>> iter = (x for x in range(4))
>>> iter
<generator object <genexpr> at 0x01300468>
>>> a,b,c,d = iter
>>> print a,b,c,d
0 1 2 3
>>>  color="#0000ff"># You could also do the (fancy), but non-* stuff too
>>> a,b,(c,d) = (1,2,(3,4))
>>>  color="#804040">print a,b,c,d
1 2 3 4
>>>  color="#0000ff"># But not
>>> a,*b = (1,2,3,4)
SyntaxError: invalid syntax
>>>

In Python 3.0 you get all of the above plus the ability to 'soak up'
the rest of the values using style="font-weight: bold;">* (which I pronounce
as 'star' in this
context).

It is similar to how you can use *args in a function definition to
assign extra positional arguments to a name. In an assignment however
it means 'take
whatever's left of the what is being assigned and assign it to me as a
list'.

Here are some examples:

>>> import platform
>>> platform.python_version()
' color="#ff00ff">3.0.1'
>>> a,*b = [1,2,3,4]
>>>  color="#804040">print(a,b, sep=' color="#ff00ff">    ')
1    [2, 3, 4]
>>> a = [1,2,3,4]
>>>  color="#0000ff"># Note the comma after the *b to make it a tuple
>>> *b, = [1,2,3,4]
>>> a == b
True
>>>  color="#0000ff"># Only one * per iterable (or sub-iterable) being assigned
>>> *a,*b = [1,2,3,4]
 style="color: red;">SyntaxError: two starred expressions  style="color: red;" color="#804040">in style="color: red;"> assignment (<pyshell style="color: red;" color="#0000ff">#87>, line 1)
>>>
>>>  color="#0000ff"># but note:
>>> a,*b,(c,*d) = [1,2,3,[4,5,6]]
>>>  color="#804040">print(a,b,c,d, sep=' color="#ff00ff">    ')
1    [2, 3]    4    [5, 6]
>>>  color="#0000ff"># Notice that  style="font-weight: bold;">*  style="text-decoration: underline;">always creates a list
>>> *a, = (1,2,3)
>>> a
 style="font-weight: bold; color: rgb(0, 153, 0);">[1, 2, 3 style="font-weight: bold; color: rgb(0, 153, 0);">]
>>> a,b,(*c,) = (x  color="#804040">if x<2  color="#804040">else range(5)  color="#804040">for x  color="#804040">in range(3))
>>>  color="#804040">print(a,b,c, sep=' color="#ff00ff">    ')
0    1    [0, 1, 2, 3, 4]
>>>
>>>  color="#0000ff"># Although you can have only one  style="font-weight: bold;">* per (sub)iterable, its position is flexible
>>> a,*b,c = [1,2,3,4]
>>> b
[2, 3]
>>> *b,a,c = [1,2,3,4]
>>> b
[1, 2]
>>> a,c,*b = [1,2,3,4]
>>>
>>>  color="#0000ff"># And look!
>>> a,c,*b,d,e = [1,2,3,4]
>>> b
[]
>>> a,b,c,d,e
(1, [], 2, 3, 4)
>>>  color="#0000ff"># So * can accumulate  style="text-decoration: underline;">no values in some cases

I'm looking forward to using this feature in 3.1.

Critique of pseudocode explanations of the Closest Pair Algorithm

2009-05-12T06:30:00.000+01:00

I watched this
task on Rosetta Code from its inception on the 6^th of
May, but had problems understanding the pseudo-code given, and some
of the stuff I googled too. after a couple of days though I was more
familiar with attempts at O(nlog(n)) divide-n-conquer
algorithms having seen a few and happened to find this
explanation (.ppt file) where it seemed to click for me.

After problems with the other references I thought I might go
through and identify what I found problematic with them.

The Rosetta Code Pseudo-code:

closestPair of P(1), P(2), ... P(N)
if N ≤ 3 then
  return closest points of P using brute-force algorithm
else
  xP ← P ordered by the x coordinate, in ascending order
  P_L ← points of xP from 1 to ⌈N/2⌉
  P_R ← points of xP from ⌈N/2⌉+1 to N
  (d_L, pair_L) ← closestPair of P_L

  (d_R, pair_R) ← closestPair of P_R
  (d_min, pair_min) ← (d_R, pair_R)
  if d_L < d_R then

    (d_min, pair_min) ← (d_L, pair_L)
  endif
  x_M ← xP(⌈N/2⌉)_x

  S ← { p ∈ xP : |x_M - p_x| < d_min }
  yP ← S ordered by the y coordinate, in ascending order
  nP ← number of points in yP
  (closest, closestPair) ← (d_min, pair_min)
  for i from 1 to nP - 1
    k ← i + 1
    while k ≤ nP and yP(k)_y - yP(k)_y < d_min

      if |yP(k) - yP(i)| < closest then
        (closest, closestPair) ← (|yP(k) - yP(i)|, {yP(k), yP(i)})
      endif
      k ← k + 1
    endwhile
  return closest, closestPair

endif

The author had noted that the above could have problems so I was
warned, but still started by trying to follow the above. I was
blindly following the pseudo-code until I got the the statement: x_M
← xP(⌈N/2⌉)_x which I couldn't understand.

This prompted my wider search for a better explanation. After
finding it (linked in the first paragraph), I modified my Python
program to follow the new explanation.

The RC references

There were several references quoted.

The Wikipedia article
Did not give pseudo-code
Closest
Pair (McGill)
Section 3.4 on improving their previous
O(nlog²n) into an O(nlogn)
algorithm I found to be very vague. Just how do you “returning
the points in each set in sorted order by y-coordinate “
without doing an nlogn sort leading again to an O(nlog²n)
solution?
Closest
Pair (UCBS)
This seems to concur with my preferred reference
but is given at a higher level than pseudo-code, and so would be
harder to implement from the description given.
Closest
pair (WUStL) by Dr. Jeremy Buhler
Their seems to be a
problem with having N, the number of points as an argument of the
function. What if N is odd? The pseudo-code seems to rely on other
material probably presented in a class and is at a higher level than
my preferred solution as it seems to be set as an exercise for
students.

My Python Solution

The only tests done is to generate random points and compare the
brute force algorithms solution with that of the divide and conquer.

'''

  Compute nearest pair of points using two algorithms
  
  First algorithm is 'brute force' comparison of every possible pair.
  Second, 'divide and conquer', is based on:
    www.cs.iupui.edu/~xkzou/teaching/CS580/Divide-and-conquer-closestPair.ppt 
'''

from random import randint

infinity = float('inf')


# Note the use of complex numbers to represent 2D points making distance == abs(P1-P2)

def bruteForceClosestPair(point):
    numPoints = len(point)
    if numPoints < 2:
        return infinity, (None, None)
    return min( ((abs(point[i] - point[j]), (point[i], point[j]))
                 for i in range(numPoints-1)
                 for j in range(i+1,numPoints)),
                key=lambda x: x[0])


def closestPair(point):
    xP = sorted(point, key= lambda p: p.real)
    yP = sorted(point, key= lambda p: p.imag)
    return _closestPair(xP, yP)

def _closestPair(xP, yP):
    numPoints = len(xP)
    if numPoints <= 3:
        return bruteForceClosestPair(xP)
    Pl = xP[:numPoints/2]
    Pr = xP[numPoints/2:]
    Yl, Yr = [], []
    xDivider = Pl[-1].real
    for p in yP:
        if p.real <= xDivider:
            Yl.append(p)
        else:
            Yr.append(p)
    dl, pairl = _closestPair(Pl, Yl)
    dr, pairr = _closestPair(Pr, Yr)
    dm, pairm = (dl, pairl) if dl < dr else (dr, pairr)
    # Points within dm of xDivider sorted by Y coord

    closeY = [p for p in yP  if abs(p.real - xDivider) < dm]
    numCloseY = len(closeY)
    if numCloseY > 1:
        # There is a proof that you only need compare a max of 7 other points

        closestY = min( ((abs(closeY[i] - closeY[j]), (closeY[i], closeY[j]))
                         for i in range(numCloseY-1)
                         for j in range(i+1,min(i+8, numCloseY))),
                        key=lambda x: x[0])
        return (dm, pairm) if dm <= closestY[0] else closestY
    else:
        return dm, pairm


def times():
    ''' Time the different functions
    '''
    import timeit

    functions = [bruteForceClosestPair, closestPair]
    for f in functions:
        print 'Time for', f.__name__, timeit.Timer(
            '%s(pointList)' % f.__name__,
            'from closestpair import %s, pointList' % f.__name__).timeit(number=1)



pointList = [randint(0,1000)+1j*randint(0,1000) for i in range(2000)]


if __name__ == '__main__':
    pointList = [(5+9j), (9+3j), (2+0j), (8+4j), (7+4j), (9+10j), (1+9j), (8+2j), 10j, (9+6j)]
    print pointList
    print '  bruteForceClosestPair:', bruteForceClosestPair(pointList)
    print '            closestPair:', closestPair(pointList)
    for i in range(10):
        pointList = [randint(0,10)+1j*randint(0,10) for i in range(10)]
        print '\n', pointList
        print ' bruteForceClosestPair:', bruteForceClosestPair(pointList)
        print '           closestPair:', closestPair(pointList)
    print '\n'
    times()
    times()
    times()

Sample output:

[(5+9j), (9+3j), (2+0j), (8+4j), (7+4j), (9+10j), (1+9j), (8+2j), 10j, (9+6j)]
  bruteForceClosestPair: (1.0, ((8+4j), (7+4j)))
            closestPair: (1.0, ((8+4j), (7+4j)))

[(6+8j), (4+1j), 7j, (4+10j), (10+6j), (1+9j), (10+10j), (5+0j), (3+1j), (6+4j)]
 bruteForceClosestPair: (1.0, ((4+1j), (3+1j)))
           closestPair: (1.0, ((3+1j), (4+1j)))

[(5+5j), (10+9j), 3j, (6+1j), (9+5j), 1j, (8+0j), 8j, (3+2j), (10+4j)]
 bruteForceClosestPair: (1.4142135623730951, ((9+5j), (10+4j)))
           closestPair: (1.4142135623730951, ((9+5j), (10+4j)))

[(9+2j), (1+2j), 5j, (10+0j), (6+6j), (6+4j), (5+1j), (2+5j), (8+6j), (3+2j)]
 bruteForceClosestPair: (2.0, ((1+2j), (3+2j)))
           closestPair: (2.0, ((6+6j), (6+4j)))

[(2+3j), (9+0j), (9+7j), 8j, (4+5j), (6+7j), (6+2j), (3+3j), (3+10j), (7+5j)]
 bruteForceClosestPair: (1.0, ((2+3j), (3+3j)))
           closestPair: (1.0, ((2+3j), (3+3j)))

[(7+2j), (2+10j), (3+7j), (8+6j), (2+1j), 3j, (8+5j), (6+9j), (1+0j), 0j]
 bruteForceClosestPair: (1.0, ((8+6j), (8+5j)))
           closestPair: (1.0, ((8+6j), (8+5j)))

[(10+7j), (6+3j), (2+7j), (3+6j), (9+3j), (8+5j), (5+7j), (3+10j), 5j, (3+7j)]
 bruteForceClosestPair: (1.0, ((2+7j), (3+7j)))
           closestPair: (1.0, ((3+6j), (3+7j)))

[(4+6j), (5+4j), 2j, (7+0j), 4j, (6+5j), (7+5j), (8+9j), (10+5j), (8+2j)]
 bruteForceClosestPair: (1.0, ((6+5j), (7+5j)))
           closestPair: (1.0, ((6+5j), (7+5j)))

[(4+0j), (7+9j), 7j, 8j, (7+1j), (2+5j), (7+3j), (1+3j), (3+9j), (10+7j)]
 bruteForceClosestPair: (1.0, (7j, 8j))
           closestPair: (1.0, (7j, 8j))

[(9+1j), (5+0j), (8+3j), (6+9j), (3+7j), (4+6j), (7+6j), (8+6j), (10+8j), (2+5j)]
 bruteForceClosestPair: (1.0, ((7+6j), (8+6j)))
           closestPair: (1.0, ((7+6j), (8+6j)))

[(2+2j), (10+8j), (3+0j), (8+0j), (1+1j), (6+5j), 10j, (6+8j), (10+4j), (4+10j)]
 bruteForceClosestPair: (1.4142135623730951, ((2+2j), (1+1j)))
           closestPair: (1.4142135623730951, ((1+1j), (2+2j)))


Time for bruteForceClosestPair 4.59288093878
Time for closestPair 0.120782948671
Time for bruteForceClosestPair 5.00096353028
Time for closestPair 0.120047939054
Time for bruteForceClosestPair 4.86709875138
Time for closestPair 0.123363723602

Note how much slower the brute force algorithm is for only 2000
points.

(Oh, I guess I should add that the above is all my own work –
students might want to say the same).

On "RailsConf: What killed Smalltalk could kill Ruby, too"

2009-05-10T01:54:00.001+01:00

Browsing Reddit introduced me to RailsConf: What killed Smalltalk could kill Ruby, too The keynote speech by Robert Martin. I found it to be quite interesting as I had been monitoring odd posts on the Rails community and its general attitude.

Some of the points made by Robert on why Smalltalk contracted in popularity included:

It was too easy to make a mess.
And the mess was found late in the project timeline, leading to failed projects.
Smalltakers didn't want to deal with the boring 'corporate' tasks.
Leading to the rise of programmers using other languages to earn a corporate living.
When Smalltalk was 'riding high', a tendency to look down on other languages.
This would alienate the general population of programmers.

Robert put in a big, recurring, plug for both TDD and refactoring IDEs as ways to keep Ruby alive.

One comment he made was that he thought that in the dynamic vs static typing, language wars: Dynamic won and cited as evidence that most mainstream 'static' languages such as Java, C++ and C# have ways to do dynamic type checking – although badly. I wouldn't go that far, as I suspect that little of the mindset of programming in a dynamic language such as Python would be adopted by a C# or Java programmer. I suspect that they will use the dynamism available to them as a last resort as part of a design pattern to solve a particular sub-problem.

Tech Issues For The Euro Elections

2009-05-03T06:17:00.000+01:00

Elections for members of the European Parliament are href="http://news.bbc.co.uk/1/hi/world/europe/7819889.stm"
target="_blank">coming up soon. At the moment, it
seems certain companies are lobbying to have href="http://www.theregister.co.uk/2009/05/02/microsoft_open_source_eu_patents/">software
patents extended into Europe which is something I
don't agree with, as on algorithms - how can you distinguish between
the maths and the code? On something like look-and-feel of a GUI - how
can you distinguish between the idea and its implementation? If you see
that person A used a pinching motion of the fingers to
zoom-out centered on the pinch then I don't think that should stop
person B from coding his own implementation using different code. I do
think that the screen technology that allows a particular
implementation of multi-touch should be patentable however.

There is another tech issue that I need to talk to prospective
candidates about: their stance on href="http://www.theregister.co.uk/2008/04/04/ooxml_ec_investigation_iso/"
target="_blank"> ODF vs OOX ML, and more widely:
their stance on the promotion of the use of open style="text-decoration: underline;">standards
within government.

Can't think of any more tech issues at the moment but if you can think
of any...

P.S. Happy birthday to me. Happy birthday to me,...

Fight for the community you want

2009-04-29T08:20:00.002+01:00

Someone took offence at a puerile and sexist presentation given at what is supposed to be a technical event in the Rails community. Good on them. The community is what you make of it.

"Using sex to sell spanners" doesn't make me think they are good spanners.

- Paddy.

Monitoring a linux process as it reads files

2009-04-21T22:40:00.001+01:00

I am processing 20 thousand files with some proprietary software and
needed to monitor how far it got in reading the files. In my own Python
version of the utility the reading of data was ten times faster than
the subsequent processing and i wanted to find out if this proprietary
solution, which was havinfg performance problems , was equally spending
most of its time reading data.

The proprietary program took a list of 20000 files to process as its
first argument and I remembered that on Linux, the /proc
directory had info on running processes and sure enough , the
/proc/<process id>/fd directory had info on all the file
descriptors currently open by the process as links. So by opening the
list of files to in my editor and searching within it for the file name
shown on one of the file descriptors, I could gauge how many
files been read so far.

I decided to automate the checking and wrote a shell script using
cat/fgrep/gawk/... that then told me what line in the list of files to
process the program was currently at.

Now I've had time to refine things to use mainly python but to
demonstrate its use I also have to generate a test environment

First create some test files to process

style="font-family: monospace;">bash$ style="font-weight: bold;">mkdir -p /tmp/test style="font-family: monospace;">
bash$ style="font-weight: bold;">for ((i=0; i < 100; i++))
do touch /tmp/test/file$i ;done style="font-family: monospace;">
bash$ style="font-weight: bold;">/bin/ls /tmp/test/file* >
/tmp/test/all_files.lst style="font-family: monospace;">
bash$ style="font-weight: bold;">head /tmp/test/all_files.lst style="font-family: monospace;">
/tmp/test/file0 style="font-family: monospace;">
/tmp/test/file1 style="font-family: monospace;">
/tmp/test/file10 style="font-family: monospace;">
/tmp/test/file11 style="font-family: monospace;">
/tmp/test/file12 style="font-family: monospace;">
/tmp/test/file13 style="font-family: monospace;">
/tmp/test/file14 style="font-family: monospace;">
/tmp/test/file15 style="font-family: monospace;">
/tmp/test/file16 style="font-family: monospace;">
/tmp/test/file17 style="font-family: monospace;">
bash$

Now lets create a test executable to monitor

This script just holds each file open for reading for twenty seconds
before closing the file.

style="font-family: monospace;">bash$ style="font-weight: bold;">python -c ' style="color: rgb(0, 0, 153);">import sys,time style="font-family: monospace; color: rgb(0, 0, 153);">
for
name in file(sys.argv[1]): style="font-family: monospace; color: rgb(0, 0, 153);">

f = file(name.strip()) style="font-family: monospace; color: rgb(0, 0, 153);">

time.sleep(45) style="font-family: monospace; color: rgb(0, 0, 153);">

f.close()

' style="font-weight: bold;">/tmp/test/all_files.lst
& style="font-family: monospace;">
[2] style="font-weight: bold;"> style="font-family: monospace; color: rgb(0, 0, 153);"> style="font-weight: bold; color: red;">7984 style="font-family: monospace;">
bash$

here is what the fd directory looks like

style="font-family: monospace;">bash$ ls -l /proc/ style="font-family: monospace; color: rgb(0, 0, 153);"> style="font-weight: bold; color: red;">7984 style="font-family: monospace;">/fd style="font-family: monospace;">
total 0 style="font-family: monospace;">
lrwxrwxrwx 1 HP
DV8025EA None 0 Apr 21 22:17 0 -> /dev/tty1 style="font-family: monospace;">
lrwxrwxrwx 1 HP
DV8025EA None 0 Apr 21 22:17 1 -> /dev/tty1 style="font-family: monospace;">
lrwxrwxrwx 1 HP
DV8025EA None 0 Apr 21 22:17 2 -> /dev/tty1 style="font-family: monospace;">
lrwxrwxrwx 1 HP
DV8025EA None 0 Apr 21 22:17 3 -> /tmp/test/all_files.lst style="font-family: monospace;">
lrwxrwxrwx 1 HP
DV8025EA None 0 Apr 21 22:17 4 -> /tmp/test/file0 style="font-family: monospace;">
bash$

And here is a python script to monitor that fd directories
link number 4 periodically

style="font-family: monospace;">bash$ style="font-weight: bold;">python -c ' style="color: rgb(0, 0, 153);">import sys,time,os,datetime style="color: rgb(0, 0, 153);">
name2index =
dict((name.strip(), index) for index,name in
enumerate(file(sys.argv[1]))) style="color: rgb(0, 0, 153);">
all = len(name2index) style="color: rgb(0, 0, 153);">
while True: style="color: rgb(0, 0, 153);">
path =
os.path.realpath("/proc/7984/fd/ style="font-weight: bold; color: red;">4").strip() style="color: rgb(0, 0, 153);">
print
name2index[path],"/",all, path, datetime.datetime.now().isoformat() style="color: rgb(0, 0, 153);">

time.sleep(30)

' /tmp/test/all_files.lst

22 / 100 /tmp/test/file29 2009-04-21T22:34:07.817750

23 / 100 /tmp/test/file3 2009-04-21T22:34:37.820750

24 / 100 /tmp/test/file30 2009-04-21T22:35:07.825750

24 / 100 /tmp/test/file30 2009-04-21T22:35:37.834750

I watched the monitor output over the next couple of hours and found
out when file reading ended and processing of read data started.

END.

Bimap. Bi-directional mapping/dictionary for Python?

2009-04-16T19:08:00.001+01:00

I had a chance encounter with the boost C++ library documentation and
came across a description of their href="http://cablemodem.fibertel.com.ar/mcape/boost/libs/bimap/index.html">bimap
(that is bimap without a 't') .

It seems to be a bi-directional mapping.

In a Python dict, keys map to values. In a bimap, values would also map
to keys. It seems they put keys/values on an equal footing by renaming
them as left and right members, and having methods that work on either
the right to left or the left to right mapping.

I had just re-visited a program of mine that dealt with what I would
consider large sets of values. I was mapping test names to test
coverage where their were tens of trhousands of test names, each of
~100 characters long and for each file, I had many thousands of code
coverage points, which themselves were strings of ~150 characters each.
I new I had to work with a dict mapping test_names to sets of covered
points for each test and do a lot of set arithmatic to rank the tests
in order of contribution to the overall coverage figure, so I decided
to use indices; replacing files by negative numbers and coverage points
by positive numbers. My program works, and is very fast but I note that
I create both a testname2index mapping dict and the inverse
index2testname mapping, (and coverpoint2index and index2coverpoint
mapping dicts).

What I had been constructing, unawares, was a bimap by using two dicts.

I then spent time googling for Python bimaps to no avail, (although I
think there are haskel and Java implementations out there - and Google
'helps' by sugesting you meant bitmap with a 't').

So my questions are:

Is their a Python bimap implementation out there?

Would anyone want one? (I already code around its
absense).

What would its methods be?

Methods

On my third question, I guess you could have all the normal methods of
a dict but preceeded by either
left_ or right_
with calling of a method without those prefixes raising an exception,
as it seems important to not favour one mapping direction over another.

This would give:

style="font-family: monospace; font-weight: bold; text-decoration: underline;">
DICT      BIMOD
(LEFT) BIMOD (RIGHT) style="font-family: monospace; font-weight: bold;">
style="font-family: monospace;">
clear
left_clear right_clear


copy
left_copy right_copy


fromkeys
left_fromkeys right_fromkeys style="font-family: monospace;">

get
left_get
right_get      style="font-family: monospace;">

has_key
left_has_key right_has_key style="font-family: monospace;">

items
left_items right_items

iteritems    left_iteritems
right_iteritems


iterkeys
left_iterkeys right_iterkeys style="font-family: monospace;">
itervalues
left_itervalues right_itervalues style="font-family: monospace;">

keys
left_keys right_keys


pop
left_pop
right_pop      style="font-family: monospace;">

popitem
left_popitem right_popitem style="font-family: monospace;">
setdefault
left_setdefault right_setdefault style="font-family: monospace;">

update
left_update right_update   style="font-family: monospace;">

values
left_values right_values

The above couldn't be extended for indexing however, unless you did
something like:

style="font-family: monospace;">bimod_instance.left[leftindex]
and bimod_instance..right[rightindex]

but then wouldn't it be better to use the dotted form for all the
methods and have:

style="font-family: monospace;">
bimod_instance.left.haskey # (with a dot after left), ...

I note that their are redundancies above, for example, left_clear ==
right_clear ...

Strewth, I've only scratched the surface :-)

- Paddy.

Knapsack solution by OO Calc

2009-04-10T22:56:00.001+01:00

I had previously solved the href="http://www.rosettacode.org/wiki/Knapsack_Problem"
target="_blank">knapsack problem in Python. I came
across the Solver function of OpenOffice Calc and so tried to use it to
solve knapsack.

A quick recap of the problem

A traveller gets diverted and has to make an unscheduled stop
in
what turns out to be Shangri La. Opting to leave, he is allowed to take
as much as he likes of the following items, so long as it will fit in
his knapsack, and he can carry it.
He knows that he can carry no more than 25 'weights' in total; and that
the capacity of his knapsack is 0.25 'cubic lengths'.

Looking just above the bar codes on the items he finds their
weights and volumes. He digs out his recent copy of a financial paper
and gets the value of each item.

cellpadding="2" cellspacing="2">

nowrap="nowrap" valign="middle">Item	nowrap="nowrap" valign="middle">Explanation	nowrap="nowrap" valign="middle">Value (each)	nowrap="nowrap" valign="middle">weight	nowrap="nowrap" valign="middle">Volume (each)
panacea (vials of)	Incredible healing properties	3000	0.3	0.025
ichor (ampules of)	Vampires blood	1800	0.2	0.015
gold (bars)	Shiney shiney	2500	2.0	0.002
align="left" nowrap="nowrap" valign="middle">Knapsack	align="left" nowrap="nowrap" valign="middle">For the carrying of	align="left" nowrap="nowrap" valign="middle">-	align="left" nowrap="nowrap" valign="middle"><=25	align="left" nowrap="nowrap" valign="middle"><=0.25

He can only take whole units of any item,, but there is much
more of any item than he could ever carry

How many of each item does he take to maximise the
value of items he is carrying away with him?

Note:

There are four solutions that maximise the value taken.
Only one need be given.

OO Calc setup

I am using version 3.0.1. I opened a new spreadsheat and literally
cut-n-pasted the above table into the sheet. The solver needs to be
able to change some cells (the variables) that affect one cell
containing the value to optimise.. I added the table in blue, at the
right of the figure below to calculate the value, weight and
volume of what is in the knapsack based on the amount of each item
chosen, so you can change what is in the style="font-weight: bold;">Number column and
get new totals calculated in the Totals
row

style="width: 914px; height: 444px;"
alt="Knapsack problem formulae" title="(Showing cell formulae)"
src="http://lh3.ggpht.com/_ueV7OWnroeY/Sd-9LlNhkmI/AAAAAAAAAD8/4oXkLjx0xzU/Knapsack_problem_fomulae.PNG">

I show the formulae in the cells in the image above.

The Solver

Menu Tools-> Solver shows a form for filling in. from top to
bottom:

The target cell is the cell I want to maximise, which is
the total value of all items in cell H8

I want to Maximise.

I want to change the Number
of each item i.e. vary cells G4:G6

My limiting conditions are:

The weight is <= 25

The volume is <= 0.25

The Number of each item cannot be negative.

style="width: 455px; height: 373px;" alt="Solver Menu"
title="Solver menu"
src="http://lh6.ggpht.com/_ueV7OWnroeY/Sd-9L2VfwsI/AAAAAAAAAEU/8mlo9NrxqIg/Knapsack_problem_solver.PNG">

You also need to click the Solvers Options
button and ensure the options are set like this:

style="width: 431px; height: 286px;" alt="Solver Options"
title="Solver Options"
src="http://lh6.ggpht.com/_ueV7OWnroeY/Sd-9L2P_QgI/AAAAAAAAAEM/Af-j8UvPRMg/Knapsack_problem_solve_optionsr.PNG">

Ok the Options, then hit Solve on the Solver window and eventually you
will get the following result:

style="width: 914px; height: 444px;"
alt="Result of Knapsack constraints problem"
title="Result of Knapsack constraints problem"
src="http://lh5.ggpht.com/_ueV7OWnroeY/Sd-9Lx1yvmI/AAAAAAAAAEE/j9ZFL0JLe_k/Knapsack_problem_result.PNG">

The solver has correctly determined that carrying no panacea, fifteen
ampules of ichor, and eleven gold bars will maximise the value, subject
to the given constraints..

Comment

Their are actually other values that give the same total value under
the same constraints but I can't find an easy way for the solver to
give them all.

END.

Memoization and stack use.

2009-04-10T09:11:00.000+01:00

I was caught out today when investigating how to apply memoization to
allow computation of larger values of a recursive function.

Mutual Recursion

I had thought up a new task for Rosetta Code called href="http://www.rosettacode.org/wiki/Mutual_Recursion">Mutual
Recursion, as I had remembered that some languages needed
at least the signature of a function to be defined before it
could be called, and so creating mutually recursive functions would
allow you to compare the languages in an interesting way - which is
what R.C. is all about.

For the example for implementation I used href="http://en.wikipedia.org/wiki/Hofstadter_sequence#Hofstadter_Female_and_Male_sequences"
target="_blank">Hoffstadters Female and Male
sequences, and produced the following href="http://www.rosettacode.org/wiki/Mutual_Recursion#Python">Python
definitions:

def  color="#008080">F(n): return 1  color="#804040">if n == 0  color="#804040">else n - M(F(n-1))
 color="#804040">def  color="#008080">M(n): return 0  color="#804040">if n == 0  color="#804040">else n - F(M(n-1))

It seems that other members of R.C. liked the task and soon added many
examples in other languages.

Apart from the computation of the series for 0<=n<20, I
did try to compute F(200) but it blew the stack so I went to bed and
decided to memoize the functions another day.

target="_blank">Wot, no memoize?

Memoization came very easily to mind, so I thought that Python would
have such available as a decorator waiting to be applied. No such luck,
I would have to craft my own . I did remember correctly however that
there was some helper function that made wrapped functions
look a lot more like what had been wrapped, and so used functools.wrap:

from functools  color="#a020f0">import wraps

 color="#804040">def  color="#008080">memoize0(func):
    ' color="#ff00ff"> Adds cache as attribute to wrapped function for ease of access'
     color="#a020f0">@wraps(func)
     color="#804040">def  color="#008080">wrapper(*args):
        argsl = tuple(args)
         color="#804040">if argsl  color="#804040">in wrapper._cache:
             color="#804040">return wrapper._cache[argsl]
         color="#804040">else:
             color="#804040">return wrapper._cache.setdefault(argsl, func(*args))
    wrapper._cache = dict()
     color="#804040">return wrapper

 color="#804040">def  color="#008080">memoize(func):
    ' color="#ff00ff"> Creates closure around locally created cache'
    cache = dict()
     color="#a020f0">@wraps(func)
     color="#804040">def  color="#008080">wrapper(*args):
        argsl = tuple(args)
         color="#804040">if argsl  color="#804040">in cache:
             color="#804040">return cache[argsl]
         color="#804040">else:
             color="#804040">return cache.setdefault(argsl, func(*args))
     color="#804040">return wrapper

I haven't got over my (probably infantile) lust for attaching
attributes to functions so justified memoize0 as it allows you to look
at the contents of the cache. Combined with the fact that you cannot do
a memoization decorator without applying it to a href="http://en.wikipedia.org/wiki/Fibonacci_function">Fibonacci
number function:

if __name__ == ' color="#ff00ff">__main__':
     color="#a020f0">@memoize0
     color="#804040">def  color="#008080">fib(n): return n  color="#804040">if n  color="#804040">in (0,1)  color="#804040">else fib(n - 1) + fib(n - 2)
     color="#a020f0">@memoize0
     color="#804040">def  color="#008080">fib2(n):  color="#804040">return n  color="#804040">if n  color="#804040">in (0,1)  color="#804040">else fib2(n - 1) + fib2(n - 2)
     color="#804040">assert fib._cache  color="#804040">is  color="#804040">not fib2._cache
    fib(50); fib2(50)
     color="#804040">assert fib._cache == fib2._cache  color="#804040">and (fib._cache  color="#804040">is  color="#804040">not fib2._cache)

It is important that caches are not shared for mutually recursive
functions.

A blown stack

I had initial success with the Male and Female sequence and was easily
able to compute F(200) with:


    @ color="#008080">memoize
     color="#804040">def  color="#008080">F(n):
        ' color="#ff00ff"> Hofstadters Female Male sequences'
         color="#804040">return 1  color="#804040">if n == 0  color="#804040">else n - M(F(n-1))

     color="#a020f0">@memoize
     color="#804040">def  color="#008080">M(n):
        ' color="#ff00ff"> Hofstadters Female Male sequences'
         color="#804040">return 0  color="#804040">if n == 0  color="#804040">else n - F(M(n-1))

     color="#804040">assert F(200) == 124

But I got cocky, as F(2000) ran out of stack.

A little thought made me realise that I only had a cached results for F
and M up to around 200 and so a call of F(2000) was going to eat up
huge amounts of stack before recursing down to cached values.

With that knowledge in mind I thought I might gradually jump to higher
values of n and, indeed, the following worked:

    print ([F(n)  color="#804040">for n  color="#804040">in range(0,20000+1,100)])

In real-world applicatrions, you might have to have some scheme to
reduce or limit cache size if memory rather than the stack becomes an
issue, and I have a great solution that I have just completed in the
margin :-)

- Paddy.

(Ab)use of Pythons in-built data structures

2009-04-02T22:07:00.005+01:00

I like to use Pythons in-built data structures quit a lot, and tend to
force myself to ask whether I should create my own classes,
which allows you to use meaningful names for fields and add
comments/docstrings to the datastructure but usually at the cost of
adding more lines of text.

After writing about Huffman codes and posting solutions to href="http://www.rosettacode.org/w/index.php?title=Huffman_codes&oldid=27804#Python">Rosetta
code and href="http://paddy3118.blogspot.com/2009/03/huffman-encoding-in-python.html">my
blog:

 29  color="#804040">def  color="#008080">codecreate(symbol2weights, tutor= False):
 color="#804040"> 30     ''' color="#ff00ff"> Huffman encode the given dict mapping symbols to weights '''
 color="#804040"> 31     global decode
 color="#804040"> 32 
 color="#804040"> 33     heap = [ [float(wt), [[sym, []]], repr(sym)]  color="#804040">for sym, wt  color="#804040">in symbol2weights.iteritems() ]
 color="#804040"> 34     heapify(heap)
 color="#804040"> 35     if tutor:  color="#804040">print " color="#ff00ff">ENCODING:", sorted(symbol2weights.iteritems())
 color="#804040"> 36     while len(heap) >1:
 color="#804040"> 37         lo = heappop(heap)
 color="#804040"> 38         hi = heappop(heap)
 color="#804040"> 39          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  COMBINING:", lo, ' color="#6a5acd">\n        AND:', hi
 color="#804040"> 40          color="#804040">for i  color="#804040">in lo[1]: i[1].insert(0, ' color="#ff00ff">0')
 color="#804040"> 41          color="#804040">for i  color="#804040">in hi[1]: i[1].insert(0, ' color="#ff00ff">1')
 color="#804040"> 42         lohi = [ lo[0] + hi[0] ] + [lo[1] + hi[1]]
 color="#804040"> 43         lohi.append(' color="#ff00ff">(%s if nextbit() else %s)' % (hi[2], lo[2]))
 color="#804040"> 44          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  PRODUCING:", lohi, ' color="#6a5acd">\n'
 color="#804040"> 45         heappush(heap, lohi)
 color="#804040"> 46     wt, codes, decoder = heappop(heap)
 color="#804040"> 47     decode = eval(' color="#ff00ff">lambda :' + decoder, globals())
 color="#804040"> 48     decode.__doc__ = decoder
 color="#804040"> 49     for i  color="#804040">in codes: i[1] = ''.join(i[1])
 color="#804040"> 50     #for i in codes: i[::] = i[:2]
 color="#804040"> 51     return sorted(codes, key= color="#804040">lambda x: (len(x[-1]), x))
 color="#804040"> 52

Someone posted a href="http://www.rosettacode.org/w/index.php?title=Huffman_codes&oldid=27804#Java">
Java solution, that was over three times as long. Instead of
just counting and comparing line counts, I took a closer
look to see what the extra code was for.

import java.util.*;

 color="#2e8b57">abstract  color="#2e8b57">class HuffmanTree  color="#2e8b57">implements Comparable<HuffmanTree> {
     color="#2e8b57">public  color="#2e8b57">int frequency;  color="#0000ff">// the frequency of this tree
     color="#2e8b57">public HuffmanTree( color="#2e8b57">int freq) { frequency = freq; }

     color="#0000ff">// compares on the frequency
     color="#2e8b57">public  color="#2e8b57">int compareTo(HuffmanTree tree) {
         color="#804040">return frequency - tree.frequency;
    }
}

 color="#2e8b57">class HuffmanLeaf  color="#2e8b57">extends HuffmanTree {
     color="#2e8b57">public  color="#2e8b57">char value;  color="#0000ff">// the character this leaf represents

     color="#2e8b57">public HuffmanLeaf( color="#2e8b57">int freq,  color="#2e8b57">char val) {
         color="#2e8b57">super(freq);
        value = val;
    }
}

 color="#2e8b57">class HuffmanNode  color="#2e8b57">extends HuffmanTree {
     color="#2e8b57">public HuffmanTree left, right;  color="#0000ff">// subtrees

     color="#2e8b57">public HuffmanNode(HuffmanTree l, HuffmanTree r) {
         color="#2e8b57">super(l.frequency + r.frequency);
        left = l;
        right = r;
    }
}

 color="#2e8b57">public  color="#2e8b57">class HuffmanCode {
     color="#0000ff">// input is an array of frequencies, indexed by character code
     color="#2e8b57">public  color="#2e8b57">static HuffmanTree buildTree( color="#2e8b57">int[] charFreqs) {
        PriorityQueue<HuffmanTree> trees =  color="#804040">new PriorityQueue<HuffmanTree>();
         color="#0000ff">// initially, we have a forest of leaves
         color="#0000ff">// one for each non-empty character
         color="#804040">for ( color="#2e8b57">int i =  color="#ff00ff">0; i < charFreqs.length; i++)
             color="#804040">if (charFreqs[i] >  color="#ff00ff">0)
                trees.offer( color="#804040">new HuffmanLeaf(charFreqs[i], ( color="#2e8b57">char)i));

         color="#804040">assert trees.size() >  color="#ff00ff">0;
         color="#0000ff">// loop until there is only one tree left
         color="#804040">while (trees.size() >  color="#ff00ff">1) {
             color="#0000ff">// two trees with least frequency
            HuffmanTree a = trees.poll();
            HuffmanTree b = trees.poll();

             color="#0000ff">// put into new node and re-insert into queue
            trees.offer( color="#804040">new HuffmanNode(a, b));
        }
         color="#804040">return trees.poll();
    }

     color="#2e8b57">public  color="#2e8b57">static  color="#2e8b57">void printCodes(HuffmanTree tree, Stack<Character> prefix) {
         color="#804040">assert tree !=  color="#ff00ff">null;
         color="#804040">if (tree  color="#804040">instanceof HuffmanLeaf) {
            HuffmanLeaf leaf = (HuffmanLeaf)tree;

             color="#0000ff">// print out character and frequency
            System.out.print(leaf.value +  color="#ff00ff">"\t color="#ff00ff">" + leaf.frequency +  color="#ff00ff">"\t color="#ff00ff">");

             color="#0000ff">// print out code for this leaf, which is just the prefix
             color="#804040">for ( color="#2e8b57">char bit : prefix)
                System.out.print(bit);
            System.out.println();

        }  color="#804040">else  color="#804040">if (tree  color="#804040">instanceof HuffmanNode) {
            HuffmanNode node = (HuffmanNode)tree;

             color="#0000ff">// traverse left
            prefix.push( color="#ff00ff">'0');
            printCodes(node.left, prefix);
            prefix.pop();

             color="#0000ff">// traverse right
            prefix.push( color="#ff00ff">'1');
            printCodes(node.right, prefix);
            prefix.pop();
        }
    }

     color="#2e8b57">public  color="#2e8b57">static  color="#2e8b57">void main(String[] args) {
        String test =  color="#ff00ff">"this is an example for huffman encoding";

         color="#0000ff">// we will assume that all our characters will have
         color="#0000ff">// code less than 256, for simplicity
         color="#2e8b57">int[] charFreqs =  color="#804040">new  color="#2e8b57">int[ color="#ff00ff">256];
         color="#0000ff">// read each character and record the frequencies
         color="#804040">for ( color="#2e8b57">char c : test.toCharArray())
            charFreqs[c]++;

         color="#0000ff">// build tree
        HuffmanTree tree = buildTree(charFreqs);

         color="#0000ff">// print out results
        System.out.println( color="#ff00ff">"SYMBOL\t color="#ff00ff">WEIGHT\t color="#ff00ff">HUFFMAN CODE");
        printCodes(tree,  color="#804040">new Stack<Character>());
    }
}

I
noticed that the Java was tidy, and that they had defined their own
classes for a HuffmanLeaf structure used when constructing a
HuffmanTree.

In my original Python, I used a nested list as
the equivalent to the HuffmanLeaf of the Java example. It worked, but I
had to introduce 'magic constants' to access the fields of the Python
leaf , which is also un-named as a structure in the program (lline 33
of
the Python code creates a list of Leaf structures):

33  heap = [ [float(wt), [sym, []]]  color="#804040">for sym, wt  color="#804040">in symbol2weights.iteritems() ]

Now I didn't want to go to the trouble of creating my own
classes, but that's were the new href="http://docs.python.org/dev/py3k/library/collections.html#collections.namedtuple">namedtuple
class factory of Python 2.6 came to the rescue.

namedtuple

namedtuple
will generate a subtype of the tuple class that allows the fields of
the tuple to be named and accessed via subscription, [], or as if the
fields are instance variable names.

I modified my Python program to use two named tuples:

For the Leaf structure as a whole in line 3.

For a component of the Leaf structure that holds the symbol
and the Huffman code (accumulated so far), in line 4.

 1  color="#a020f0">from collections  color="#a020f0">import namedtuple
 color="#804040"> 2 
 color="#804040"> 3 Leaf = namedtuple(' color="#ff00ff">Leaf', 'weight, symbols')
 color="#804040"> 4 SH = namedtuple(' color="#ff00ff">SH', 'sym, huff')
 color="#804040"> 5 
 color="#804040"> 6 def  color="#008080">codecreate2(symbol2weights, tutor= False):
 color="#804040"> 7     ''' color="#ff00ff"> Huffman codecreate2 the given dict mapping symbols to weights '''
 color="#804040"> 8     heap = [ Leaf(weight=float(wt), symbols=[ SH(sym, []) ])
 color="#804040"> 9               color="#804040">for sym, wt  color="#804040">in symbol2weights.iteritems() ]
 color="#804040">10     heapify(heap)
 color="#804040">11     if tutor:  color="#804040">print " color="#ff00ff">ENCODING:", sorted(symbol2weights.iteritems())
 color="#804040">12     while len(heap) >1:
 color="#804040">13         lo = heappop(heap)
 color="#804040">14         hi = heappop(heap)
 color="#804040">15          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  COMBINING:", lo, ' color="#6a5acd">\n        AND:', hi
 color="#804040">16          color="#804040">for sh  color="#804040">in lo.symbols: sh.huff.insert(0, ' color="#ff00ff">0')
 color="#804040">17          color="#804040">for sh  color="#804040">in hi.symbols: sh.huff.insert(0, ' color="#ff00ff">1')
 color="#804040">18         lohi = Leaf(weight  = lo.weight + hi.weight,
 color="#804040">19                     symbols = lo.symbols + hi.symbols)
 color="#804040">20          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  PRODUCING:", lohi, ' color="#6a5acd">\n'
 color="#804040">21         heappush(heap, lohi)
 color="#804040">22     symbols = heappop(heap).symbols
 color="#804040">23     symbols = [SH(sym, ''.join(huff))  color="#804040">for sym, huff  color="#804040">in symbols]
 color="#804040">24     return sorted(symbols, key= color="#804040">lambda sh: (len(sh.huff), sh))
 color="#804040">25

The
class generation is succinct in lines 3 and 4, and I use the field
names for clarity for example when creating the original list of Leaves
that will form the heap in line 8

Printing namedtuples

The print statements stay the same, but before you got output like this:

style="font-family: monospace;"> style="font-family: monospace;">...

COMBINING: [2.5, ['C', []]] style="font-family: monospace;">

AND: [5.0, ['A', []]]

PRODUCING: [7.5, ['C', ['0']], ['A', ['1']]] style="font-family: monospace;">
...

Now you get the clearer:

style="font-family: monospace;"> ENCODING:
[('A', '5'), ('B', '25'), ('C', '2.5'), ('D', '12.5')]

COMBINING: Leaf(weight=2.5, symbols=[SH(sym='C', huff=[])])


AND: Leaf(weight=5.0, symbols=[SH(sym='A', huff=[])])

PRODUCING: Leaf(weight=7.5, symbols=[SH(sym='C',
huff=['0']), SH(sym='A', huff=['1'])])

COMBINING: Leaf(weight=7.5, symbols=[SH(sym='C',
huff=['0']), SH(sym='A', huff=['1'])])


AND: Leaf(weight=12.5, symbols=[SH(sym='D', huff=[])])

PRODUCING: Leaf(weight=20.0, symbols=[SH(sym='C', huff=['0', '0']),
SH(sym='A', huff=['0', '1']), SH(sym='D', huff=['1'])])

COMBINING: Leaf(weight=20.0, symbols=[SH(sym='C', huff=['0', '0']),
SH(sym='A', huff=['0', '1']), SH(sym='D', huff=['1'])])


AND: Leaf(weight=25.0, symbols=[SH(sym='B', huff=[])])

PRODUCING: Leaf(weight=45.0, symbols=[SH(sym='C', huff=['0', '0',
'0']), SH(sym='A', huff=['0', '0', '1']), SH(sym='D', huff=['0', '1']),
SH(sym='B', huff=['1'])])

Conclusion

Although
you can do so much with Python lists, if each position in the list has
a fixed meaning then you might be best to use tuples, and if you should
be using tuples, then namedtuples can make your code more readable
without bloating it with long class definitions.

- Paddy.

Huffman Encoding in Python

2009-03-28T00:11:00.000Z

Huffman encoding came up on href="http://www.rosettacode.org/wiki/Huffman_codes">Rosetta
Code.

Huffman encoding is a way to assign binary codes to symbols
that
reduces the overall number of bits used to encode a typical string of
of those symbols.

For example, if you use letters as symbols and have details of
the frequency of occurence of those letters in typical strings, then
you could just encode each letter with a fixed number of bits, such as
in ASCII codes. You can do better than this by encoding more frequently
occurring letters such as e and a, with smaller bit strings; and less
frequently occurring letters such as q and x with longer bit strings.

Any string of letters will be encoded as a string of bits that
are no-longer of the same length per letter. To successfully decode
such as string, the smaller codes assigned to letters such as 'e'
cannot occur as a prefix in the larger codes such as that for 'x'.

If you were to assign a code 01 for 'e' and code 011 for
'x',
then if the bits to decode started as 011... then you would not know
iif you should decode an 'e' or an 'x'.

The Huffman coding scheme takes each symbol and its weight (or
frequency of occurrence), and generates proper encodings for each
symbol taking account of the weights of each symbol, so that higher
weighted symbols have less bits in their encoding. (See the href="http://en.wikipedia.org/wiki/Huffman_coding" class="extiw"
title="wp:Huffman_coding">WP article for more
information).

A Huffman encoding can be computed by first creating a tree of
nodes:

cellpadding="1" cellspacing="1">

Create a leaf node for each symbol and add it to the
priority queue.

While there is more than one node in the queue:
1. Remove the node of highest priority (lowest
  probability) twice to get two nodes.
2. Create a new internal node with these two nodes
  as children
  and with probability equal to the sum of the two nodes' probabilities.
3. Add the new node to the queue.

The remaining node is the root node and the tree is
complete.

Traverse the constructed binary tree from root to leaves
assigning
and accumulating a '0' for one branch and a '1' for the other at each
node. The accumulated zeroes and ones at each leaf constitute a Huffman
encoding for those symbols and weights.

href="http://www.rosettacode.org/wiki/Image:Huffman_coding_example.jpg"
class="image" title="Huffman coding example.jpg"> alt=""
src="http://www.rosettacode.org/w/images/thumb/8/8b/Huffman_coding_example.jpg/250px-Huffman_coding_example.jpg"
border="0" height="120" width="250">

In Python

I orginally gave an an example that matched a definition that
was later found to be insufficient, so substituted my own definition
above.. My href="http://www.rosettacode.org/w/index.php?title=Huffman_codes&oldid=26720#Python">first
Python solution on RC to the wrong definition, did
have the advantage, (as I saw it), of not having to traverse a tree.

My 'true' Huffman code creator assembles each symbol and its
weight into the following structure initially (the style="font-weight: bold;">leaf structure):

style="font-family: monospace;">[ weight, [ symbol, []]]

The weight applies to every (in this case only one), of the style="font-family: monospace;"> [symbol, []]
pairs after it in the same list.

The empty list is used to accumulate the Huffman code for the symbol as
we manipulate the heap, without having to walk a constructed tree
structure.

There are two types of input to the program that I am running examples
with:

A string of space separated symbol, weight pairs, as used
in small examples.

A sample of text for which letters and letter frequencies
are extracted.

The if statement at line 23
allows me to switch between the two types of input whilst
exploring the algorithm.

The tutor argument to the encode function shows what is happening in
the loop around the heap pops

 1 
 color="#804040"> 2 from heapq  color="#a020f0">import heappush, heappop, heapify
 color="#804040"> 3 
 color="#804040"> 4 def  color="#008080">codecreate(symbol2weights, tutor= False):
 color="#804040"> 5     ''' color="#ff00ff"> Huffman encode the given dict mapping symbols to weights '''
 color="#804040"> 6     heap = [ [float(wt), [sym, []]]  color="#804040">for sym, wt  color="#804040">in symbol2weights.iteritems() ]
 color="#804040"> 7     heapify(heap)
 color="#804040"> 8     if tutor:  color="#804040">print " color="#ff00ff">ENCODING:", sorted(symbol2weights.iteritems())
 color="#804040"> 9     while len(heap) >1:
 color="#804040">10         lo = heappop(heap)
 color="#804040">11         hi = heappop(heap)
 color="#804040">12          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  COMBINING:", lo, ' color="#6a5acd">\n        AND:', hi
 color="#804040">13          color="#804040">for i  color="#804040">in lo[1:]: i[1].insert(0, ' color="#ff00ff">0')
 color="#804040">14          color="#804040">for i  color="#804040">in hi[1:]: i[1].insert(0, ' color="#ff00ff">1')
 color="#804040">15         lohi = [ lo[0] + hi[0] ] + lo[1:] + hi[1:]
 color="#804040">16          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  PRODUCING:", lohi, ' color="#6a5acd">\n'
 color="#804040">17         heappush(heap, lohi)
 color="#804040">18     codes = heappop(heap)[1:]
 color="#804040">19     for i  color="#804040">in codes: i[1] = ''.join(i[1])
 color="#804040">20     return sorted(codes, key= color="#804040">lambda x: (len(x[-1]), x))
 color="#804040">21 
 color="#804040">22 # Input types
 color="#804040"> style="font-weight: bold; background-color: rgb(255, 255, 102);">23  color="#804040">if 1:
 color="#804040">24     readin = " color="#ff00ff">B 25   C 2.5 D  12.5 A 5  color="#6a5acd">\n"
 color="#804040">25     #readin = "a .1 b .15 c .3 d .16 e .29" # Wikipedia sample
 color="#804040">26     #readin = "a1 .4 a2 .35 a3 .2 a4 .05" # Wikipedia sample
 color="#804040">27     #readin = "A 50 B 25 C 12.5 D 12.5" # RC example
 color="#804040">28 
 color="#804040">29     cleaned = readin.strip().split()
 color="#804040">30     symbol2weights = dict((symbol, wt)
 color="#804040">31                           color="#804040">for symbol, wt  color="#804040">in zip(cleaned[0::2], cleaned[1::2]) )
 color="#804040">32 else:
 color="#804040">33     astring = " color="#ff00ff">this is an example for huffman encoding"
 color="#804040">34     symbol2weights = dict((ch, astring.count(ch))  color="#804040">for ch  color="#804040">in set(astring))  color="#0000ff"># for astring
 color="#804040">35 
 color="#804040">36 huff = codecreate(symbol2weights, True)
 color="#804040">37 print " color="#6a5acd">\nSYMBOL color="#6a5acd">\tWEIGHT color="#6a5acd">\tHUFFMAN CODE"
 color="#804040">38 for h  color="#804040">in huff:
 color="#804040">39     print " color="#ff00ff">%s\t color="#ff00ff">%s\t color="#ff00ff">%s" % (h[0], symbol2weights[h[0]], h[1])

A run, with the tutor enabled gives the following output:

style="font-family: monospace;">ENCODING: [('A', '5'), ('B',
'25'), ('C', '2.5'), ('D', '12.5')] style="font-family: monospace;">

COMBINING: [2.5, ['C', []]] style="font-family: monospace;">

AND: [5.0, ['A', []]]

PRODUCING: [7.5, ['C', ['0']], ['A', ['1']]] style="font-family: monospace;">

COMBINING: [7.5, ['C', ['0']], ['A', ['1']]] style="font-family: monospace;">

AND: [12.5, ['D', []]] style="font-family: monospace;">

PRODUCING: [20.0, ['C', ['0', '0']], ['A', ['0', '1']], ['D', ['1']]] style="font-family: monospace;">

COMBINING: [20.0, ['C', ['0', '0']], ['A', ['0', '1']], ['D', ['1']]] style="font-family: monospace;">

AND: [25.0, ['B', []]] style="font-family: monospace;">

PRODUCING: [45.0, ['C', ['0', '0', '0']], ['A', ['0', '0', '1']], ['D',
['0', '1']], ['B', ['1']]] style="font-family: monospace;">

SYMBOL
WEIGHT    HUFFMAN CODE style="font-family: monospace;">
B
25
1 style="font-family: monospace;">
D
12.5    01 style="font-family: monospace;">
A
5
001 style="font-family: monospace;">
C
2.5
000

Encode/Decode Round-tripping

I realised that I could use a method similar to how I accumulate the
codes in the heap loop, to generate a single function that can
recognise a single symbol from the beginning of the encoded symbols. By
using the function in a loop, I could regenerate the symbol list.

In the codecreate function, the leaf
structure is modified:

style="font-family: monospace;">[weight, [ [symbol, []] ],
repr(sym)]

Their is an extra level of list around the [symbol, code accumulation
list pair] as well as a new item: 'repr(sym)' it is the third
item in the outer list and will always be the function to
generate the item as accumulated so far.. This function starts off by
returning just the symbol and an outer if/then/else expression is added
as we go around the heap loop (see line style="font-weight: bold;">43)

After the outer expression is accumulated, it is turned into a lambda
expression, the string eval'd, and assigned to a global variable (see
line 47)

I use function probchoice (from earlier RC work), to create an
arbitrary sequence of symbols to in the given weighting then encode and
decode it as well as giving some stats on space saving.

  1 
 color="#804040">  2 from heapq  color="#a020f0">import heappush, heappop, heapify
 color="#804040">  3 import random, bisect
 color="#804040">  4 
 color="#804040">  5 
 color="#804040">  6 # Helper routine for generating test sequences
 color="#804040"> style="background-color: rgb(255, 255, 102);">  7  color="#804040">def  color="#008080">probchoice(items, probs):
 color="#804040">  8   ''' color="#6a5acd">\
 color="#804040">  9   Splits the interval 0.0-1.0 in proportion to probs
 color="#804040"> 10   then finds where each random.random() choice lies
 color="#804040"> 11   (This routine, probchoice, was released under the
 color="#804040"> 12   GNU Free Documentation License 1.2)
 color="#804040"> 13   '''
 color="#804040"> 14 
 color="#804040"> 15   prob_accumulator = 0
 color="#804040"> 16   accumulator = []
 color="#804040"> 17   for p  color="#804040">in probs:
 color="#804040"> 18     prob_accumulator += p
 color="#804040"> 19     accumulator.append(prob_accumulator)
 color="#804040"> 20 
 color="#804040"> 21   while True:
 color="#804040"> 22     r = random.random()
 color="#804040"> 23     yield items[bisect.bisect(accumulator, r)]
 color="#804040"> 24 
 color="#804040"> 25 
 color="#804040"> 26 # placeholder
 color="#804040"> 27 decode =  color="#804040">lambda : None
 color="#804040"> 28 
 color="#804040"> 29 def  color="#008080">codecreate(symbol2weights, tutor= False):
 color="#804040"> 30     ''' color="#ff00ff"> Huffman encode the given dict mapping symbols to weights '''
 color="#804040"> 31     global decode
 color="#804040"> 32 
 color="#804040"> 33     heap = [ [float(wt), [[sym, []]], repr(sym)]  color="#804040">for sym, wt  color="#804040">in symbol2weights.iteritems() ]
 color="#804040"> 34     heapify(heap)
 color="#804040"> 35     if tutor:  color="#804040">print " color="#ff00ff">ENCODING:", sorted(symbol2weights.iteritems())
 color="#804040"> 36     while len(heap) >1:
 color="#804040"> 37         lo = heappop(heap)
 color="#804040"> 38         hi = heappop(heap)
 color="#804040"> 39          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  COMBINING:", lo, ' color="#6a5acd">\n        AND:', hi
 color="#804040"> 40          color="#804040">for i  color="#804040">in lo[1]: i[1].insert(0, ' color="#ff00ff">0')
 color="#804040"> 41          color="#804040">for i  color="#804040">in hi[1]: i[1].insert(0, ' color="#ff00ff">1')
 color="#804040"> 42         lohi = [ lo[0] + hi[0] ] + [lo[1] + hi[1]]
 color="#804040"> style="background-color: rgb(255, 255, 102);"> 43         lohi.append(' color="#ff00ff">(%s if nextbit() else %s)' % (hi[2], lo[2]))
 color="#804040"> 44          color="#804040">if tutor:  color="#804040">print " color="#ff00ff">  PRODUCING:", lohi, ' color="#6a5acd">\n'
 color="#804040"> 45         heappush(heap, lohi)
 color="#804040"> 46     wt, codes, decoder = heappop(heap)
 color="#804040"> style="background-color: rgb(255, 255, 102);"> 47     decode = eval(' color="#ff00ff">lambda :' + decoder, globals())
 color="#804040"> 48     decode.__doc__ = decoder
 color="#804040"> 49     for i  color="#804040">in codes: i[1] = ''.join(i[1])
 color="#804040"> 50     #for i in codes: i[::] = i[:2]
 color="#804040"> 51     return sorted(codes, key= color="#804040">lambda x: (len(x[-1]), x))
 color="#804040"> 52 
 color="#804040"> 53 # Input types
 color="#804040"> style="background-color: rgb(255, 255, 102);"> 54  color="#804040">if 1:
 color="#804040"> 55     tutor = True
 color="#804040"> 56     sequencecount = 50
 color="#804040"> 57     readin = " color="#ff00ff">B 25   C 2.5 D  12.5 A 5  color="#6a5acd">\n"
 color="#804040"> 58     #readin = "a .1 b .15 c .3 d .16 e .29" # Wikipedia sample
 color="#804040"> 59     #readin = "a1 .4 a2 .35 a3 .2 a4 .05" # Wikipedia sample
 color="#804040"> 60     #readin = "A 50 B 25 C 12.5 D 12.5" # RC example
 color="#804040"> 61 
 color="#804040"> 62     cleaned = readin.strip().split()
 color="#804040"> 63     symbol2weights = dict((symbol, wt)
 color="#804040"> 64                           color="#804040">for symbol, wt  color="#804040">in zip(cleaned[0::2], cleaned[1::2]) )
 color="#804040"> 65 else:
 color="#804040"> 66     tutor = False
 color="#804040"> 67     sequencecount = 500
 color="#804040"> 68     astring = " color="#ff00ff">this is an example for huffman encoding"
 color="#804040"> 69     symbol2weights = dict((ch, astring.count(ch))  color="#804040">for ch  color="#804040">in set(astring))  color="#0000ff"># for astring
 color="#804040"> 70 
 color="#804040"> 71 huff = codecreate(symbol2weights, tutor= tutor)
 color="#804040"> 72 print " color="#6a5acd">\nSYMBOL color="#6a5acd">\tWEIGHT color="#6a5acd">\tHUFFMAN CODE"
 color="#804040"> 73 for h  color="#804040">in huff:
 color="#804040"> 74     print " color="#ff00ff">%s\t color="#ff00ff">%s\t color="#ff00ff">%s" % (h[0], symbol2weights[h[0]], h[1])
 color="#804040"> 75 
 color="#804040"> 76 ##
 color="#804040"> 77 ## encode-decode check
 color="#804040"> 78 ##
 color="#804040"> 79 symbol2code = dict(huff)
 color="#804040"> 80 symbols, weights = zip(*symbol2weights.iteritems())
 color="#804040"> 81 # normalize weights
 color="#804040"> 82 weights = [float(wt)  color="#804040">for wt  color="#804040">in weights]
 color="#804040"> 83 tot = sum(weights)
 color="#804040"> 84 weights = [wt/tot  color="#804040">for wt  color="#804040">in weights]
 color="#804040"> 85 # Generate a sequence
 color="#804040"> 86 nxt = probchoice(symbols, weights).next
 color="#804040"> 87 symbolsequence = [nxt()  color="#804040">for i  color="#804040">in range(sequencecount)]
 color="#804040"> 88 # encode it
 color="#804040"> 89 bitsequence = ''.join(symbol2code[sym]  color="#804040">for sym  color="#804040">in symbolsequence)
 color="#804040"> 90 
 color="#804040"> 91 sslen, slen, blen = len(symbolsequence), len(symbols), len(bitsequence)
 color="#804040"> 92 countlen = len(bin(slen-1)[2:])
 color="#804040"> 93 print '''
 color="#804040"> 94 
 color="#804040"> 95 
 color="#804040"> 96 ROUND-TRIPPING
 color="#804040"> 97 ==============
 color="#804040"> 98 I have generated a random sequence of %i symbols to the given weights.
 color="#804040"> 99 If I use a binary count to encode each of the %i symbols I would need
 color="#804040">100 %i * %i = %i bits to encode the sequence.
 color="#804040">101 Using the Huffman code, I need only %i bits.
 color="#804040">102 ''' % (sslen, slen, sslen, countlen, sslen * countlen, blen )
 color="#804040">103 
 color="#804040">104 ## decoding
 color="#804040">105 nextbit = (bit==' color="#ff00ff">1' for bit  color="#804040">in bitsequence).next
 color="#804040">106 
 color="#804040">107 decoded = []
 color="#804040">108 try:
 color="#804040">109     while 1:
 color="#804040">110         decoded.append(decode())
 color="#804040">111 except StopIteration:
 color="#804040">112     pass
 color="#804040">113 
 color="#804040">114 print " color="#ff00ff">Comparing the decoded sequence with the original I get:", decoded == symbolsequence

This short run is in tutor mode, so you can track the accumulation of
the decode function:

style="font-family: monospace;">ENCODING: [('A', '5'), ('B',
'25'), ('C', '2.5'), ('D', '12.5')] style="font-family: monospace;">

COMBINING: [2.5, [['C', []]], "'C'"]


AND: [5.0, [['A', []]], "'A'"] style="font-family: monospace;">

PRODUCING: [7.5, [['C', ['0']], ['A', ['1']]], style="font-weight: bold;">"('A' if nextbit() else 'C')"]

COMBINING: [7.5, [['C', ['0']], ['A', ['1']]], "('A' if nextbit() else
'C')"]


AND: [12.5, [['D', []]], "'D'"] style="font-family: monospace;">

PRODUCING: [20.0, [['C', ['0', '0']], ['A', ['0', '1']], ['D', ['1']]],
"('D' if nextbit() else
('A' if nextbit() else 'C'))"] style="font-family: monospace;">

COMBINING: [20.0, [['C', ['0', '0']], ['A', ['0', '1']], ['D', ['1']]],
"('D' if nextbit() else ('A' if nextbit() else 'C'))"] style="font-family: monospace;">

AND: [25.0, [['B', []]], "'B'"] style="font-family: monospace;">

PRODUCING: [45.0, [['C', ['0', '0', '0']], ['A', ['0', '0', '1']],
['D', ['0', '1']], ['B', ['1']]], "('B'
if nextbit() else ('D' if nextbit() else ('A' if nextbit() else 'C')))"]

SYMBOL
WEIGHT    HUFFMAN CODE style="font-family: monospace;">
B
25    style="font-family: monospace;">
   style="font-family: monospace;">1 style="font-family: monospace;">
D

style="font-family: monospace;">12.5
   01 style="font-family: monospace;">
A

style="font-family: monospace;">5
001 style="font-family: monospace;">
C

style="font-family: monospace;">2.5
000 style="font-family: monospace;">

ROUND-TRIPPING style="font-family: monospace;">
============== style="font-family: monospace;">
I have generated a
random sequence of 50 symbols to the given weights. style="font-family: monospace;">
If I use a binary
count to encode each of the 4 symbols I would need style="font-family: monospace;">
50 * 2 = style="font-weight: bold;">100 bits to encode
the sequence.

Using the Huffman
code, I need only 90
bits.

Comparing the
decoded sequence with the original I get: style="font-weight: bold;">True

And if I change line 54
to be False, I get the following:

style="font-family: monospace;">SYMBOL
WEIGHT    HUFFMAN CODE style="font-family: monospace;">

6    101 style="font-family: monospace;">
n
4    010 style="font-family: monospace;">
a
3    1001 style="font-family: monospace;">
e
3    1100 style="font-family: monospace;">
f
3    1101 style="font-family: monospace;">
h
2    0001 style="font-family: monospace;">
i
3    1110 style="font-family: monospace;">
m
2    0010 style="font-family: monospace;">
o
2    0011 style="font-family: monospace;">
s
2    0111 style="font-family: monospace;">
g
1    00000 style="font-family: monospace;">
l
1    00001 style="font-family: monospace;">
p
1    01100 style="font-family: monospace;">
r
1    01101 style="font-family: monospace;">
t
1    10000 style="font-family: monospace;">
u
1    10001 style="font-family: monospace;">
x
1    11110 style="font-family: monospace;">
c
1    111110 style="font-family: monospace;">
d
1    111111 style="font-family: monospace;">

ROUND-TRIPPING style="font-family: monospace;">
============== style="font-family: monospace;">
I have generated a
random sequence of 500 symbols to the given weights. style="font-family: monospace;">
If I use a binary
count to encode each of the 19 symbols I would need style="font-family: monospace;">
500 * 5 = 2500 bits
to encode the sequence. style="font-family: monospace;">
Using the Huffman
code, I need only 2012 bits. style="font-family: monospace;">

Comparing the
decoded sequence with the original I get: True style="font-family: monospace;">

Note: The purpose of the program is to teach me more about Huffman
coding and is not an exercise in speed!

I am the author of all the Python on this page. The diagram is from
Wikipedia. Please refrain from passing-off my code as your own (that
one is mainly for students).