 |
 |
 |
SAGE Programming Guide |  |
 |
 |
|
|
|
SAGE Programming Guide | |
|
|
This chapter assumes that you have a C compiler suite (e.g., GCC) installed on your computer, since you should not write compiled extension code on a computer without a C compiler. If you are using OS X, download and install the latest version of XCode.5.1 For Linux, use your package manager to install GCC. If you're using Microsoft Windows, then you should be using SAGE via colinux, which includes GCC by default.
Interpreted code: way easier to debug, trace, document, structure (multiple inheritence), etc.; good when few for loops; lots of structure and calls to compiled code.
Compiled code: two reasons to use--make fast code and make optimal use of libraries written in a compiled language; hard to write and debug; takes longer to write. The methods explained in this chapter can produce code where speed is only limited by your choice of algorithm.
TODO:
The primary goal of this chapter is to provide excellent documentation about how to use Pyrex to write very fast code in the context of SAGE. This chapter is self contained in that the reader should not have to look at any of the documentation distributed with Pyrex. The reader is urged to at some point read the Python/C API Reference manual.
There is also much to be learned from the Pyrex mailing list [#!pyrex:list!#].
This fusion of interpreted and compiled code is extremely natural for mathematics software; some algorithms are much better implemented in an interpreter because all time-critical steps involve low level arithmetic -- other algorithms, e.g., matrix multiplication, must be implemented in a compiled language in order to perform optimally. Also, existing compiled libraries are sometimes of very high quality, and a compiled language is needed to create the best possible interface to such libraries. It's crucial that both approaches to programming be fully supported. When deciding how to implement SAGE, I searched for months for an environment that could support both approaches. Python and Pyrex provide exactly this in a way that I believe is much easier conceptually than the implementation models of any of the other systems mentioned above.
In 2004, Stein considered using SWIG for SAGE extensively, and implemented several basic types using this model. The idea was to write code in C++ for SAGE that needed to be fast, then wrap it in SWIG. This ground to a halt, because the result was not sufficiently fast. First, there is overhead when writing code in C++ in the first place. Second, SWIG generates several layers of code between Python and the code that does the actual work; for implementing any of SAGE's numerous basic arithmetic type, e.g., integers, the overhead is completely unacceptable, at least if one wants code that is as fast as Magma and PARI. For example, in my initial implementation of integer arithmetic using SWIG and GMP, there were the following six levels between the SAGE interpreter and the GMP C library:
| Python code to provide a clean interface |
| SWIG Autogenerated Python code |
| SWIG Autogenerated C++ extension code |
| Handcode C++ Integer class |
| GMP's C++ Interface |
| GMP C Library |
| Pyrex generated C code that provides everything needed |
| GMP C Library |
Moreover, in my own experience with SAGE, the speedups from Psycho were not significant, except on code that had an obvious C equivalent. Such code can be easily rewritten in Pyrex, which always results in something that is faster, more robust, and more reliable than what Pyscho produces, and which can be debugged using gdb.
Though Psycho is not included with or used by SAGE, it is possible to try out in Psycho by installing it and putting the appropriate line at the top of a file.
def hello(name):
"""
Print hello with the given name.
"""
print("Hello %s"%name)
sage: load "hello.spyx"
Compiling hello.spyx...
sage: hello('World')
Hello World
sage: attach "hello.spyx"
SAGE integers. For example, you want to
do the following computation (timings on a Macbook Pro 2Ghz
under OS X):
sage: a = 939393; b = 19393; c = 39390283 sage: time for _ in range(10000): a = a*b + c CPU times: user 0.33 s, sys: 0.21 s, total: 0.53 s
There is no ``multiply then add'' operation built into SAGE, and you're very frustrated because you really want this to be as fast as possible. There is possibly significant overhead in the computation of a*b + c in Python.
First create a new compiled ``multiply and add'' command in SAGE.
cimport sage.ext.integer
def muladd(sage.ext.integer.Integer a,
sage.ext.integer.Integer b,
sage.ext.integer.Integer c):
"""
Compute a*b + c
"""
cdef sage.ext.integer.Integer d
d = sage.ext.integer.Integer()
mpz_mul(d.value, a.value, b.value)
mpz_add(d.value, d.value, c.value)
return d
Do this either by putting the above code in a spyx file or pasting it into the notebook in a cell that starts with %spyx. (See Section above.)
Now we can do the above computation more quickly:
sage: a = 939393; b = 19393; c = 39390283 sage: time for _ in range(10000): a = muladd(a,b,c) CPU times: user 0.11 s, sys: 0.05 s, total: 0.15 s
To squeeze out even more performance we can put the entire loop in single function.
def muladdloop(sage.ext.integer.Integer a,
sage.ext.integer.Integer b,
sage.ext.integer.Integer c,
int n):
cdef int i
cdef mpz_t t
cdef sage.ext.integer.Integer aa
aa = sage.ext.integer.Integer(a)
mpz_init(t)
for i from 0 <= i < n:
mpz_mul(t, aa.value, b.value)
mpz_add(aa.value, t, c.value)
mpz_clear(t)
return aa
We then have
sage: a = 939393; b = 19393; c = 39390283 sage: time z=muladdloop(a,b,c,10000) CPU times: user 0.11 s, sys: 0.00 s, total: 0.11 s
%pyrex
two = int(2)
def sumsquarespy(n):
return sum(i**two for i in xrange(1,n+1))
#<
time v=[sumsquarespy(100) for _ in xrange(10000)]
#> CPU time: 0.49 s, Wall time: 0.49 s
def sumsquares(int n):
cdef int i, j
j = 0
for i from 1 <= i <= n:
j = j + i*i
return j
#<
time v=[sumsquares(100) for _ in xrange(10000)]
#> CPU time: 0.06 s, Wall time: 0.06 s
Test Drive wrote: > I have tried "cdef > object j", however in that case users will not be able to use the j when > user will create an instance of the object Test. You want 'cdef public int j' to expose it to Python code. Or 'cdef readonly int j' if you don't want them to be able to write to it. Like builtin types, Pyrex extension types don't have a dict, so you can't assign arbitrary attributes to them.
def allows multiple inheritance, but no c-level data structures. Maybe creation is faster. What exactly is the deal here?!
very powerful, but use with care.
evidently there is bad overhead with them on OS X. not so bad in linux.
If we can find the module name where foo is defined, then check if the appropriate foo.pyx file is available. Load it. Find the function. Figure out which if it is defined twice. Show all if we can't figure it out.
If totally screwed and we can't even find the file, then modify pyrex so it embeds extra info in the docstring. Strip that out on displaying.
type checking
Etc.
How to do explicit casts
def __getslice__(self, Py_ssize_t i, Py_ssize_t j):
...
The object type.
Pyrex code compiles with C++ compiler. Using pyrexembed.
def addit(v, s):
"""
EXAMPLE:
sage: load a.spyx
Compiling a.spyx...
sage: addit([1,2,3], 5)
[5, 10, 15]
"""
return eval('[x*s for x in v]', {'v':v, 's':s})
|
|
|
SAGE Programming Guide | |
|
|