1.16. Weave (scipy.weave)¶

1.16.7.2.1. Binary search¶

Lets look at the example of searching a sorted list of integers for a value. For inspiration, we’ll use Kalle Svensson’s binary_search() algorithm from the Python Cookbook. His recipe follows:

def binary_search(seq, t):
    min = 0; max = len(seq) - 1
    while 1:
        if max < min:
            return -1
        m = (min  + max)  / 2
        if seq[m] < t:
            min = m  + 1
        elif seq[m] > t:
            max = m  - 1
        else:
            return m

This Python version works for arbitrary Python data types. The C version below is specialized to handle integer values. There is a little type checking done in Python to assure that we’re working with the correct data types before heading into C. The variables seq and t don’t need to be declared because weave handles converting and declaring them in the C code. All other temporary variables such as min, max, etc. must be declared – it is C after all. Here’s the new mixed Python/C function:

def c_int_binary_search(seq,t):
    # do a little type checking in Python
    assert(type(t) == type(1))
    assert(type(seq) == type([]))

    # now the C code
    code = """
           #line 29 "binary_search.py"
           int val, m, min = 0;
           int max = seq.length() - 1;
           PyObject *py_val;
           for(;;)
           {
               if (max < min  )
               {
                   return_val =  Py::new_reference_to(Py::Int(-1));
                   break;
               }
               m =  (min + max) /2;
               val = py_to_int(PyList_GetItem(seq.ptr(),m),"val");
               if (val  < t)
                   min = m  + 1;
               else if (val >  t)
                   max = m - 1;
               else
               {
                   return_val = Py::new_reference_to(Py::Int(m));
                   break;
               }
           }
           """
    return inline(code,['seq','t'])

We have two variables seq and t passed in. t is guaranteed (by the assert) to be an integer. Python integers are converted to C int types in the transition from Python to C. seq is a Python list. By default, it is translated to a CXX list object. Full documentation for the CXX library can be found at its website. The basics are that the CXX provides C++ class equivalents for Python objects that simplify, or at least object orientify, working with Python objects in C/C++. For example, seq.length() returns the length of the list. A little more about CXX and its class methods, etc. is in the Type Conversions section.

Most of the algorithm above looks similar in C to the original Python code. There are two main differences. The first is the setting of return_val instead of directly returning from the C code with a return statement. return_val is an automatically defined variable of type PyObject* that is returned from the C code back to Python. You’ll have to handle reference counting issues when setting this variable. In this example, CXX classes and functions handle the dirty work. All CXX functions and classes live in the namespace Py::. The following code converts the integer m to a CXX Int() object and then to a PyObject* with an incremented reference count using Py::new_reference_to().

return_val = Py::new_reference_to(Py::Int(m));

The second big differences shows up in the retrieval of integer values from the Python list. The simple Python seq[i] call balloons into a C Python API call to grab the value out of the list and then a separate call to py_to_int() that converts the PyObject* to an integer. py_to_int() includes both a NULL check and a PyInt_Check() call as well as the conversion call. If either of the checks fail, an exception is raised. The entire C++ code block is executed with in a try/catch block that handles exceptions much like Python does. This removes the need for most error checking code.

It is worth note that CXX lists do have indexing operators that result in code that looks much like Python. However, the overhead in using them appears to be relatively high, so the standard Python API was used on the seq.ptr() which is the underlying PyObject* of the List object.

The #line directive that is the first line of the C code block isn’t necessary, but it’s nice for debugging. If the compilation fails because of the syntax error in the code, the error will be reported as an error in the Python file “binary_search.py” with an offset from the given line number (29 here).

So what was all our effort worth in terms of efficiency? Well not a lot in this case. The examples/binary_search.py file runs both Python and C versions of the functions As well as using the standard bisect module. If we run it on a 1 million element list and run the search 3000 times (for 0- 2999), here are the results we get:

C:\home\ej\wrk\scipy\weave\examples> python binary_search.py
Binary search for 3000 items in 1000000 length list of integers:
speed in python: 0.159999966621
speed of bisect: 0.121000051498
speed up: 1.32
speed in c: 0.110000014305
speed up: 1.45
speed in c(no asserts): 0.0900000333786
speed up: 1.78

So, we get roughly a 50-75% improvement depending on whether we use the Python asserts in our C version. If we move down to searching a 10000 element list, the advantage evaporates. Even smaller lists might result in the Python version being faster. I’d like to say that moving to NumPy lists (and getting rid of the GetItem() call) offers a substantial speed up, but my preliminary efforts didn’t produce one. I think the log(N) algorithm is to blame. Because the algorithm is nice, there just isn’t much time spent computing things, so moving to C isn’t that big of a win. If there are ways to reduce conversion overhead of values, this may improve the C/Python speed up. Anyone have other explanations or faster code, please let me know.

1.16.7.2.2. Dictionary Sort¶

The demo in examples/dict_sort.py is another example from the Python CookBook. This submission, by Alex Martelli, demonstrates how to return the values from a dictionary sorted by their keys:

def sortedDictValues3(adict):
    keys = adict.keys()
    keys.sort()
    return map(adict.get, keys)

Alex provides 3 algorithms and this is the 3rd and fastest of the set. The C version of this same algorithm follows:

def c_sort(adict):
    assert(type(adict) == type({}))
    code = """
    #line 21 "dict_sort.py"
    Py::List keys = adict.keys();
    Py::List items(keys.length()); keys.sort();
    PyObject* item = NULL;
    for(int i = 0;  i < keys.length();i++)
    {
        item = PyList_GET_ITEM(keys.ptr(),i);
        item = PyDict_GetItem(adict.ptr(),item);
        Py_XINCREF(item);
        PyList_SetItem(items.ptr(),i,item);
    }
    return_val = Py::new_reference_to(items);
    """
    return inline_tools.inline(code,['adict'],verbose=1)

Like the original Python function, the C++ version can handle any Python dictionary regardless of the key/value pair types. It uses CXX objects for the most part to declare python types in C++, but uses Python API calls to manipulate their contents. Again, this choice is made for speed. The C++ version, while more complicated, is about a factor of 2 faster than Python.

C:\home\ej\wrk\scipy\weave\examples> python dict_sort.py
Dict sort of 1000 items for 300 iterations:
 speed in python: 0.319999933243
[0, 1, 2, 3, 4]
 speed in c: 0.151000022888
 speed up: 2.12
[0, 1, 2, 3, 4]

1.16.7.2.3. NumPy – cast/copy/transpose¶

CastCopyTranspose is a function called quite heavily by Linear Algebra routines in the NumPy library. Its needed in part because of the row-major memory layout of multi-dimensional Python (and C) arrays vs. the col-major order of the underlying Fortran algorithms. For small matrices (say 100x100 or less), a significant portion of the common routines such as LU decomposition or singular value decomposition are spent in this setup routine. This shouldn’t happen. Here is the Python version of the function using standard NumPy operations.

def _castCopyAndTranspose(type, array):
    if a.typecode() == type:
        cast_array = copy.copy(NumPy.transpose(a))
    else:
        cast_array = copy.copy(NumPy.transpose(a).astype(type))
    return cast_array

And the following is a inline C version of the same function:

from weave.blitz_tools import blitz_type_factories
from weave import scalar_spec
from weave import inline
def _cast_copy_transpose(type,a_2d):
    assert(len(shape(a_2d)) == 2)
    new_array = zeros(shape(a_2d),type)
    NumPy_type = scalar_spec.NumPy_to_blitz_type_mapping[type]
    code = \
    """
    for(int i = 0;i < _Na_2d[0]; i++)
        for(int j = 0;  j < _Na_2d[1]; j++)
            new_array(i,j) = (%s) a_2d(j,i);
    """ % NumPy_type
    inline(code,['new_array','a_2d'],
           type_factories = blitz_type_factories,compiler='gcc')
    return new_array

This example uses blitz++ arrays instead of the standard representation of NumPy arrays so that indexing is simpler to write. This is accomplished by passing in the blitz++ “type factories” to override the standard Python to C++ type conversions. Blitz++ arrays allow you to write clean, fast code, but they also are sloooow to compile (20 seconds or more for this snippet). This is why they aren’t the default type used for Numeric arrays (and also because most compilers can’t compile blitz arrays...). inline() is also forced to use ‘gcc’ as the compiler because the default compiler on Windows (MSVC) will not compile blitz code. (‘gcc’ I think will use the standard compiler on Unix machine instead of explicitly forcing gcc (check this)) Comparisons of the Python vs inline C++ code show a factor of 3 speed up. Also shown are the results of an “inplace” transpose routine that can be used if the output of the linear algebra routine can overwrite the original matrix (this is often appropriate). This provides another factor of 2 improvement.

#C:\home\ej\wrk\scipy\weave\examples> python cast_copy_transpose.py
# Cast/Copy/Transposing (150,150)array 1 times
#  speed in python: 0.870999932289
#  speed in c: 0.25
#  speed up: 3.48
#  inplace transpose c: 0.129999995232
#  speed up: 6.70

1.16.7.2.4. wxPython¶

inline knows how to handle wxPython objects. That’s nice in and of itself, but it also demonstrates that the type conversion mechanism is reasonably flexible. Chances are, it won’t take a ton of effort to support special types you might have. The examples/wx_example.py borrows the scrolled window example from the wxPython demo, accept that it mixes inline C code in the middle of the drawing function.

def DoDrawing(self, dc):

    red = wxNamedColour("RED");
    blue = wxNamedColour("BLUE");
    grey_brush = wxLIGHT_GREY_BRUSH;
    code = \
    """
    #line 108 "wx_example.py"
    dc->BeginDrawing();
    dc->SetPen(wxPen(*red,4,wxSOLID));
    dc->DrawRectangle(5,5,50,50);
    dc->SetBrush(*grey_brush);
    dc->SetPen(wxPen(*blue,4,wxSOLID));
    dc->DrawRectangle(15, 15, 50, 50);
    """
    inline(code,['dc','red','blue','grey_brush'])

    dc.SetFont(wxFont(14, wxSWISS, wxNORMAL, wxNORMAL))
    dc.SetTextForeground(wxColour(0xFF, 0x20, 0xFF))
    te = dc.GetTextExtent("Hello World")
    dc.DrawText("Hello World", 60, 65)

    dc.SetPen(wxPen(wxNamedColour('VIOLET'), 4))
    dc.DrawLine(5, 65+te[1], 60+te[0], 65+te[1])
    ...

Here, some of the Python calls to wx objects were just converted to C++ calls. There isn’t any benefit, it just demonstrates the capabilities. You might want to use this if you have a computationally intensive loop in your drawing code that you want to speed up. On windows, you’ll have to use the MSVC compiler if you use the standard wxPython DLLs distributed by Robin Dunn. That’s because MSVC and gcc, while binary compatible in C, are not binary compatible for C++. In fact, its probably best, no matter what platform you’re on, to specify that inline use the same compiler that was used to build wxPython to be on the safe side. There isn’t currently a way to learn this info from the library – you just have to know. Also, at least on the windows platform, you’ll need to install the wxWindows libraries and link to them. I think there is a way around this, but I haven’t found it yet – I get some linking errors dealing with wxString. One final note. You’ll probably have to tweak weave/wx_spec.py or weave/wx_info.py for your machine’s configuration to point at the correct directories etc. There. That should sufficiently scare people into not even looking at this... :)

1.16.7.5.1. Keyword Option Examples¶

We’ll walk through several examples here to demonstrate the behavior of inline and also how the various arguments are used. In the simplest (most) cases, code and arg_names are the only arguments that need to be specified. Here’s a simple example run on Windows machine that has Microsoft VC++ installed.

>>> from weave import inline
>>> a = 'string'
>>> code = """
...        int l = a.length();
...        return_val = Py::new_reference_to(Py::Int(l));
...        """
>>> inline(code,['a'])
 sc_86e98826b65b047ffd2cd5f479c627f12.cpp
Creating
   library C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\sc_86e98826b65b047ffd2cd5f479c627f12.lib
and object C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\sc_86e98826b65b047ff
d2cd5f479c627f12.exp
6
>>> inline(code,['a'])
6

When inline is first run, you’ll notice that pause and some trash printed to the screen. The “trash” is actually part of the compiler’s output that distutils does not suppress. The name of the extension file, sc_bighonkingnumber.cpp, is generated from the SHA-256 check sum of the C/C++ code fragment. On Unix or windows machines with only gcc installed, the trash will not appear. On the second call, the code fragment is not compiled since it already exists, and only the answer is returned. Now kill the interpreter and restart, and run the same code with a different string.

>>> from weave import inline
>>> a = 'a longer string'
>>> code = """
...        int l = a.length();
...        return_val = Py::new_reference_to(Py::Int(l));
...        """
>>> inline(code,['a'])
15

Notice this time, inline() did not recompile the code because it found the compiled function in the persistent catalog of functions. There is a short pause as it looks up and loads the function, but it is much shorter than compiling would require.

You can specify the local and global dictionaries if you’d like (much like exec or eval() in Python), but if they aren’t specified, the “expected” ones are used – i.e. the ones from the function that called inline(). This is accomplished through a little call frame trickery. Here is an example where the local_dict is specified using the same code example from above:

>>> a = 'a longer string'
>>> b = 'an even  longer string'
>>> my_dict = {'a':b}
>>> inline(code,['a'])
15
>>> inline(code,['a'],my_dict)
21

Every time the code is changed, inline does a recompile. However, changing any of the other options in inline does not force a recompile. The force option was added so that one could force a recompile when tinkering with other variables. In practice, it is just as easy to change the code by a single character (like adding a space some place) to force the recompile.

I use verbose sometimes for debugging. When set to 2, it’ll output all the information (including the name of the .cpp file) that you’d expect from running a make file. This is nice if you need to examine the generated code to see where things are going haywire. Note that error messages from failed compiles are printed to the screen even if verbose is set to 0.

The following example demonstrates using gcc instead of the standard msvc compiler on windows using same code fragment as above. Because the example has already been compiled, the force=1 flag is needed to make inline() ignore the previously compiled version and recompile using gcc. The verbose flag is added to show what is printed out:

>>>inline(code,['a'],compiler='gcc',verbose=2,force=1)
running build_ext
building 'sc_86e98826b65b047ffd2cd5f479c627f13' extension
c:\gcc-2.95.2\bin\g++.exe -mno-cygwin -mdll -O2 -w -Wstrict-prototypes -IC:
\home\ej\wrk\scipy\weave -IC:\Python21\Include -c C:\DOCUME~1\eric\LOCAL
S~1\Temp\python21_compiled\sc_86e98826b65b047ffd2cd5f479c627f13.cpp
-o C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\sc_86e98826b65b04ffd2cd5f479c627f13.o
skipping C:\home\ej\wrk\scipy\weave\CXX\cxxextensions.c
(C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\cxxextensions.o up-to-date)
skipping C:\home\ej\wrk\scipy\weave\CXX\cxxsupport.cxx
(C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\cxxsupport.o up-to-date)
skipping C:\home\ej\wrk\scipy\weave\CXX\IndirectPythonInterface.cxx
(C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\indirectpythoninterface.o up-to-date)
skipping C:\home\ej\wrk\scipy\weave\CXX\cxx_extensions.cxx
(C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\cxx_extensions.o
up-to-date)
writing C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\sc_86e98826b65b047ffd2cd5f479c627f13.def
c:\gcc-2.95.2\bin\dllwrap.exe --driver-name g++ -mno-cygwin
-mdll -static --output-lib
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\libsc_86e98826b65b047ffd2cd5f479c627f13.a --def
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\sc_86e98826b65b047ffd2cd5f479c627f13.def
-sC:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\sc_86e98826b65b047ffd2cd5f479c627f13.o
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\cxxextensions.o
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\cxxsupport.o
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\indirectpythoninterface.o
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\temp\Release\cxx_extensions.o -LC:\Python21\libs
-lpython21 -o
C:\DOCUME~1\eric\LOCALS~1\Temp\python21_compiled\sc_86e98826b65b047ffd2cd5f479c627f13.pyd
15

That’s quite a bit of output. verbose=1 just prints the compile time.

>>>inline(code,['a'],compiler='gcc',verbose=1,force=1)
Compiling code...
finished compiling (sec):  6.00800001621
15

The support_code argument is likely to be used a lot. It allows you to specify extra code fragments such as function, structure or class definitions that you want to use in the code string. Note that changes to support_code do not force a recompile. The catalog only relies on code (for performance reasons) to determine whether recompiling is necessary. So, if you make a change to support_code, you’ll need to alter code in some way or use the force argument to get the code to recompile. I usually just add some innocuous whitespace to the end of one of the lines in code somewhere. Here’s an example of defining a separate method for calculating the string length:

>>> from weave import inline
>>> a = 'a longer string'
>>> support_code = """
...                PyObject* length(Py::String a)
...                {
...                    int l = a.length();
...                    return Py::new_reference_to(Py::Int(l));
...                }
...                """
>>> inline("return_val = length(a);",['a'],
...        support_code = support_code)
15

customize is a left over from a previous way of specifying compiler options. It is a custom_info object that can specify quite a bit of information about how a file is compiled. These info objects are the standard way of defining compile information for type conversion classes. However, I don’t think they are as handy here, especially since we’ve exposed all the keyword arguments that distutils can handle. Between these keywords, and the support_code option, I think customize may be obsolete. We’ll see if anyone cares to use it. If not, it’ll get axed in the next version.

The type_factories variable is important to people who want to customize the way arguments are converted from Python to C. We’ll talk about this in the next chapter xx of this document when we discuss type conversions.

auto_downcast handles one of the big type conversion issues that is common when using NumPy arrays in conjunction with Python scalar values. If you have an array of single precision values and multiply that array by a Python scalar, the result is upcast to a double precision array because the scalar value is double precision. This is not usually the desired behavior because it can double your memory usage. auto_downcast goes some distance towards changing the casting precedence of arrays and scalars. If your only using single precision arrays, it will automatically downcast all scalar values from double to single precision when they are passed into the C++ code. This is the default behavior. If you want all values to keep there default type, set auto_downcast to 0.

1.16.7.5.2. Returning Values¶

Python variables in the local and global scope transfer seamlessly from Python into the C++ snippets. And, if inline were to completely live up to its name, any modifications to variables in the C++ code would be reflected in the Python variables when control was passed back to Python. For example, the desired behavior would be something like:

# THIS DOES NOT WORK
>>> a = 1
>>> weave.inline("a++;",['a'])
>>> a
2

Instead you get:

>>> a = 1
>>> weave.inline("a++;",['a'])
>>> a
1

Variables are passed into C++ as if you are calling a Python function. Python’s calling convention is sometimes called “pass by assignment”. This means its as if a c_a = a assignment is made right before inline call is made and the c_a variable is used within the C++ code. Thus, any changes made to c_a are not reflected in Python’s a variable. Things do get a little more confusing, however, when looking at variables with mutable types. Changes made in C++ to the contents of mutable types are reflected in the Python variables.

>>> a= [1,2]
>>> weave.inline("PyList_SetItem(a.ptr(),0,PyInt_FromLong(3));",['a'])
>>> print a
[3, 2]

So modifications to the contents of mutable types in C++ are seen when control is returned to Python. Modifications to immutable types such as tuples, strings, and numbers do not alter the Python variables. If you need to make changes to an immutable variable, you’ll need to assign the new value to the “magic” variable return_val in C++. This value is returned by the inline() function:

>>> a = 1
>>> a = weave.inline("return_val = Py::new_reference_to(Py::Int(a+1));",['a'])
>>> a
2

The return_val variable can also be used to return newly created values. This is possible by returning a tuple. The following trivial example illustrates how this can be done:

# python version
def multi_return():
    return 1, '2nd'

# C version.
def c_multi_return():
    code =  """
              py::tuple results(2);
              results[0] = 1;
              results[1] = "2nd";
              return_val = results;
            """
    return inline_tools.inline(code)

The example is available in examples/tuple_return.py. It also has the dubious honor of demonstrating how much inline() can slow things down. The C version here is about 7-10 times slower than the Python version. Of course, something so trivial has no reason to be written in C anyway.

1.16.7.5.3. A quick look at the code¶

weave generates a C++ file holding an extension function for each inline code snippet. These file names are generated using from the SHA-256 signature of the code snippet and saved to a location specified by the PYTHONCOMPILED environment variable (discussed later). The cpp files are generally about 200-400 lines long and include quite a few functions to support type conversions, etc. However, the actual compiled function is pretty simple. Below is the familiar printf example:

>>> import weave
>>> a = 1
>>> weave.inline('printf("%d\\n",a);',['a'])
1

And here is the extension function generated by inline:

static PyObject* compiled_func(PyObject*self, PyObject* args)
{
    py::object return_val;
    int exception_occurred = 0;
    PyObject *py__locals = NULL;
    PyObject *py__globals = NULL;
    PyObject *py_a;
    py_a = NULL;

    if(!PyArg_ParseTuple(args,"OO:compiled_func",&py__locals,&py__globals))
        return NULL;
    try
    {
        PyObject* raw_locals = py_to_raw_dict(py__locals,"_locals");
        PyObject* raw_globals = py_to_raw_dict(py__globals,"_globals");
        /* argument conversion code */
        py_a = get_variable("a",raw_locals,raw_globals);
        int a = convert_to_int(py_a,"a");
        /* inline code */
        /* NDARRAY API VERSION 90907 */
        printf("%d\n",a);    /*I would like to fill in changed locals and globals here...*/
    }
    catch(...)
    {
        return_val =  py::object();
        exception_occurred = 1;
    }
    /* cleanup code */
    if(!(PyObject*)return_val && !exception_occurred)
    {
        return_val = Py_None;
    }
    return return_val.disown();
}

Every inline function takes exactly two arguments – the local and global dictionaries for the current scope. All variable values are looked up out of these dictionaries. The lookups, along with all inline code execution, are done within a C++ try block. If the variables aren’t found, or there is an error converting a Python variable to the appropriate type in C++, an exception is raised. The C++ exception is automatically converted to a Python exception by SCXX and returned to Python. The py_to_int() function illustrates how the conversions and exception handling works. py_to_int first checks that the given PyObject* pointer is not NULL and is a Python integer. If all is well, it calls the Python API to convert the value to an int. Otherwise, it calls handle_bad_type() which gathers information about what went wrong and then raises a SCXX TypeError which returns to Python as a TypeError.

int py_to_int(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyInt_Check(py_obj))
        handle_bad_type(py_obj,"int", name);
    return (int) PyInt_AsLong(py_obj);
}

void handle_bad_type(PyObject* py_obj, char* good_type, char* var_name)
{
    char msg[500];
    sprintf(msg,"received '%s' type instead of '%s' for variable '%s'",
            find_type(py_obj),good_type,var_name);
    throw Py::TypeError(msg);
}

char* find_type(PyObject* py_obj)
{
    if(py_obj == NULL) return "C NULL value";
    if(PyCallable_Check(py_obj)) return "callable";
    if(PyString_Check(py_obj)) return "string";
    if(PyInt_Check(py_obj)) return "int";
    if(PyFloat_Check(py_obj)) return "float";
    if(PyDict_Check(py_obj)) return "dict";
    if(PyList_Check(py_obj)) return "list";
    if(PyTuple_Check(py_obj)) return "tuple";
    if(PyFile_Check(py_obj)) return "file";
    if(PyModule_Check(py_obj)) return "module";

    //should probably do more interrogation (and thinking) on these.
    if(PyCallable_Check(py_obj) && PyInstance_Check(py_obj)) return "callable";
    if(PyInstance_Check(py_obj)) return "instance";
    if(PyCallable_Check(py_obj)) return "callable";
    return "unknown type";
}

Since the inline is also executed within the try/catch block, you can use CXX exceptions within your code. It is usually a bad idea to directly return from your code, even if an error occurs. This skips the clean up section of the extension function. In this simple example, there isn’t any clean up code, but in more complicated examples, there may be some reference counting that needs to be taken care of here on converted variables. To avoid this, either uses exceptions or set return_val to NULL and use if/then's to skip code after errors.

1.16.7.8.1. NumPy Argument Conversion¶

Integer, floating point, and complex arguments are handled in a very similar fashion. Consider the following inline function that has a single integer variable passed in:

>>> a = 1
>>> inline("",['a'])

The argument conversion code inserted for a is:

/* argument conversion code */
int a = py_to_int (get_variable("a",raw_locals,raw_globals),"a");

get_variable() reads the variable a from the local and global namespaces. py_to_int() has the following form:

static int py_to_int(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyInt_Check(py_obj))
        handle_bad_type(py_obj,"int", name);
    return (int) PyInt_AsLong(py_obj);
}

Similarly, the float and complex conversion routines look like:

static double py_to_float(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyFloat_Check(py_obj))
        handle_bad_type(py_obj,"float", name);
    return PyFloat_AsDouble(py_obj);
}

static std::complex py_to_complex(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyComplex_Check(py_obj))
        handle_bad_type(py_obj,"complex", name);
    return std::complex(PyComplex_RealAsDouble(py_obj),
                                PyComplex_ImagAsDouble(py_obj));
}

NumPy conversions do not require any clean up code.

1.16.7.8.2. String, List, Tuple, and Dictionary Conversion¶

Strings, Lists, Tuples and Dictionary conversions are all converted to SCXX types by default. For the following code,

>>> a = [1]
>>> inline("",['a'])

The argument conversion code inserted for a is:

/* argument conversion code */
Py::List a = py_to_list(get_variable("a",raw_locals,raw_globals),"a");

get_variable() reads the variable a from the local and global namespaces. py_to_list() and its friends have the following form:

static Py::List py_to_list(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyList_Check(py_obj))
        handle_bad_type(py_obj,"list", name);
    return Py::List(py_obj);
}

static Py::String py_to_string(PyObject* py_obj,char* name)
{
    if (!PyString_Check(py_obj))
        handle_bad_type(py_obj,"string", name);
    return Py::String(py_obj);
}

static Py::Dict py_to_dict(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyDict_Check(py_obj))
        handle_bad_type(py_obj,"dict", name);
    return Py::Dict(py_obj);
}

static Py::Tuple py_to_tuple(PyObject* py_obj,char* name)
{
    if (!py_obj || !PyTuple_Check(py_obj))
        handle_bad_type(py_obj,"tuple", name);
    return Py::Tuple(py_obj);
}

SCXX handles reference counts on for strings, lists, tuples, and dictionaries, so clean up code isn’t necessary.

1.16.7.8.3. File Conversion¶

For the following code,

>>> a = open("bob",'w')
>>> inline("",['a'])

The argument conversion code is:

/* argument conversion code */
PyObject* py_a = get_variable("a",raw_locals,raw_globals);
FILE* a = py_to_file(py_a,"a");

get_variable() reads the variable a from the local and global namespaces. py_to_file() converts PyObject* to a FILE* and increments the reference count of the PyObject*:

FILE* py_to_file(PyObject* py_obj, char* name)
{
    if (!py_obj || !PyFile_Check(py_obj))
        handle_bad_type(py_obj,"file", name);

    Py_INCREF(py_obj);
    return PyFile_AsFile(py_obj);
}

Because the PyObject* was incremented, the clean up code needs to decrement the counter

/* cleanup code */
Py_XDECREF(py_a);

Its important to understand that file conversion only works on actual files – i.e. ones created using the open() command in Python. It does not support converting arbitrary objects that support the file interface into C FILE* pointers. This can affect many things. For example, in initial printf() examples, one might be tempted to solve the problem of C and Python IDE’s (PythonWin, PyCrust, etc.) writing to different stdout and stderr by using fprintf() and passing in sys.stdout and sys.stderr. For example, instead of

>>> weave.inline('printf("hello\\n");')

You might try:

>>> buf = sys.stdout
>>> weave.inline('fprintf(buf,"hello\\n");',['buf'])

This will work as expected from a standard python interpreter, but in PythonWin, the following occurs:

>>> buf = sys.stdout
>>> weave.inline('fprintf(buf,"hello\\n");',['buf'])
Traceback (most recent call last):
    File "", line 1, in ?
    File "C:\Python21\weave\inline_tools.py", line 315, in inline
        auto_downcast = auto_downcast,
    File "C:\Python21\weave\inline_tools.py", line 386, in compile_function
        type_factories = type_factories)
    File "C:\Python21\weave\ext_tools.py", line 197, in __init__
        auto_downcast, type_factories)
    File "C:\Python21\weave\ext_tools.py", line 390, in assign_variable_types
        raise TypeError, format_error_msg(errors)
    TypeError: {'buf': "Unable to convert variable 'buf' to a C++ type."}

The traceback tells us that inline() was unable to convert ‘buf’ to a C++ type (If instance conversion was implemented, the error would have occurred at runtime instead). Why is this? Let’s look at what the buf object really is:

>>> buf
pywin.framework.interact.InteractiveView instance at 00EAD014

PythonWin has reassigned sys.stdout to a special object that implements the Python file interface. This works great in Python, but since the special object doesn’t have a FILE* pointer underlying it, fprintf doesn’t know what to do with it (well this will be the problem when instance conversion is implemented...).

1.16.7.8.4. Callable, Instance, and Module Conversion¶

>>> def a():
    pass
>>> inline("",['a'])

Callable and instance variables are converted to PyObject*. Nothing is done to their reference counts.

/* argument conversion code */
PyObject* a = py_to_callable(get_variable("a",raw_locals,raw_globals),"a");

get_variable() reads the variable a from the local and global namespaces. The py_to_callable() and py_to_instance() don’t currently increment the ref count.

PyObject* py_to_callable(PyObject* py_obj, char* name)
{
    if (!py_obj || !PyCallable_Check(py_obj))
        handle_bad_type(py_obj,"callable", name);
    return py_obj;
}

PyObject* py_to_instance(PyObject* py_obj, char* name)
{
    if (!py_obj || !PyFile_Check(py_obj))
        handle_bad_type(py_obj,"instance", name);
    return py_obj;
}

There is no cleanup code for callables, modules, or instances.

1.16.7.8.5. Customizing Conversions¶

Converting from Python to C++ types is handled by xxx_specification classes. A type specification class actually serve in two related but different roles. The first is in determining whether a Python variable that needs to be converted should be represented by the given class. The second is as a code generator that generates C++ code needed to convert from Python to C++ types for a specific variable.

When

>>> a = 1
>>> weave.inline('printf("%d",a);',['a'])

is called for the first time, the code snippet has to be compiled. In this process, the variable ‘a’ is tested against a list of type specifications (the default list is stored in weave/ext_tools.py). The first specification in the list is used to represent the variable.

Examples of xxx_specification are scattered throughout numerous “xxx_spec.py” files in the weave package. Closely related to the xxx_specification classes are yyy_info classes. These classes contain compiler, header, and support code information necessary for including a certain set of capabilities (such as blitz++ or CXX support) in a compiled module. xxx_specification classes have one or more yyy_info classes associated with them. If you’d like to define your own set of type specifications, the current best route is to examine some of the existing spec and info files. Maybe looking over sequence_spec.py and cxx_info.py are a good place to start. After defining specification classes, you’ll need to pass them into inline using the type_factories argument. A lot of times you may just want to change how a specific variable type is represented. Say you’d rather have Python strings converted to std::string or maybe char* instead of using the CXX string object, but would like all other type conversions to have default behavior. This requires that a new specification class that handles strings is written and then prepended to a list of the default type specifications. Since it is closer to the front of the list, it effectively overrides the default string specification. The following code demonstrates how this is done: ...

1.16.7.9.1. Function Storage¶

Function caches are stored as dictionaries where the key is the entire C++ code string and the value is either a single function (as in the “level 1” cache) or a list of functions (as in the main catalog cache). On disk catalogs are stored in the same manor using standard Python shelves.

Early on, there was a question as to whether md5 checksums of the C++ code strings should be used instead of the actual code strings. I think this is the route inline Perl took. Some (admittedly quick) tests of the md5 vs. the entire string showed that using the entire string was at least a factor of 3 or 4 faster for Python. I think this is because it is more time consuming to compute the md5 value than it is to do look-ups of long strings in the dictionary. Look at the examples/md5_speed.py file for the test run.

1.16.7.9.2. Catalog search paths and the PYTHONCOMPILED variable¶

The default location for catalog files on Unix is ~/.pythonXX_compiled where XX is version of Python being used. If this directory doesn’t exist, it is created the first time a catalog is used. The directory must be writable. If, for any reason it isn’t, then the catalog attempts to create a directory based on your user id in the /tmp directory. The directory permissions are set so that only you have access to the directory. If this fails, I think you’re out of luck. I don’t think either of these should ever fail though. On Windows, a directory called pythonXX_compiled is created in the user’s temporary directory.

The actual catalog file that lives in this directory is a Python shelf with a platform specific name such as “nt21compiled_catalog” so that multiple OSes can share the same file systems without trampling on each other. Along with the catalog file, the .cpp and .so or .pyd files created by inline will live in this directory. The catalog file simply contains keys which are the C++ code strings with values that are lists of functions. The function lists point at functions within these compiled modules. Each function in the lists executes the same C++ code string, but compiled for different input variables.

You can use the PYTHONCOMPILED environment variable to specify alternative locations for compiled functions. On Unix this is a colon (‘:’) separated list of directories. On windows, it is a (‘;’) separated list of directories. These directories will be searched prior to the default directory for a compiled function catalog. Also, the first writable directory in the list is where all new compiled function catalogs, .cpp and .so or .pyd files are written. Relative directory paths (‘.’ and ‘..’) should work fine in the PYTHONCOMPILED variable as should environment variables.

There is a “special” path variable called MODULE that can be placed in the PYTHONCOMPILED variable. It specifies that the compiled catalog should reside in the same directory as the module that called it. This is useful if an admin wants to build a lot of compiled functions during the build of a package and then install them in site-packages along with the package. User’s who specify MODULE in their PYTHONCOMPILED variable will have access to these compiled functions. Note, however, that if they call the function with a set of argument types that it hasn’t previously been built for, the new function will be stored in their default directory (or some other writable directory in the PYTHONCOMPILED path) because the user will not have write access to the site-packages directory.

An example of using the PYTHONCOMPILED path on bash follows:

PYTHONCOMPILED=MODULE:/some/path;export PYTHONCOMPILED;

If you are using python21 on linux, and the module bob.py in site-packages has a compiled function in it, then the catalog search order when calling that function for the first time in a python session would be:

/usr/lib/python21/site-packages/linuxpython_compiled
/some/path/linuxpython_compiled
~/.python21_compiled/linuxpython_compiled

The default location is always included in the search path.

The in-module cache (in weave.inline_tools reduces the overhead of calling inline functions by about a factor of 2. It can be reduced a little more for type loop calls where the same function is called over and over again if the cache was a single value instead of a dictionary, but the benefit is very small (less than 5%) and the utility is quite a bit less. So, we’ll stick with a dictionary as the cache.

1.16.8.4.1. Parser¶

Tearing NumPy expressions apart, examining the pieces, and then rebuilding them as C++ (blitz) expressions requires a parser of some sort. I can imagine someone attacking this problem with regular expressions, but it’d likely be ugly and fragile. Amazingly, Python solves this problem for us. It actually exposes its parsing engine to the world through the parser module. The following fragment creates an Abstract Syntax Tree (AST) object for the expression and then converts to a (rather unpleasant looking) deeply nested list representation of the tree.

>>> import parser
>>> import scipy.weave.misc
>>> ast = parser.suite("a = b * c + d")
>>> ast_list = ast.tolist()
>>> sym_list = scipy.weave.misc.translate_symbols(ast_list)
>>> pprint.pprint(sym_list)
['file_input',
 ['stmt',
  ['simple_stmt',
   ['small_stmt',
    ['expr_stmt',
     ['testlist',
      ['test',
       ['and_test',
        ['not_test',
         ['comparison',
          ['expr',
           ['xor_expr',
            ['and_expr',
             ['shift_expr',
              ['arith_expr',
               ['term',
                ['factor', ['power', ['atom', ['NAME', 'a']]]]]]]]]]]]]]],
     ['EQUAL', '='],
     ['testlist',
      ['test',
       ['and_test',
        ['not_test',
         ['comparison',
          ['expr',
           ['xor_expr',
            ['and_expr',
             ['shift_expr',
              ['arith_expr',
               ['term',
                ['factor', ['power', ['atom', ['NAME', 'b']]]],
                ['STAR', '*'],
                ['factor', ['power', ['atom', ['NAME', 'c']]]]],
               ['PLUS', '+'],
               ['term',
                ['factor', ['power', ['atom', ['NAME', 'd']]]]]]]]]]]]]]]]],
   ['NEWLINE', '']]],
 ['ENDMARKER', '']]

Despite its looks, with some tools developed by Jeremy H., it’s possible to search these trees for specific patterns (sub-trees), extract the sub-tree, manipulate them converting python specific code fragments to blitz code fragments, and then re-insert it in the parse tree. The parser module documentation has some details on how to do this. Traversing the new blitzified tree, writing out the terminal symbols as you go, creates our new blitz++ expression string.

1.16.8.4.2. Blitz and NumPy¶

The other nice discovery in the project is that the data structure used for NumPy arrays and blitz arrays is nearly identical. NumPy stores “strides” as byte offsets and blitz stores them as element offsets, but other than that, they are the same. Further, most of the concept and capabilities of the two libraries are remarkably similar. It is satisfying that two completely different implementations solved the problem with similar basic architectures. It is also fortuitous. The work involved in converting NumPy expressions to blitz expressions was greatly diminished. As an example, consider the code for slicing an array in Python with a stride:

>>> a = b[0:4:2] + c
>>> a
[0,2,4]

In Blitz it is as follows:

Array<2,int> b(10);
Array<2,int> c(3);
// ...
Array<2,int> a = b(Range(0,3,2)) + c;

Here the range object works exactly like Python slice objects with the exception that the top index (3) is inclusive where as Python’s (4) is exclusive. Other differences include the type declarations in C++ and parentheses instead of brackets for indexing arrays. Currently, weave.blitz handles the inclusive/exclusive issue by subtracting one from upper indices during the translation. An alternative that is likely more robust/maintainable in the long run is to write a PyRange class that behaves like Python’s range. This is likely very easy.

The stock blitz also doesn’t handle negative indices in ranges. The current implementation of the blitz() has a partial solution to this problem. It calculates and index that starts with a ‘-‘ sign by subtracting it from the maximum index in the array so that:

                upper index limit
                    /-----\
b[:-1] -> b(Range(0,Nb[0]-1-1))

This approach fails, however, when the top index is calculated from other values. In the following scenario, if i+j evaluates to a negative value, the compiled code will produce incorrect results and could even core-dump. Right now, all calculated indices are assumed to be positive.

b[:i-j] -> b(Range(0,i+j))

A solution is to calculate all indices up front using if/then to handle the +/- cases. This is a little work and results in more code, so it hasn’t been done. I’m holding out to see if blitz++ can be modified to handle negative indexing, but haven’t looked into how much effort is involved yet. While it needs fixin’, I don’t think there is a ton of code where this is an issue.

The actual translation of the Python expressions to blitz expressions is currently a two part process. First, all x:y:z slicing expression are removed from the AST, converted to slice(x,y,z) and re-inserted into the tree. Any math needed on these expressions (subtracting from the maximum index, etc.) are also performed here. _beg and _end are used as special variables that are defined as blitz::fromBegin and blitz::toEnd.

a[i+j:i+j+1,:] = b[2:3,:]

becomes a more verbose:

a[slice(i+j,i+j+1),slice(_beg,_end)] = b[slice(2,3),slice(_beg,_end)]

The second part does a simple string search/replace to convert to a blitz expression with the following translations:

slice(_beg,_end) -> _all  # not strictly needed, but cuts down on code.
slice            -> blitz::Range
[                -> (
]                -> )
_stp             -> 1

_all is defined in the compiled function as blitz::Range.all(). These translations could of course happen directly in the syntax tree. But the string replacement is slightly easier. Note that namespaces are maintained in the C++ code to lessen the likelihood of name clashes. Currently no effort is made to detect name clashes. A good rule of thumb is don’t use values that start with ‘_’ or ‘py_’ in compiled expressions and you’ll be fine.