Day 88

Day 88: The moment of truth

Benchmarking. The time has come. This evening I will test the performance of my implementation and update this post as I get results with each test. The machine I'm going to use is a shared-memory multicore with 8 Intel(R) Xeon(R) X5450 CPU's, 3GHz, 6MB cache, 4 cores each. This little monster is called Jupiter III and together with Jupiters I and II forms a 3-piece cluster.

So, at the moment I am testing the test-bench and soon will be ready to launch. Don't switch the channel=)

[Update]

Phew! Preparation for testing took me way more time than I expected... However, look at this beautiful script! (./topographica line's were created by Chris, everything else is mine)

#!/bin/bash
# Argument =  -c cortex density  -n number of cpu's  -i number of iterations

usage()
{
cat << EOF
usage: $0 options

This script tests the performance of parallel Topographica implementation 
    
OPTIONS:
    -c        cortex density
    -n        number of cpu's
    -i        number of iterations
EOF
}
        
C_DENSITY=
N_CPU=
ITERATIONS=
while getopts “c:n:i:” OPTION
do
    case $OPTION in
        c)
        C_DENSITY=$OPTARG
        ;;
        n)
        N_CPU=$OPTARG
        ;;
        i)
        ITERATIONS=$OPTARG
        ;;
        ?)
        usage
        exit
        ;;
    esac
done

if [[ -z $C_DENSITY ]] || [[ -z $N_CPU ]] || [[ -z $ITERATIONS ]]
then
    usage
    exit 1
fi

# Mock run to make sure that all C code has compiled before proceeding to actual testing
./topographica -p cortex_density=1 examples/tiny.ty -c "topo.sim.run(1)"

# Serial run
./topographica -p cortex_density=$C_DENSITY examples/tiny.ty -c "import timeit; print 'Serial:', timeit.Timer('topo.sim.run($ITERATIONS)','import topo').timeit(number=1)" -c 'import pickle; pickle.dump(topo.sim["V1"].activity,open("results.pickle","wb"))'


# Parallel run
mpirun -n $N_CPU ./topographica -p mpi=True -p cortex_density=$C_DENSITY examples/tiny.ty -c "import timeit; print 'MPI:', timeit.Timer('topo.sim.run($ITERATIONS)','import topo').timeit(number=1)" -c 'import pickle; from numpy.testing import assert_array_equal; previous_result=pickle.load(open("results.pickle","r")); assert_array_equal(topo.sim["V1"].activity,previous_result)'

Now to do testing you just need to run this script with 3 parameters: -c for cortex density, i for the number of iterations and -n for the number of processes to use. I already tried this out on my machine and so far the results are not that inspiring:

./test_script -c 25 -n 2 -i 100
Serial: 5.36015510559
MPI: 6.15542793274

On the bright side, the results computed in MPI are still correct. Now I'm going to try this out on Jupiter.

[update 2] Ok, found a bug in the script: it ran tiny.ty instead of tiny_mpi.ty. How the hell did it get correct results? Anyway, new output (my laptop):
./test_script -c 25 -n 2 -i 100
Not using MPI
Not using MPI
Serial: 1.02229189873
Using MPI
NODE 0 initialized
NODE 1 initialized
MPI: 3.37540221214

[update 3] Ok, I need to make a break now. For some reason Jupiter won't let me run the tests, so I'll have to deal with that somehow. But later.

Day 83

Day 83: ssh me!

Hi, Mr Internet. Yesterday I was finally able to run my first experiment with PMI implementation! This means that now I can start benchmarking and see where I really am. At the moment I am working on setting up my first cluster network - something I started over Christmas Holidays, spent a couple of days working on, hadn't managed to make it and decided to shift my attention to more important things. However, now is the time. I want to see for myself that it is possible to make a scalable network of arbitrary PC's (not specialized hardware) with different architectures to work as a single unit doing computations in parallel and then gathering results on a single node. Currently setting up ssh-agent to allow password-less access from my computer to DICE (university) network, which is absolutely necessary for running parallel jobs between these machines, is giving me the most headache. Anyway, should be manageable.

Also, since it's finally woring now, here's the (part of) actual implementation of MPI_CFProjection:

class MPI_CFProjection(CFProjection):
    
    def __init__(self,initialize_cfs=True, **params):
        pmi.execfile_('topo/base/mpi_cf.py')
        self.pmiobj = pmi.create('MPI_CFProjection_node')
        super(MPI_CFProjection,self).__init__(initialize_cfs=True,**params)
        self.allow_skip_non_responding_units = True # = self.dest.allow_skip_non_responding_units
        
        self.mask = self.dest.mask
            
    
    """>>>>>>>>>>>>>>>>>>>>>>>>> PROPERTIES >>>>>>>>>>>>>>>>>>>>>>>>>"""
    
    def __set_flatcfs(self,flatcfs):
        pmi.invoke(self.pmiobj,'_set_flatcfs_chunk',flatcfs)
    def __get_flatcfs(self):
        flatcfs_list = pmi.invoke(self.pmiobj,'_get_flatcfs_chunk')
        flatcfs = []
        for flatcfs_row in flatcfs_list:
            flatcfs.extend(flatcfs_row)
        return flatcfs
    def __del_flatcfs(self):
        pmi.invoke(self.pmiobj,'_set_flatcfs_chunk',None)
    flatcfs = property(__get_flatcfs,__set_flatcfs,__del_flatcfs)
    
    
    def __set_strength(self,strength):
        pmi.invoke(self.pmiobj,'_set_strength',strength)
    def __get_strength(self):
        strength = pmi.invoke(self.pmiobj,'_get_strength')
        return strength[0]
    def __del_strength(self):
        pmi.invoke(self.pmiobj,'_set_strength',None)
    strength = property(__get_strength,__set_strength,__del_strength)
    

    def __set_activity(self,activity):
        #flattening activity matrix
        self.activity_shape = activity.shape
        #reshaping into one-dimensional matrix (2d matrix that has only one row)
        activity = activity.reshape(1,self.activity_shape[0] * self.activity_shape[1])
        pmi.invoke(self.pmiobj,'_set_activity',list(activity[0]))
    def __get_activity(self):
        activity_list = pmi.invoke(self.pmiobj,'_get_activity')
        activity = []
        for activity_row in activity_list:
            activity.extend(activity_row)
        activity = numpy.array(activity)
        return activity.reshape(self.activity_shape[0],self.activity_shape[1])
    def __del_activity(self):
        pmi.invoke(self.pmiobj,'_set_activity', None)
    activity = property(__get_activity,__set_activity,__del_activity)

This is the serial part of implementation. This means that this code is run as if mpirun command wasn't issued. Members of the super-class (CFProjection) that have to be distributed in order to do computations are implemented as properties. For each property, get and set methods call pmi.invoke which, in turn, calls the MPI method specified as the parameter to mpi.invoke call. If something else is passed as parameter to pmi.invoke, it goes through to parallel method as a parameter, one copy per node. Thus, if you do pmi.invoke(self.pmiobj, "some_method", "abc"), then some_method will be called in the parallel mode and an instance of "abc" string will be passed to it. Obviously, some_method(self, some_string) has to exist somewhere.

Here's the (part of) parallel code that is being run every time pmi.invoke call occurs:

class MPI_CFProjection_node(CFProjection):
    def __init__(self):
        self.comm = MPI.COMM_WORLD
        self.rank = self.comm.Get_rank()
        self.size = self.comm.Get_size()
        print "NODE", MPI.COMM_WORLD.Get_rank(), "initialized"
    
    def _set_flatcfs_chunk(self, flatcfs):
        if flatcfs==None:
            self.flatcfs = None
        else:
            cfs_per_node = int(round(len(flatcfs)/self.size))
            if self.rank+1<self.size:
                self.flatcfs = flatcfs[self.rank * cfs_per_node : (self.rank+1) * cfs_per_node]
            else:
                self.flatcfs = flatcfs[self.rank * cfs_per_node : len(flatcfs)]
    def _get_flatcfs_chunk(self):
        return self.flatcfs

    
    def _set_activity(self, activity):
        if activity==None:
            self.activity = None
        else:
            items_per_node = int(round(len(activity) / self.size))
            if self.rank+1<self.size:
                self.activity = activity[self.rank * items_per_node : (self.rank+1) * items_per_node]
            else:
                self.activity = activity[self.rank * items_per_node : ]
            self.activity = numpy.array([self.activity])
    def _get_activity(self):
        if self.activity==None:
            return None
        else:
            return list(self.activity[0])

That's how it works.

Day 80

Day 80: out of sync


Hi! It's been a really busy week, had to submit 3 assignments. Now some virus has struck me down, I asked for extension for this week's assignment and have at least a couple of days to work on the project. So, the status is: I am really close to my first (real) distributed simulation with pmi. There is one thing that stands between me and this target: tight coupling of one of the components I need to distribute. Hopefully, as soon as I deal with it, I'll be able to topo.sim.run(1) *hwang! rock'n'roll*... ah, whatever.

Day 73

Day 73

Ok, now it's only two...

Day 72

Day 72.

3 deadlines coming up, sorry guys.

Day 68

Day 68.

A quick update. PMI has proven to be very useful, the code now looks a lot better. I've rewritten the MPI_CFProjection so that now it is entirely serial and all the content of this class that has to be parallelized now has properties that link to the MPI_CFProjection_node parallel implementation of CFProjection. This is what I'm doing my best to debug at the moment. I'll post some code as soon as it starts working.

Day 67

Day 67. "Wow, PMI rules!"

Yes, indeed I was wrong. Lets have a look at this file (my edited version of hello_class.py, which comes with the PMI archive):

if __name__ != 'pmi':
    ##################################################
    ## Serial code
    ##################################################
    import pmi

    pmi.setup()
    pmi.execfile_(__file__)

    create a frontend class

    class Hello(object):
        def __init__(self, name):
            self.pmiobj = pmi.create('HelloLocal', name)
        def __call__(self):
            return pmi.invoke(self.pmiobj, '__call__')

        def func(self):
            return pmi.invoke(self.pmiobj,'func')

    # use the class
    hello = Hello('Olaf')
    print('\n'.join(hello()))

    a = hello.func()

    print "Wow, PMI rules!"

else:
    ##################################################
    ## Parallel code
    ##################################################
    from mpi4py import MPI

    class HelloLocal(object):
        def __init__(self, name):
            self.name = name
        def __call__(self):
            return 'Hello %s, this is MPI task %d!' % (self.name, MPI.COMM_WORLD.rank)
        def func(self):
            if MPI.COMM_WORLD.rank == 0:
                return raw_input()
            else:
                return "Slave node" 

Things to pay attention to:

1) First, name object is passed as a parameter from the serial code into the parallel class, and than a separate call to parallel code uses this object. This means that all objects persist between separate calls to parallel code.

2) While calling the func function, only the master node (rank 0) is waiting for user input. However, as opposed to using just mpi, all other nodes are waiting for the master node to complete operation before all nodes join again in the serial code.

3) Notice that "Wow, PMI rules!" is printed out only once, no matter how many nodes you create with mpirun. Now comment out everything else, leave just this print statement, and mpirun this file once again, creating say 4 processes. What will you get?

"Wow, PMI rules!"
"Wow, PMI rules!"
"Wow, PMI rules!"
"Wow, PMI rules!"

Indeed, it does.

Day 66

Day 66

Yesterday, one thought wouldn't let me sleep - I might be wrong about PMI. I think, I'll have another look at it.

Day 66

Day 66

It is a sad, sad day. My favorite forge www.kenai.com has officially announced that it's life is coming to an end.

Day 65

Day 65

I just thought I haven't posted for a while, so this is a quick update on what's currently going on.

First, this Monday my laptop has refused to switch on and now seems to be pretty much dead, but hopefully this is just coma and anyways my warranty hasn't expired yet, so this Wednesday it's going to be picked up and hopefully returned in a week, fully functional. Not a big issue, just a bit annoying.

Second, some time ago Chris suggested I have a look at PMI - "a pure python module that allows libraries to provide functions that are parallelized using MPI but that can nonetheless be called from serial Python scripts." I'll write about it if I find the way to utilize this thing, but at the moment it seems to be quite a handy module. However, here's what I want to get rid of in my implementation by using PMI:

'''
Created on 20 Jan 2010

@author: megatelevizor
'''
from mpi4py import MPI
import time

from MPI_client import client


""" Vocabulary

a   = assign flatcfs
c   = compute activity
i   = Initialize
i_a = Get a copy of input_activity matrix
k   = Kill
p   = Print
r   = Return activity
s   = Set strength


"""

class communicator(object):

   def listen(self):
       self.comm = MPI.COMM_WORLD
       self.rank = self.comm.Get_rank()
      
       self.client_node = client()
      
       while True:
           while not self.comm.Iprobe(source=0, tag=0):
               time.sleep(0.1)
           command_prefix = self.comm.recv(source=0, tag=0)
          
           """ Print """
           if command_prefix == "p":
               i_command_body = self.comm.recv(source=0,tag=1)
               self.client_node.print_command(command = i_command_body)
          
           """ Kill """
           if command_prefix == "k":
               exit()
          
           """ Initialize MPI_CFProjection_slave"""
           if command_prefix == "i":
               self.client_node.init_mpi_cfprojection_slave(node_id=self.rank)
          
           """ Assign flatcfs"""
           if command_prefix == "a":
               new_flatcfs = self.comm.recv(source=0,tag=1)
               print "Node",self.rank,": received", len(new_flatcfs),"flatcfs"
               self.client_node.assign_cfs(new_flatcfs)
          
           """ Set a copy of input activity matrix """
           if command_prefix == "i_a":
               input_activity = self.comm.recv(source=0,tag=1)
               self.client_node.set_input_activity(input_activity)
                      
           """ Compute activity """
           if command_prefix == "c":
               self.client_node.compute_activity()
              
           """ Return activity """
           if command_prefix == "r":
               activity = self.client_node.get_activity()
               print "Node",self.rank,"activity:",activity
               self.comm.send(obj=activity, dest=0, tag=2)
          
           """ Set strength"""
           if command_prefix == "s":
               strength = self.comm.recv(source=0,tag=1)
               self.client_node.set_strength(strength)
              

              
              

""" TAGS:
0 - command prefix
1 - input data
2 - output data
"""

This class allows master node to communicate with slave nodes. Potentially, PMI will eliminate all need to have such class.



[update] No, it will not. What PMI allows you to do is hide the parallel implementation, so that your programme looks like it is a serial (single-processed) code. However, it does not enhance the underlying MPI implementation in any way, so it still needs to be implicit. What does that mean in my case? Well, if I was to introduce PMI into my code as it is now, then there would have been a class - something like CFProjection_proxy, linked to the real (i.e. parallel) implementation of CFProjection_slave using PMI mechanism. CFPRojection_slave's would still have to communicate with the master node CFProjection using MPI_io (the code listed above) because PMI does not influence or change the way different nodes of an MPI application communicate. In simple words, it just makes your code look like it's not using MPI, when it actually does.

Well, so what's so bad about it? Why not use PMI anyway, since from all that it looks like it's actually doing a good job? Actually, there's nothing wrong with it and it might be a good idea to use it. No, it won't add any new functionality, but it will definitelly make things look a lot nicer. However, I prefer dealing with things one at a time, and also finishing what I started, so I think I'll concentrate on looking for the solution to my initial problem.

To summarise what I'm looking for, have a look at this: Java RMI . Without going into too much detail, RMI is a mechanism that allows Java objects calling each other's methods accross several different JVM's (which can be run on different hosts!) and it stands for Remote Method Invocation. Combined with serialization (something similar to Python's pickling), it allows to hide all explicit communication-level code away from the user, which is great! Thus, my assumption is that if something exists in Java, it has a chance to exist in Python. As simple as that.

Searching goes on.