Ticket #5093: trac5093-fast-callable.patch

File trac5093-fast-callable.patch, 216.1 KB (added by cwitty, 14 years ago)
  • .hgignore

    # HG changeset patch
    # User Carl Witty <cwitty@newtonlabs.com>
    # Date 1235880780 28800
    # Node ID 0949e5d7cc4cc9bade2b54973a1fedf2005b49c9
    # Parent  5ba8b66a0ae4da7d663e961da9e50f2fcbd2b829
    Add fast_callable (rewritten clone of fast_eval)
    
    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c .hgignore
    a b  
    110110^c_lib/config.status$
    111111^c_lib/libtool$
    112112^c_lib/libcsage.la$
     113^sage/ext/interpreters/
    113114^sage/libs/ntl/misc.h$
    114115^sage/libs/ntl/ntl.h$
    115116^sage/rings/integer.h$
  • module_list.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c module_list.py
    a b  
    151151    ##
    152152    ################################
    153153
     154    Extension('sage.ext.fast_callable',
     155              sources = ['sage/ext/fast_callable.pyx']),
     156   
    154157    Extension('sage.ext.fast_eval',
    155158              sources = ['sage/ext/fast_eval.pyx']),
    156159   
  • sage/all.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/all.py
    a b  
    126126
    127127from sage.parallel.all   import *
    128128
     129from sage.ext.fast_callable  import fast_callable
     130from sage.ext.fast_eval      import fast_float
     131
    129132from copy import copy, deepcopy
    130133
    131134# The code executed here uses a large amount of Sage components
  • sage/calculus/calculus.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/calculus/calculus.py
    a b  
    47414741            return fast_float.fast_float_constant(float(self))
    47424742        except:
    47434743            raise NotImplementedError # return lambda x: float(self(x))
    4744        
     4744
    47454745
    47464746class Symbolic_object(SymbolicExpression):
    47474747    r"""
     
    47684768    def _fast_float_(self, *vars):
    47694769        return fast_float.fast_float_constant(float(self))
    47704770
     4771    def _fast_callable_(self, etb):
     4772        r"""
     4773        Given an ExpressionTreeBuilder, return an Expression representing
     4774        this value.
     4775
     4776        EXAMPLES::
     4777
     4778            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     4779            sage: etb = ExpressionTreeBuilder(vars=['x'])
     4780            sage: SR(pi)._fast_callable_(etb)
     4781            pi
     4782            sage: etb = ExpressionTreeBuilder(vars=['x'], domain=RDF)
     4783            sage: SR(pi)._fast_callable_(etb)
     4784            3.14159265359
     4785        """
     4786        return etb.constant(self._obj)
     4787
    47714788    def __float__(self):
    47724789        """
    47734790        EXAMPLES:
     
    55925609        fops = [op._fast_float_(*vars) for op in self._operands]
    55935610        return self._operator(*fops)
    55945611
     5612    def _fast_callable_(self, etb):
     5613        r"""
     5614        Given an ExpressionTreeBuilder, return an Expression representing
     5615        this value.
     5616
     5617        EXAMPLES::
     5618
     5619            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     5620            sage: etb = ExpressionTreeBuilder(vars=['x','y'])
     5621            sage: var('x,y')
     5622            (x, y)
     5623            sage: (x+y)._fast_callable_(etb)
     5624            add(v_0, v_1)
     5625            sage: (-x)._fast_callable_(etb)
     5626            neg(v_0)
     5627        """
     5628        fops = [op._fast_callable_(etb) for op in self._operands]
     5629        return self._operator(*fops)
     5630
    55955631    def _convert(self, typ):
    55965632        """
    55975633        Convert self to the given type by converting each of the operands
     
    61686204        except TypeError:
    61696205            raise ValueError, "free variable: %s" % self._name
    61706206
     6207    def _fast_callable_(self, etb):
     6208        r"""
     6209        Given an ExpressionTreeBuilder, return an Expression representing
     6210        this value.
     6211
     6212        EXAMPLES::
     6213
     6214            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     6215            sage: etb = ExpressionTreeBuilder(vars=['x','y'])
     6216            sage: var('x,y,z')
     6217            (x, y, z)
     6218            sage: x._fast_callable_(etb)
     6219            v_0
     6220            sage: y._fast_callable_(etb)
     6221            v_1
     6222            sage: z._fast_callable_(etb)
     6223            Traceback (most recent call last):
     6224            ...
     6225            ValueError: list.index(x): x not in list
     6226        """
     6227        return etb.var(self)
     6228
    61716229    def _recursive_sub(self, kwds):
    61726230        # do the replacement if needed
    61736231        if kwds.has_key(self):
     
    66756733            vars = self.arguments()
    66766734        return self._expr._fast_float_(*vars)
    66776735
     6736    def _fast_callable_(self, etb):
     6737        r"""
     6738        Given an ExpressionTreeBuilder, return an Expression representing
     6739        this value.
     6740
     6741        EXAMPLES::
     6742
     6743            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     6744            sage: etb = ExpressionTreeBuilder(vars=['x','y'])
     6745            sage: g(x) = sin(x) + 2
     6746            sage: g._fast_callable_(etb)
     6747            add(sin(v_0), 2)
     6748        """
     6749        return self._expr._fast_callable_(etb)
     6750
    66786751    def __float__(self):
    66796752        return float(self._expr)
    66806753
     
    72257298            else:
    72267299                return fast_float.fast_float_func(f, g)
    72277300
     7301    def _fast_callable_(self, etb):
     7302        r"""
     7303        Given an ExpressionTreeBuilder, return an Expression representing
     7304        this value.
     7305
     7306        EXAMPLES::
     7307
     7308            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     7309            sage: etb = ExpressionTreeBuilder(vars=['x','y'])
     7310            sage: var('x,y')
     7311            (x, y)
     7312            sage: sin(sqrt(x+y))._fast_callable_(etb)
     7313            sin(sqrt(add(v_0, v_1)))
     7314        """
     7315        return etb.call(self._operands[0], etb(self._operands[1]))
     7316
    72287317    def __complex__(self):
    72297318        """
    72307319        Convert this symbolic composition to a Python complex number.
     
    75067595       
    75077596            sage: from sage.ext.fast_eval import fast_float
    75087597            sage: fast_float(sin)
    7509             <sage.ext.fast_eval.FastDoubleFunc object at 0x...>
     7598            <sage.ext... object at 0x...>
    75107599            sage: sin._fast_float_()
    75117600            <sage.ext.fast_eval.FastDoubleFunc object at 0x...>
    75127601            sage: sin._fast_float_()(0)
     
    75377626        except TypeError:
    75387627            return fast_float.fast_float_func(self, *args)
    75397628
     7629    def _fast_callable_(self, etb):
     7630        r"""
     7631        Given an ExpressionTreeBuilder, return an Expression representing
     7632        this value.
     7633
     7634        EXAMPLES::
     7635
     7636            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     7637            sage: etb = ExpressionTreeBuilder(vars=['x','y'])
     7638            sage: sin._fast_callable_(etb)
     7639            sin(v_0)
     7640            sage: erf._fast_callable_(etb)
     7641            {erf}(v_0)
     7642        """
     7643        args = [etb._var_number(n) for n in range(self.number_of_arguments())]
     7644        return etb.call(self, *args)
     7645
    75407646_syms = {}
    75417647
    75427648class Function_erf(PrimitiveFunction):
  • new file sage/ext/fast_callable.pxd

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/ext/fast_callable.pxd
    - +  
     1cdef class Wrapper:
     2    cdef readonly object _orig_args
     3    cdef readonly object _metadata
  • new file sage/ext/fast_callable.pyx

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/ext/fast_callable.pyx
    - +  
     1r"""
     2Fast Expression Evaluation.
     3
     4For many applications such as numerical integration, differential
     5equation approximation, plotting a 3d surface, optimization problems,
     6monte-carlo simulations, etc., one wishes to pass around and evaluate
     7a single algebraic expression many, many times at various floating
     8point values.  Other applications may need to evaluate an expression
     9many times in interval arithmetic, or in a finite field.  Doing this
     10via recursive calls over a python representation of the object (even
     11if maxima or other outside packages are not involved) is extremely
     12inefficient.
     13
     14This module provides a function, \function{fast_callable}, to
     15transform such expressions into a form where they can be evaluated
     16quickly.
     17
     18sage: f = sin(x) + 3*x^2
     19sage: ff = fast_callable(f)
     20sage: ff(3.5)
     2136.3992167723104
     22sage: ff(RIF(3.5))
     2336.39921677231038?
     24
     25By default, \function{fast_callable} only removes some interpretive
     26overhead from the evaluation, but all of the individual arithmetic
     27operations are done using standard \sage arithmetic.  This is still a
     28huge win over sage.calculus, which evidently has a lot of overhead.
     29Compare the cost of evaluating Wilkinson's polynomial (in unexpanded
     30form) at x=30:
     31
     32sage: wilk = prod((x-i) for i in [1 .. 20]); wilk
     33(x - 20)*(x - 19)*(x - 18)*(x - 17)*(x - 16)*(x - 15)*(x - 14)*(x - 13)*(x - 12)*(x - 11)*(x - 10)*(x - 9)*(x - 8)*(x - 7)*(x - 6)*(x - 5)*(x - 4)*(x - 3)*(x - 2)*(x - 1)
     34sage: timeit('wilk.subs(x=30)') # random, long time
     35625 loops, best of 3: 1.46 ms per loop
     36sage: fc_wilk = fast_callable(wilk)
     37sage: timeit('fc_wilk(30)') # random, long time
     38625 loops, best of 3: 10.4 us per loop
     39
     40You can specify a particular domain for the evaluation using
     41\code{domain=}:
     42
     43sage: fc_wilk_zz = fast_callable(wilk, domain=ZZ)
     44
     45The meaning of domain=D is that each intermediate and final result
     46is converted to type D.  For instance, the previous example of
     47``sin(x) + 3*x^2`` with domain=D would be equivalent to
     48``D(D(sin(D(x))) + D(D(3)*D(D(x)^2)))``.  (This example also
     49demonstrates the one exception to the general rule: if an exponent is an
     50integral constant, then it is not wrapped with D().)
     51
     52At first glance, this seems like a very bad idea if you want to
     53compute quickly.  And it is a bad idea, for types where we don't
     54have a special interpreter.  It's not too bad of a slowdown, though.
     55To mitigate the costs, we check whether the value already has
     56the correct parent before we call D.
     57
     58We don't yet have a special interpreter with domain ZZ, so we can see
     59how that compares to the generic fc_wilk example above:
     60
     61sage: timeit('fc_wilk_zz(30)') # random, long time
     62625 loops, best of 3: 15.6 us per loop
     63
     64However, for other types, using domain=D will get a large speedup,
     65because we have special-purpose interpreters for those types.  One
     66example is RDF.  Since with domain=RDF we know that every single
     67operation will be floating-point, we can just execute the
     68floating-point operations directly and skip all the Python object
     69creations that you would get from actually using RDF objects.
     70
     71sage: fc_wilk_rdf = fast_callable(wilk, domain=RDF)
     72sage: timeit('fc_wilk_rdf(30.0)') # random, long time
     73625 loops, best of 3: 5.13 us per loop
     74
     75We also have support for RR:
     76
     77sage: fc_wilk_rr = fast_callable(wilk, domain=RR)
     78sage: timeit('fc_wilk_rr(30.0)') # random, long time
     79625 loops, best of 3: 12.9 us per loop
     80
     81By default, \function{fast_callable} uses the same variable names in the
     82same order that the \method{__call__} method on its argument would use;
     83for instance, on symbolic expressions, the variables are used in alphabetical
     84order.
     85
     86sage: var('y,z,x')
     87(y, z, x)
     88sage: f = 10*y + 100*z + x
     89sage: f(1,2,3)
     90321
     91sage: ff = fast_callable(f)
     92sage: ff(1,2,3)
     93321
     94
     95However, this can be overridden with \code{vars=}:
     96
     97sage: ff = fast_callable(f, vars=('x','w','z','y'))
     98sage: ff(1,2,3,4)
     99341
     100
     101This should be enough for normal use of \function{fast_callable}; let's
     102discuss some more advanced topics.
     103
     104Sometimes it may be useful to create a fast version of an expression
     105without going through symbolic expressions or polynomials; perhaps
     106because you want to describe to \function{fast_callable} an expression
     107with common subexpressions.
     108
     109Internally, \function{fast_callable} works in two stages: it constructs
     110an expression tree from its argument, and then it builds a
     111fast evaluator from that expression tree.  You can bypass the first phase
     112by building your own expression tree and passing that directly to
     113\function{fast_callable}, using an \class{ExpressionTreeBuilder}.
     114
     115sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     116sage: etb = ExpressionTreeBuilder(vars=('x','y','z'))
     117
     118An \class{ExpressionTreeBuilder} has three interesting methods:
     119\method{constant}, \method{var}, and \method{call}.
     120All of these methods return \class{ExpressionTree} objects.
     121
     122The \method{var} method takes a string, and returns an expression tree
     123for the corresponding variable.
     124
     125sage: x = etb.var('x')
     126sage: y = etb.var('y')
     127sage: z = etb.var('y')
     128
     129Expression trees support Python's numeric operators, so you can easily
     130build expression trees representing arithmetic expressions.
     131
     132sage: v1 = (x+y)*(y+z) + (y//z)
     133
     134The \method{constant} method takes a \sage value, and returns an
     135expression tree representing that value. 
     136
     137sage: v2 = etb.constant(3.14159) * x + etb.constant(1729) * y
     138
     139The \method{call} method takes a \sage/Python function and zero or more
     140expression trees, and returns an expression tree representing
     141the function call.
     142
     143sage: v3 = etb.call(sin, v1+v2)
     144sage: v3
     145sin(add(add(mul(add(v_0, v_1), add(v_1, v_1)), floordiv(v_1, v_1)), add(mul(3.14159000000000, v_0), mul(1729, v_1))))
     146
     147Many \sage/Python built-in functions are specially handled; for instance,
     148when evaluating an expression involving \function{sin} over \code{RDF},
     149the C math library function \function{sin} is called.  Arbitrary functions
     150are allowed, but will be much slower since they will call back to
     151Python code on every call; for example, the following will work.
     152
     153sage: def my_sqrt(x): return pow(x, 0.5)
     154sage: e = etb.call(my_sqrt, v1); e
     155{my_sqrt}(add(mul(add(v_0, v_1), add(v_1, v_1)), floordiv(v_1, v_1)))
     156sage: fast_callable(e)(1, 2, 3)
     1573.60555127546399
     158
     159To provide \function{fast_callable} for your own class (so that
     160\code{fast_callable(x)} works when \variable{x} is an instance of your
     161class), implement a method \code{_fast_callable_(self, etb)} for your class.
     162This method takes an \class{ExpressionTreeBuilder}, and returns an
     163expression tree built up using the methods described above.
     164
     165EXAMPLES:
     166    sage: var('x')
     167    x
     168    sage: f = fast_callable(sqrt(x^7+1), domain=RDF)
     169
     170    sage: f(1)
     171    1.4142135623730951
     172    sage: f.op_list()
     173    [('load_arg', 0), ('load_const', 7.0), 'pow', ('load_const', 1.0), 'add', 'sqrt', 'return']
     174   
     175    To interpret that last line, we load argument 0 ('x' in this case) onto
     176    the stack, push the constant 7.0 onto the stack, call the pow function
     177    (which takes 2 arguments from the stack), push the constant 1.0, add the
     178    top two arguments of the stack, and then call sqrt.
     179
     180Here we take sin of the first argument and add it to f:
     181
     182    sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     183    sage: etb = ExpressionTreeBuilder('x')
     184    sage: x = etb.var('x')
     185    sage: f = etb.call(sqrt, x^7 + 1)
     186    sage: g = etb.call(sin, x)
     187    sage: fast_callable(f+g).op_list()
     188    [('load_arg', 0), ('load_const', 7), 'pow', ('load_const', 1), 'add', ('py_call', sqrt, 1), ('load_arg', 0), ('py_call', sin, 1), 'add', 'return']
     189
     190
     191AUTHOR:
     192    -- Carl Witty (2009-02): initial version (heavily inspired by
     193       Robert Bradshaw's fast_eval.pyx)
     194"""
     195
     196
     197#*****************************************************************************
     198#       Copyright (C) 2008 Robert Bradshaw <robertwb@math.washington.edu>
     199#       Copyright (C) 2009 Carl Witty <Carl.Witty@gmail.com>
     200#
     201#  Distributed under the terms of the GNU General Public License (GPL)
     202#
     203#    This code is distributed in the hope that it will be useful,
     204#    but WITHOUT ANY WARRANTY; without even the implied warranty of
     205#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
     206#    General Public License for more details.
     207#
     208#  The full text of the GPL is available at:
     209#
     210#                  http://www.gnu.org/licenses/
     211#*****************************************************************************
     212
     213# The following bits of text were written for the module docstring.
     214# They are not true yet, but I hope they will be true someday, at
     215# which point I will move them into the docstring.
     216#------------------------------ WRONG (for now) docs follow
     217# The final interesting method of \class{ExpressionTreeBuilder} is
     218# \method{choice}.  This produces conditional expressions, like the C
     219# \code{COND ? T : F} expression or the Python {T if COND else F}.
     220# This lets you define piecewise functions using \function{fast_callable}.
     221
     222# sage: v4 = etb.choice(v3 >= etb.constant(0), v1, v2)
     223
     224# The arguments are \code{(COND, T, F)} (the same order as in C), so the
     225# above means that if \variable{v3} evaluates to a nonnegative number,
     226# then \variable{v4} will evaluate to the result of \variable{v1};
     227# otherwise, \variable{v4} will evaluate to the result of \variable{v2}.
     228
     229# Let's see an example where we see that \function{fast_callable} does not
     230# evaluate common subexpressions more than once.  We'll make a
     231# \function{fast_callable} expression that gives the result
     232# of 16 iterations of the Mandelbrot function.
     233
     234# sage: etb = ExpressionTreeBuilder('c')
     235# sage: z = etb.constant(0)
     236# sage: c = etb.var('c')
     237# sage: for i in range(16):
     238# ...       z = z*z + c
     239# sage: mand = fast_callable(z, domain=CDF) # not tested
     240
     241# Now \variable{ff} does 32 complex arithmetic operations on each call
     242# (16 additions and 16 multiplications).  However, if \code{z*z} produced
     243# code that evaluated \variable{z} twice, then this would do many
     244# thousands of arithmetic operations instead.
     245
     246# Note that the handling for common subexpressions only checks whether
     247# expression trees are the same Python object; for instance, the following
     248# code will evaluate \code{x+1} twice:
     249
     250# sage: etb = ExpressionTreeBuilder('x')
     251# sage: x = etb.var('x')
     252# sage: (x+1)*(x+1)
     253# *(+(v_0, 1), +(v_0, 1))
     254
     255# but this code will only evaluate \code{x+1} once:
     256
     257# sage: v = x+1; v*v
     258# *(+(v_0, 1), +(v_0, 1))
     259#------------------------------ done with WRONG (for now) docs
     260
     261
     262import operator
     263from copy import copy
     264from sage.rings.real_mpfr cimport RealField, RealNumber
     265from sage.structure.element cimport Element
     266from sage.rings.all import RDF
     267from sage.libs.mpfr cimport mpfr_t, mpfr_ptr, mpfr_init2, mpfr_set, GMP_RNDN
     268
     269include "stdsage.pxi"
     270
     271def fast_callable(x, domain=None, vars=None):
     272    r"""
     273    Given an expression x, compiles it into a form that can be quickly
     274    evaluated, given values for the variables in x.
     275
     276    x can be an expression object, or anything that has a .variables()
     277    method and a ._fast_callable_() method (this includes SR, univariate
     278    polynomials, and multivariate polynomials).
     279
     280    By default, x is evaluated the same way that a Python function would
     281    evaluate it -- addition maps to PyNumber_Add, etc.  However, you
     282    can specify domain=D where D is some Sage parent; in this case,
     283    all arithmetic is done in that parent.  If we have a special-purpose
     284    interpreter for that parent (like RDF), domain=RDF will trigger the
     285    use of that interpreter.
     286
     287    If vars is None, then we will attempt to determine the set of
     288    variables from x; otherwise, we will use the given set.
     289
     290    EXAMPLES:
     291        sage: var('x')
     292        x
     293        sage: expr = sin(x) + 3*x^2
     294        sage: f = fast_callable(expr)
     295        sage: f(2)
     296        sin(2) + 12
     297        sage: f(2.0)
     298        12.9092974268257
     299        sage: f_rdf = fast_callable(expr, domain=RDF)
     300        sage: f_rdf(2)
     301        12.909297426825681
     302        sage: f = fast_callable(expr, vars=('z','x','y'))
     303        sage: f(1, 2, 3)
     304        sin(2) + 12
     305        sage: K.<x> = QQ[]
     306        sage: p = K.random_element(6); p
     307        -x^6 - 12*x^5 + 1/2*x^4 - 1/95*x^3 - 1/2*x^2 - 4
     308        sage: fp = fast_callable(p, domain=RDF)
     309        sage: fp.op_list()
     310        [('load_arg', 0), ('load_const', -1.0), 'mul', ('load_const', -12.0), 'add', ('load_arg', 0), 'mul', ('load_const', 0.5), 'add', ('load_arg', 0), 'mul', ('load_const', -0.0105263157895), 'add', ('load_arg', 0), 'mul', ('load_const', -0.5), 'add', ('load_arg', 0), 'mul', ('load_arg', 0), 'mul', ('load_const', -4.0), 'add', 'return']
     311        sage: fp(3.14159)
     312        -4594.1618236401764
     313        sage: K.<x,y,z> = QQ[]
     314        sage: p = K.random_element(degree=3, terms=5); p
     315        -x*y^2 - x*z^2 - 6*x^2 - y^2 - 3*x*z
     316        sage: fp = fast_callable(p, domain=RDF)
     317        sage: fp.op_list()
     318        [('load_const', 0.0), ('load_const', -3.0), ('load_arg', 0), ('load_const', 1.0), 'pow', ('load_arg', 2), ('load_const', 1.0), 'pow', 'mul', 'mul', 'add', ('load_const', -1.0), ('load_arg', 0), ('load_const', 1.0), 'pow', ('load_arg', 1), ('load_const', 2.0), 'pow', 'mul', 'mul', 'add', ('load_const', -6.0), ('load_arg', 0), ('load_const', 2.0), 'pow', 'mul', 'add', ('load_const', -1.0), ('load_arg', 1), ('load_const', 2.0), 'pow', 'mul', 'add', ('load_const', -1.0), ('load_arg', 0), ('load_const', 1.0), 'pow', ('load_arg', 2), ('load_const', 2.0), 'pow', 'mul', 'mul', 'add', 'return']
     319        sage: fp(e, pi, sqrt(2))
     320        -98.001564033629322
     321        sage: symbolic_result = p(e, pi, sqrt(2)); symbolic_result
     322        -1*e*pi^2 - pi^2 - 6*e^2 - 3*sqrt(2)*e - 2*e
     323        sage: n(symbolic_result)
     324        -98.0015640336293
     325
     326        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     327        sage: etb = ExpressionTreeBuilder(vars=('x','y'), domain=RDF)
     328        sage: x = etb.var('x')
     329        sage: y = etb.var('y')
     330        sage: expr = etb.call(sin, x^2 + y); expr
     331        sin(add(pow(v_0, 2.0), v_1))
     332        sage: fc = fast_callable(expr, domain=RDF)
     333        sage: fc(5, 7)
     334        0.55142668124169059
     335    """
     336    cdef Expression et
     337    if isinstance(x, Expression):
     338        et = x
     339        vars = et._etb._vars
     340    else:
     341        if vars is None or len(vars) == 0:
     342            vars = x.variables()
     343            # XXX This is pretty gross... there should be a "callable_variables"
     344            # method that does all this.
     345            from sage.calculus.calculus import SR
     346            from sage.rings.all import is_PolynomialRing, is_MPolynomialRing
     347            if x.parent() is SR and x.number_of_arguments() > len(vars):
     348                vars = list(vars) + ['EXTRA_VAR%d' % n for n in range(len(vars), x.number_of_arguments())]
     349            if is_PolynomialRing(x.parent()) or is_MPolynomialRing(x.parent()):
     350                vars = x.parent().variable_names()
     351        etb = ExpressionTreeBuilder(vars=vars, domain=domain)
     352        et = x._fast_callable_(etb)
     353    if isinstance(domain, RealField):
     354        import sage.ext.interpreters.wrapper_rr
     355        builder = sage.ext.interpreters.wrapper_rr.Wrapper_rr
     356
     357        str = InstructionStream(sage.ext.interpreters.wrapper_rr.metadata,
     358                                len(vars),
     359                                domain)
     360    elif domain == RDF:
     361        import sage.ext.interpreters.wrapper_rdf
     362        builder = sage.ext.interpreters.wrapper_rdf.Wrapper_rdf
     363        str = InstructionStream(sage.ext.interpreters.wrapper_rdf.metadata,
     364                                len(vars))
     365    elif domain is None:
     366        import sage.ext.interpreters.wrapper_py
     367        builder = sage.ext.interpreters.wrapper_py.Wrapper_py
     368        str = InstructionStream(sage.ext.interpreters.wrapper_py.metadata,
     369                                len(vars))
     370    else:
     371        import sage.ext.interpreters.wrapper_el
     372        builder = sage.ext.interpreters.wrapper_el.Wrapper_el
     373        str = InstructionStream(sage.ext.interpreters.wrapper_el.metadata,
     374                                len(vars),
     375                                domain)
     376    generate_code(et, str)
     377    str.instr('return')
     378    return builder(str.get_current())
     379
     380def function_name(fn):
     381    r"""
     382    Given a function, returns a string giving a name for the function.
     383
     384    For functions we recognize, we use our standard opcode name for the
     385    function (so operator.add becomes 'add', and sage.all.sin becomes 'sin').
     386
     387    For functions we don't recognize, we try to come up with a name,
     388    but the name will be wrapped in braces; this is a signal that
     389    we'll definitely use a slow Python call to call this function.
     390    (We may use a slow Python call even for functions we do recognize,
     391    if we're targeting an interpreter without an opcode for the function.)
     392
     393    Only used when printing Expressions.
     394
     395    EXAMPLES:
     396        sage: from sage.ext.fast_callable import function_name
     397        sage: function_name(operator.pow)
     398        'pow'
     399        sage: function_name(cos)
     400        'cos'
     401        sage: function_name(factorial)
     402        '{factorial}'
     403    """
     404    builtins = get_builtin_functions()
     405    if fn in builtins:
     406        return builtins[fn]
     407    try:
     408        return "{%s}" % fn.__name__
     409    except AttributeError:
     410        return "{%r}" % fn
     411
     412cdef class ExpressionTreeBuilder:
     413    r"""
     414    A class with helper methods for building Expressions.
     415
     416    An instance of this class is passed to _fast_callable_ methods;
     417    you can also instantiate it yourself to create your own expressions
     418    for fast_callable, bypassing _fast_callable_.
     419
     420    EXAMPLES:
     421        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     422        sage: etb = ExpressionTreeBuilder('x')
     423        sage: x = etb.var('x')
     424        sage: (x+3)*5
     425        mul(add(v_0, 3), 5)
     426    """
     427
     428    cdef readonly object _domain
     429    cdef readonly object _vars
     430
     431    def __init__(self, vars, domain=None):
     432        r"""
     433        Initialize an instance of ExpressionTreeBuilder.  Takes
     434        a list or tuple of variable names to use, and also an optional
     435        domain.  If a domain is given, then creating an ExpressionConstant
     436        node with the __call__, make, or constant methods will convert
     437        the value into the given domain.
     438
     439        Note that this is the only effect of the domain parameter.  It
     440        is quite possible to use different domains for
     441        ExpressionTreeBuilder and for fast_callable; in that case,
     442        constants will be converted twice (once when building the
     443        Expression, and once when generating code).
     444
     445        EXAMPLES:
     446            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     447            sage: etb = ExpressionTreeBuilder('x')
     448            sage: etb(3^50)
     449            717897987691852588770249
     450            sage: etb = ExpressionTreeBuilder('x', domain=RR)
     451            sage: etb(3^50)
     452            7.17897987691853e23
     453        """
     454
     455        if isinstance(vars, tuple):
     456            vars = list(vars)
     457        elif not isinstance(vars, list):
     458            vars = [vars]
     459
     460        vars = map(self._clean_var, vars)
     461
     462        self._domain = domain
     463        self._vars = vars
     464
     465    def __call__(self, x):
     466        r"""
     467        Try to convert the given value to an Expression.  If it is already
     468        an Expression, just return it.  If it has a _fast_callable_
     469        method, then call the method with self as an argument.  Otherwise,
     470        use self.constant() to turn it into a constant.
     471
     472        EXAMPLES:
     473            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     474            sage: etb = ExpressionTreeBuilder('x')
     475            sage: v = etb(3); v, type(v)
     476            (3, <type 'sage.ext.fast_callable.ExpressionConstant'>)
     477            sage: v = etb(polygen(QQ)); v, type(v)
     478            (v_0, <type 'sage.ext.fast_callable.ExpressionVariable'>)
     479            sage: v is etb(v)
     480            True
     481        """
     482        if isinstance(x, Expression):
     483            return x
     484
     485        try:
     486            fc = x._fast_callable_
     487        except AttributeError:
     488            return self.constant(x)
     489
     490        return fc(self)
     491
     492    def _clean_var(self, v):
     493        r"""
     494        Give a variable name, given a variable.
     495
     496        EXAMPLES:
     497            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     498            sage: etb = ExpressionTreeBuilder('x')
     499            sage: var('x')
     500            x
     501            sage: etb._clean_var(x)
     502            'x'
     503            sage: x = polygen(RR); x
     504            1.00000000000000*x
     505            sage: etb._clean_var(x)
     506            'x'
     507        """
     508        # There should be a better way to do this.  (Maybe there is.)
     509        if not PY_TYPE_CHECK(v, str):
     510            v = str(v)
     511            if '*' in v:
     512                v = v[v.index('*')+1:]
     513        return v       
     514
     515    def constant(self, c):
     516        r"""
     517        Turn the argument into an ExpressionConstant, converting it to
     518        our domain if we have one.
     519
     520        EXAMPLES:
     521            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     522            sage: etb = ExpressionTreeBuilder('x')
     523            sage: etb.constant(pi)
     524            pi
     525            sage: etb = ExpressionTreeBuilder('x', domain=RealField(200))
     526            sage: etb.constant(pi)
     527            3.1415926535897932384626433832795028841971693993751058209749
     528        """
     529        if self._domain is not None:
     530            c = self._domain(c)
     531        return ExpressionConstant(self, c)
     532
     533    def var(self, v):
     534        r"""
     535        Turn the argument into an ExpressionVariable.  Looks it up in
     536        the list of variables.  (Variables are matched by name.)
     537
     538        EXAMPLES:
     539            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     540            sage: var('a,b,some_really_long_name')
     541            (a, b, some_really_long_name)
     542            sage: x = polygen(QQ)
     543            sage: etb = ExpressionTreeBuilder(vars=('a','b',some_really_long_name, x))
     544            sage: etb.var(some_really_long_name)
     545            v_2
     546            sage: etb.var('some_really_long_name')
     547            v_2
     548            sage: etb.var(x)
     549            v_3
     550            sage: etb.var('y')
     551            Traceback (most recent call last):
     552            ...
     553            ValueError: list.index(x): x not in list           
     554        """
     555        ind = self._vars.index(self._clean_var(v))
     556        return ExpressionVariable(self, ind)
     557
     558    def _var_number(self, n):
     559        r"""
     560        Given an integer, return the variable with that index.
     561
     562        EXAMPLES:
     563            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     564            sage: etb = ExpressionTreeBuilder(vars=('a','b','c','d'))
     565            sage: etb._var_number(0)
     566            v_0
     567            sage: etb._var_number(3)
     568            v_3
     569            sage: etb._var_number(4)
     570            Traceback (most recent call last):           
     571            ...
     572            ValueError: Variable number 4 out of range       
     573        """
     574        if 0 <= n < len(self._vars):
     575            return ExpressionVariable(self, n)
     576        raise ValueError, "Variable number %d out of range" % n
     577
     578    def call(self, fn, *args):
     579        r"""
     580        Construct a call node, given a function and a list of arguments.
     581        The arguments will be converted to Expressions using
     582        ExpressionTreeBuilder.__call__.
     583
     584        EXAMPLES:
     585            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     586            sage: etb = ExpressionTreeBuilder(vars=(x,))
     587            sage: etb.call(cos, x)
     588            cos(v_0)
     589            sage: etb.call(sin, 1)
     590            sin(1)
     591            sage: etb.call(sin, etb(1))
     592            sin(1)
     593            sage: etb.call(factorial, x+57)
     594            {factorial}(add(v_0, 57))
     595        """
     596        return ExpressionCall(self, fn, map(self, args))
     597
     598    def choice(self, cond, iftrue, iffalse):
     599        r"""
     600        Construct a choice node (a conditional expression), given the
     601        condition, and the values for the true and false cases.
     602
     603        (It's possible to create choice nodes, but they don't work yet.)
     604
     605        EXAMPLES:
     606            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     607            sage: etb = ExpressionTreeBuilder(vars=(x,))
     608            sage: etb.choice(etb.call(operator.eq, x, 0), 0, 1/x)
     609            (0 if {eq}(v_0, 0) else div(1, v_0))
     610        """
     611        return ExpressionChoice(self,
     612                                cond,
     613                                self(iftrue),
     614                                self(iffalse))
     615
     616# Cache these values, to make expression building a tiny bit faster
     617# (by skipping the hash-table lookup in the operator module).
     618cdef op_add = operator.add
     619cdef op_sub = operator.sub
     620cdef op_mul = operator.mul
     621cdef op_div = operator.div
     622cdef op_floordiv = operator.floordiv
     623cdef op_pow = operator.pow
     624cdef op_neg = operator.neg
     625cdef op_abs = operator.abs
     626cdef op_inv = operator.inv
     627
     628cdef class Expression:
     629    r"""
     630    Represents an expression for fast_callable.
     631
     632    Supports the standard Python arithmetic operators; if arithmetic
     633    is attempted between an Expression and a non-Expression, the
     634    non-Expression is converted to an expression (using the
     635    __call__ method of the Expression's ExpressionTreeBuilder).
     636
     637    EXAMPLES:
     638        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     639        sage: etb = ExpressionTreeBuilder(vars=(x,))
     640        sage: x = etb.var(x)
     641        sage: etb(x)
     642        v_0
     643        sage: etb(3)
     644        3
     645        sage: etb.call(sin, x)
     646        sin(v_0)
     647        sage: (x+1)/(x-1)
     648        div(add(v_0, 1), sub(v_0, 1))
     649        sage: x//5
     650        floordiv(v_0, 5)
     651        sage: -abs(~x)
     652        neg(abs(inv(v_0)))
     653    """
     654
     655    cdef ExpressionTreeBuilder _etb
     656
     657    def __init__(self, etb):
     658        r"""
     659        Initialize an Expression.  Sets the _etb member.
     660
     661        EXAMPLES:
     662            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     663            sage: etb = ExpressionTreeBuilder(vars=(x,))
     664            sage: v = etb(3); v # indirect doctest
     665            3
     666            sage: v._get_etb() is etb
     667            True
     668        """
     669        self._etb = etb
     670
     671    def _get_etb(self):
     672        r"""
     673        Returns the ExpressionTreeBuilder used to build a given expression.
     674
     675        EXAMPLES:
     676            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     677            sage: etb = ExpressionTreeBuilder(vars=(x,))
     678            sage: v = etb(3); v
     679            3
     680            sage: v._get_etb() is etb
     681            True
     682        """
     683        return self._etb
     684
     685    def __add__(s, o):
     686        r"""
     687        Compute a sum of two Expressions.
     688
     689        EXAMPLES:
     690            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     691            sage: etb = ExpressionTreeBuilder(vars=(x,))
     692            sage: x = etb(x)
     693            sage: x+x
     694            add(v_0, v_0)
     695            sage: x+1
     696            add(v_0, 1)
     697            sage: 1+x
     698            add(1, v_0)
     699            sage: x.__add__(1)
     700            add(v_0, 1)
     701            sage: x.__radd__(1)
     702            add(1, v_0)
     703        """
     704        return _expression_binop_helper(s, o, op_add)
     705
     706    def __sub__(s, o):
     707        r"""
     708        Compute a difference of two Expressions.
     709
     710        EXAMPLES:
     711            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     712            sage: etb = ExpressionTreeBuilder(vars=(x,))
     713            sage: x = etb(x)
     714            sage: x-x
     715            sub(v_0, v_0)
     716            sage: x-1
     717            sub(v_0, 1)
     718            sage: 1-x
     719            sub(1, v_0)
     720            sage: x.__sub__(1)
     721            sub(v_0, 1)
     722            sage: x.__rsub__(1)
     723            sub(1, v_0)
     724        """
     725        return _expression_binop_helper(s, o, op_sub)
     726
     727    def __mul__(s, o):
     728        r"""
     729        Compute a product of two Expressions.
     730
     731        EXAMPLES:
     732            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     733            sage: etb = ExpressionTreeBuilder(vars=(x,))
     734            sage: x = etb(x)
     735            sage: x*x
     736            mul(v_0, v_0)
     737            sage: x*1
     738            mul(v_0, 1)
     739            sage: 1*x
     740            mul(1, v_0)
     741            sage: x.__mul__(1)
     742            mul(v_0, 1)
     743            sage: x.__rmul__(1)
     744            mul(1, v_0)
     745        """
     746        return _expression_binop_helper(s, o, op_mul)
     747
     748    def __div__(s, o):
     749        r"""
     750        Compute a quotient of two Expressions.
     751
     752        EXAMPLES:
     753            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     754            sage: etb = ExpressionTreeBuilder(vars=(x,))
     755            sage: x = etb(x)
     756            sage: x/x
     757            div(v_0, v_0)
     758            sage: x/1
     759            div(v_0, 1)
     760            sage: 1/x
     761            div(1, v_0)
     762            sage: x.__div__(1)
     763            div(v_0, 1)
     764            sage: x.__rdiv__(1)
     765            div(1, v_0)
     766        """
     767        return _expression_binop_helper(s, o, op_div)
     768
     769    def __floordiv__(s, o):
     770        r"""
     771        Compute the floordiv (the floor of the quotient) of two Expressions.
     772
     773        EXAMPLES:
     774            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     775            sage: etb = ExpressionTreeBuilder(vars=(x,))
     776            sage: x = etb(x)
     777            sage: x//x
     778            floordiv(v_0, v_0)
     779            sage: x//1
     780            floordiv(v_0, 1)
     781            sage: 1//x
     782            floordiv(1, v_0)
     783            sage: x.__floordiv__(1)
     784            floordiv(v_0, 1)
     785            sage: x.__rfloordiv__(1)
     786            floordiv(1, v_0)
     787        """
     788        return _expression_binop_helper(s, o, op_floordiv)
     789
     790    def __pow__(s, o, dummy):
     791        r"""
     792        Compute a power expression from two Expressions.
     793
     794        EXAMPLES:
     795            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     796            sage: etb = ExpressionTreeBuilder(vars=(x,))
     797            sage: x = etb(x)
     798            sage: x^x
     799            pow(v_0, v_0)
     800            sage: x^1
     801            pow(v_0, 1)
     802            sage: x.__pow__(1)
     803            pow(v_0, 1)
     804            sage: x.__rpow__(1)
     805            pow(1, v_0)
     806        """
     807        # XXX There is a performance regression from the original
     808        # fast_float here; it would replace small integer powers with
     809        # multiplication.  We can't do this safely until we support
     810        # common subexpression elimination (or at least the dup instruction).
     811        # (Plus, we should consider how strict a semantics we want;
     812        # probably this sort of optimization should be controlled by a
     813        # flag.)
     814        return _expression_binop_helper(s, o, op_pow)
     815       
     816    def __neg__(self):
     817        r"""
     818        Compute the negation of an Expression.
     819
     820        EXAMPLES:
     821            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     822            sage: etb = ExpressionTreeBuilder(vars=(x,))
     823            sage: x = etb(x)
     824            sage: -x
     825            neg(v_0)
     826            sage: x.__neg__()
     827            neg(v_0)
     828        """
     829        return ExpressionCall(self._etb, op_neg, [self])
     830
     831    def __abs__(self):
     832        r"""
     833        Compute the absolute value of an Expression.
     834
     835        EXAMPLES:
     836            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     837            sage: etb = ExpressionTreeBuilder(vars=(x,))
     838            sage: x = etb(x)
     839            sage: abs(x)
     840            abs(v_0)
     841            sage: x.abs()
     842            abs(v_0)
     843            sage: x.__abs__()
     844            abs(v_0)
     845        """
     846        return ExpressionCall(self._etb, op_abs, [self])
     847
     848    def abs(self):
     849        r"""
     850        Compute the absolute value of an Expression.
     851
     852        EXAMPLES:
     853            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     854            sage: etb = ExpressionTreeBuilder(vars=(x,))
     855            sage: x = etb(x)
     856            sage: abs(x)
     857            abs(v_0)
     858            sage: x.abs()
     859            abs(v_0)
     860            sage: x.__abs__()
     861            abs(v_0)
     862        """
     863        return ExpressionCall(self._etb, op_abs, [self])
     864
     865    def __invert__(self):
     866        r"""
     867        Compute the inverse of an Expression.
     868
     869        EXAMPLES:
     870            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     871            sage: etb = ExpressionTreeBuilder(vars=(x,))
     872            sage: x = etb(x)
     873            sage: ~x
     874            inv(v_0)
     875            sage: x.__invert__()
     876            inv(v_0)
     877        """
     878        return ExpressionCall(self._etb, op_inv, [self])
     879       
     880
     881cdef class ExpressionConstant(Expression):
     882    r"""
     883    An Expression that represents an arbitrary constant.
     884
     885    EXAMPLES:
     886        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     887        sage: etb = ExpressionTreeBuilder(vars=(x,))
     888        sage: type(etb(3))
     889        <type 'sage.ext.fast_callable.ExpressionConstant'>
     890    """
     891
     892    cdef object _value
     893
     894    def __init__(self, etb, c):
     895        r"""
     896        Initialize an ExpressionConstant.
     897
     898        EXAMPLES:
     899            sage: from sage.ext.fast_callable import ExpressionTreeBuilder, ExpressionConstant
     900            sage: etb = ExpressionTreeBuilder(vars=(x,))
     901            sage: etb(3)
     902            3
     903            sage: v = ExpressionConstant(etb, 3); v
     904            3
     905            sage: v._get_etb() is etb
     906            True
     907            sage: v.value()
     908            3
     909            sage: v.value() == 3
     910            True
     911        """
     912        Expression.__init__(self, etb)
     913        self._value = c
     914
     915    def value(self):
     916        r"""
     917        Return the constant value of an ExpressionConstant.
     918
     919        EXAMPLES:
     920            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     921            sage: etb = ExpressionTreeBuilder(vars=(x,))
     922            sage: etb(3).value()
     923            3
     924        """
     925        return self._value
     926
     927    def __repr__(self):
     928        r"""
     929        Give a string representing this ExpressionConstant.
     930        (We use the repr of its value.)
     931
     932        EXAMPLES:
     933            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     934            sage: etb = ExpressionTreeBuilder(vars=(x,))
     935            sage: v = etb.constant(pi)
     936            sage: v
     937            pi
     938            sage: repr(v)
     939            'pi'
     940            sage: v.__repr__()
     941            'pi'
     942        """
     943        return repr(self._value)
     944
     945cdef class ExpressionVariable(Expression):
     946    r"""
     947    An Expression that represents a variable.
     948
     949    EXAMPLES:
     950        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     951        sage: etb = ExpressionTreeBuilder(vars=(x,))
     952        sage: type(etb.var(x))
     953        <type 'sage.ext.fast_callable.ExpressionVariable'>
     954    """
     955    cdef int _variable_index
     956
     957    def __init__(self, etb, int n):
     958        r"""
     959        Initialize an ExpressionVariable.
     960
     961        EXAMPLES:
     962            sage: from sage.ext.fast_callable import ExpressionTreeBuilder, ExpressionVariable
     963            sage: etb = ExpressionTreeBuilder(vars=(x,))
     964            sage: etb(x)
     965            v_0
     966            sage: v = ExpressionVariable(etb, 0); v
     967            v_0
     968            sage: v._get_etb() is etb
     969            True
     970            sage: v.variable_index()
     971            0
     972        """
     973        Expression.__init__(self, etb)
     974        self._variable_index = n
     975
     976    def variable_index(self):
     977        r"""
     978        Return the variable index of an ExpressionVariable.
     979
     980        EXAMPLES:
     981            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     982            sage: etb = ExpressionTreeBuilder(vars=(x,))
     983            sage: etb(x).variable_index()
     984            0
     985        """
     986        return self._variable_index
     987
     988    def __repr__(self):
     989        r"""
     990        Give a string representing this ExpressionVariable.
     991
     992        EXAMPLES:
     993            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     994            sage: etb = ExpressionTreeBuilder(vars=(x,))
     995            sage: v = etb._var_number(0)
     996            sage: v
     997            v_0
     998            sage: repr(v)
     999            'v_0'
     1000            sage: v.__repr__()       
     1001            'v_0'
     1002        """
     1003        # Should we look up the variable name in self._etb, instead?
     1004        # I think not.. I like the emphasis that we're totally removed
     1005        # from the original expression when we have an Expression.
     1006        return "v_%d" % self._variable_index
     1007
     1008cdef class ExpressionCall(Expression):
     1009    r"""
     1010    An Expression that represents a function call.
     1011
     1012    EXAMPLES:
     1013        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1014        sage: etb = ExpressionTreeBuilder(vars=(x,))
     1015        sage: type(etb.call(sin, x))   
     1016        <type 'sage.ext.fast_callable.ExpressionCall'>
     1017    """
     1018    cdef object _function
     1019    cdef object _arguments
     1020
     1021    def __init__(self, etb, fn, args):
     1022        r"""
     1023        Initialize an ExpressionCall.
     1024
     1025        EXAMPLES:
     1026            sage: from sage.ext.fast_callable import ExpressionTreeBuilder, ExpressionCall
     1027            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1028            sage: x = etb(x)
     1029            sage: etb.call(factorial, x)
     1030            {factorial}(v_0)
     1031            sage: v = ExpressionCall(etb, factorial, [x]); v
     1032            {factorial}(v_0)
     1033            sage: v._get_etb() is etb
     1034            True
     1035            sage: v.function()
     1036            factorial
     1037            sage: v.arguments()
     1038            [v_0]
     1039        """
     1040        Expression.__init__(self, etb)
     1041        self._function = fn
     1042        self._arguments = args
     1043
     1044    def function(self):
     1045        r"""
     1046        Return the function from this ExpressionCall.
     1047
     1048        EXAMPLES:
     1049            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1050            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1051            sage: etb.call(sin, x).function()
     1052            sin
     1053        """
     1054        return self._function
     1055
     1056    def arguments(self):
     1057        r"""
     1058        Return the arguments from this ExpressionCall.
     1059
     1060        EXAMPLES:
     1061            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1062            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1063            sage: etb.call(sin, x).arguments()
     1064            [v_0]
     1065        """
     1066        return copy(self._arguments)
     1067
     1068    def __repr__(self):
     1069        r"""
     1070        Give a string representing this ExpressionCall.
     1071
     1072        EXAMPLES:
     1073            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1074            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1075            sage: x = etb.var(x)
     1076            sage: etb.call(operator.add, x, 1)
     1077            add(v_0, 1)
     1078            sage: etb.call(factorial, x)
     1079            {factorial}(v_0)
     1080            sage: v = etb.call(sin, x)
     1081            sage: v
     1082            sin(v_0)
     1083            sage: repr(v)
     1084            'sin(v_0)'
     1085            sage: v.__repr__()
     1086            'sin(v_0)'
     1087        """
     1088        fn = function_name(self._function)
     1089        return '%s(%s)' % (fn, ', '.join(map(repr, self._arguments)))
     1090
     1091cdef class ExpressionChoice(Expression):
     1092    r"""
     1093    A conditional expression.
     1094
     1095    (It's possible to create choice nodes, but they don't work yet.)
     1096   
     1097    EXAMPLES:
     1098        sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1099        sage: etb = ExpressionTreeBuilder(vars=(x,))
     1100        sage: etb.choice(etb.call(operator.eq, x, 0), 0, 1/x)
     1101        (0 if {eq}(v_0, 0) else div(1, v_0))
     1102    """
     1103
     1104    cdef object _cond
     1105    cdef object _iftrue
     1106    cdef object _iffalse
     1107
     1108    def __init__(self, etb, cond, iftrue, iffalse):
     1109        r"""
     1110        Initialize an ExpressionChoice.
     1111
     1112        EXAMPLES:
     1113            sage: from sage.ext.fast_callable import ExpressionTreeBuilder, ExpressionChoice
     1114            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1115            sage: x = etb(x)
     1116            sage: etb.choice(x, ~x, 0)
     1117            (inv(v_0) if v_0 else 0)
     1118            sage: v = ExpressionChoice(etb, x, ~x, etb(0)); v
     1119            (inv(v_0) if v_0 else 0)
     1120            sage: v._get_etb() is etb
     1121            True
     1122            sage: v.condition()
     1123            v_0
     1124            sage: v.if_true()
     1125            inv(v_0)
     1126            sage: v.if_false()
     1127            0
     1128        """
     1129        Expression.__init__(self, etb)
     1130        self._cond = cond
     1131        self._iftrue = iftrue
     1132        self._iffalse = iffalse
     1133
     1134    def condition(self):
     1135        r"""
     1136        Return the condition of an ExpressionChoice.
     1137
     1138        EXAMPLES:
     1139            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1140            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1141            sage: x = etb(x)
     1142            sage: etb.choice(x, ~x, 0).condition()
     1143            v_0
     1144        """
     1145        return self._cond
     1146
     1147    def if_true(self):
     1148        r"""
     1149        Return the true branch of an ExpressionChoice.
     1150
     1151        EXAMPLES:
     1152            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1153            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1154            sage: x = etb(x)
     1155            sage: etb.choice(x, ~x, 0).if_true()
     1156            inv(v_0)
     1157        """
     1158        return self._iftrue
     1159
     1160    def if_false(self):
     1161        r"""
     1162        Return the false branch of an ExpressionChoice.
     1163
     1164        EXAMPLES:
     1165            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1166            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1167            sage: x = etb(x)
     1168            sage: etb.choice(x, ~x, 0).if_false()
     1169            0
     1170        """
     1171        return self._iffalse
     1172
     1173    def __repr__(self):
     1174        r"""
     1175        Give a string representation for this ExpressionChoice.
     1176        (Based on the syntax for Python conditional expressions.)
     1177
     1178        EXAMPLES:
     1179            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     1180            sage: etb = ExpressionTreeBuilder(vars=(x,))
     1181            sage: x = etb(x)
     1182            sage: v = etb.choice(x, ~x, 0)
     1183            sage: v
     1184            (inv(v_0) if v_0 else 0)
     1185            sage: repr(v)
     1186            '(inv(v_0) if v_0 else 0)'
     1187            sage: v.__repr__()
     1188            '(inv(v_0) if v_0 else 0)'
     1189        """
     1190        return '(%s if %s else %s)' % (repr(self._iftrue),
     1191                                       repr(self._cond),
     1192                                       repr(self._iffalse))
     1193
     1194cpdef _expression_binop_helper(s, o, op):
     1195   r"""
     1196   Makes an Expression for (s op o).  Either s or o (or both) must already
     1197   be an expression.   
     1198
     1199   EXAMPLES:
     1200       sage: from sage.ext.fast_callable import _expression_binop_helper, ExpressionTreeBuilder
     1201       sage: var('x,y')
     1202       (x, y)
     1203       sage: etb = ExpressionTreeBuilder(vars=(x,y))
     1204       sage: x = etb(x)
     1205
     1206   Now x is an Expression, but y is not.  Still, all the following
     1207   cases work.
     1208       sage: _expression_binop_helper(x, x, operator.add)
     1209       add(v_0, v_0)
     1210       sage: _expression_binop_helper(x, y, operator.add)
     1211       add(v_0, v_1)
     1212       sage: _expression_binop_helper(y, x, operator.add)
     1213       add(v_1, v_0)
     1214   """
     1215   # The Cython way of handling operator overloading on cdef classes
     1216   # (which is inherited from Python) is quite annoying.  Inside the
     1217   # code for a binary operator, you know that either the first or
     1218   # second argument (or both) is a member of your class, but you
     1219   # don't know which.
     1220
     1221   # If there is an arithmetic operator between an Expression and
     1222   # a non-Expression, I want to convert the non-Expression into
     1223   # an Expression.  But to do that, I need the ExpressionTreeBuilder
     1224   # from the Expression.
     1225
     1226   cdef Expression self
     1227   cdef Expression other
     1228
     1229   if not isinstance(o, Expression):
     1230       self = s
     1231       other = self._etb(o)
     1232   elif not isinstance(s, Expression):
     1233       other = o
     1234       self = other._etb(s)
     1235   else:
     1236       self = s
     1237       other = o
     1238       assert self._etb is other._etb
     1239
     1240   return ExpressionCall(self._etb, op, [self, other])
     1241
     1242builtin_functions = None
     1243cpdef get_builtin_functions():
     1244    r"""
     1245    To handle ExpressionCall, we need to map from Sage and
     1246    Python functions to opcode names.
     1247
     1248    This returns a dictionary which is that map.
     1249
     1250    We delay building builtin_functions to break a circular import
     1251    between sage.calculus and this file.
     1252
     1253    EXAMPLES:
     1254        sage: from sage.ext.fast_callable import get_builtin_functions
     1255        sage: builtins = get_builtin_functions()
     1256        sage: sorted(list(builtins.values()))
     1257        ['abs', 'abs', 'acos', 'acosh', 'add', 'asin', 'asinh', 'atan', 'atanh', 'ceil', 'cos', 'cosh', 'cot', 'csc', 'div', 'exp', 'floor', 'floordiv', 'inv', 'log', 'mul', 'neg', 'pow', 'sec', 'sin', 'sinh', 'sqrt', 'sub', 'tan', 'tanh']
     1258        sage: builtins[sin]
     1259        'sin'
     1260    """
     1261    # We delay building builtin_functions to break a circular import
     1262    # between sage.calculus and this file.
     1263    global builtin_functions
     1264    if builtin_functions is not None:
     1265        return builtin_functions
     1266    builtin_functions = {
     1267        operator.add: 'add',
     1268        operator.sub: 'sub',
     1269        operator.mul: 'mul',
     1270        operator.div: 'div',
     1271        operator.floordiv: 'floordiv',
     1272        operator.abs: 'abs',
     1273        operator.neg: 'neg',
     1274        operator.inv: 'inv',
     1275        operator.pow: 'pow',
     1276        }
     1277    # not handled: atan2, log2, log10
     1278    import sage.calculus.all as calc_all
     1279    for fn in ('sqrt', 'ceil', 'floor',
     1280               'sin', 'cos', 'tan', 'sec', 'csc', 'cot',
     1281               'asin', 'acos', 'atan', 'sinh', 'cosh', 'tanh',
     1282               'asinh', 'acosh', 'atanh', 'exp', 'log'):
     1283        builtin_functions[getattr(calc_all, fn)] = fn
     1284    builtin_functions[calc_all.abs_symbolic] = 'abs'
     1285    return builtin_functions
     1286
     1287cpdef generate_code(Expression expr, stream):
     1288    r"""
     1289    Generate code from an Expression tree; write the result into an
     1290    InstructionStream.
     1291
     1292    In fast_callable, first we create an Expression, either directly
     1293    with an ExpressionTreeBuilder or with _fast_callable_ methods.
     1294    Then we optimize the Expression in tree form.  (Unfortunately,
     1295    this step is currently missing -- we do no optimizations.)
     1296
     1297    Then we linearize the Expression into a sequence of instructions,
     1298    by walking the Expression and sending the corresponding stack
     1299    instructions to an InstructionStream.
     1300
     1301    EXAMPLES:
     1302        sage: from sage.ext.fast_callable import ExpressionTreeBuilder, generate_code, InstructionStream
     1303        sage: etb = ExpressionTreeBuilder('x')
     1304        sage: x = etb.var('x')
     1305        sage: expr = ((x+pi)*(x+1))
     1306        sage: from sage.ext.interpreters.wrapper_py import metadata, Wrapper_py
     1307        sage: instr_stream = InstructionStream(metadata, 1)
     1308        sage: generate_code(expr, instr_stream)
     1309        sage: instr_stream.instr('return')
     1310        sage: v = Wrapper_py(instr_stream.get_current())
     1311        sage: type(v)
     1312        <type 'sage.ext.interpreters.wrapper_py.Wrapper_py'>
     1313        sage: v(7)
     1314        8*(pi + 7)
     1315
     1316    TESTS:
     1317        sage: def my_sin(x): return sin(x)
     1318        sage: def my_norm(x, y): return x*x + y*y
     1319        sage: def my_sqrt(x): return sqrt(x, extend=False)
     1320        sage: fc = fast_callable(expr, domain=RealField(130))
     1321        sage: fc(0)
     1322        3.1415926535897932384626433832795028842
     1323        sage: fc(1)
     1324        8.2831853071795864769252867665590057684
     1325        sage: fc = fast_callable(expr, domain=RDF)
     1326        sage: fc(0)
     1327        3.1415926535897931
     1328        sage: fc(1)
     1329        8.2831853071795862
     1330        sage: fc.op_list()
     1331        [('load_arg', 0), ('load_const', pi), 'add', ('load_arg', 0), ('load_const', 1), 'add', 'mul', 'return']
     1332        sage: fc = fast_callable(etb.call(sin, x) + etb.call(sqrt, x), domain=RDF)
     1333        sage: fc(1)
     1334        1.8414709848078965
     1335        sage: fc.op_list()
     1336        [('load_arg', 0), 'sin', ('load_arg', 0), 'sqrt', 'add', 'return']
     1337        sage: fc = fast_callable(etb.call(sin, x) + etb.call(sqrt, x))
     1338        sage: fc(1)
     1339        sin(1) + 1
     1340        sage: fc.op_list()
     1341        [('load_arg', 0), ('py_call', sin, 1), ('load_arg', 0), ('py_call', sqrt, 1), 'add', 'return']
     1342        sage: fc = fast_callable(etb.call(my_sin, x), domain=RDF)
     1343        sage: fc(3)
     1344        0.14112000805986721
     1345        sage: fc = fast_callable(etb.call(my_sin, x), domain=RealField(100))
     1346        sage: fc(3)
     1347        0.14112000805986722210074480281
     1348        sage: fc.op_list()
     1349        [('load_arg', 0), ('py_call', <function my_sin at 0x...>, 1), 'return']
     1350        sage: fc = fast_callable(etb.call(my_sqrt, x), domain=RDF)
     1351        sage: fc(3)
     1352        1.7320508075688772
     1353        sage: fc(-3)
     1354        Traceback (most recent call last):
     1355        ...
     1356        ValueError: math domain error
     1357        sage: fc = fast_callable(etb.call(my_sqrt, x), domain=RR)
     1358        sage: fc(3)
     1359        1.73205080756888
     1360        sage: fc(-3)
     1361        Traceback (most recent call last):
     1362        ...
     1363        ValueError: negative number -3.00000000000000 does not have a square root in the real field
     1364        sage: etb2 = ExpressionTreeBuilder(('y','z'))
     1365        sage: y = etb2.var('y')
     1366        sage: z = etb2.var('z')
     1367        sage: fc = fast_callable(etb2.call(sqrt, etb2.call(my_norm, y, z)), domain=RDF)
     1368        sage: fc(3, 4)
     1369        5.0
     1370        sage: fc.op_list()
     1371        [('load_arg', 0), ('load_arg', 1), ('py_call', <function my_norm at 0x...>, 2), 'sqrt', 'return']
     1372        sage: fc.python_calls()
     1373        [<function my_norm at 0x...>]
     1374        sage: fc = fast_callable(etb2.call(sqrt, etb2.call(my_norm, y, z)), domain=RR)
     1375        sage: fc(3, 4)
     1376        5.00000000000000
     1377        sage: fc = fast_callable(etb2.call(my_norm, y, z), domain=ZZ)
     1378        sage: fc(3, 4)
     1379        25
     1380        sage: fc.op_list()
     1381        [('load_arg', 0), ('load_arg', 1), ('py_call', <function my_norm at 0x...>, 2), 'return']
     1382        sage: fc = fast_callable(expr)
     1383        sage: fc(3.0r)
     1384        4.0*(pi + 3.0)
     1385        sage: fc = fast_callable(x+3, domain=ZZ)
     1386        sage: fc(4)
     1387        7
     1388        sage: fc = fast_callable(x/3, domain=ZZ)
     1389        sage: fc(4)
     1390        Traceback (most recent call last):
     1391        ...
     1392        TypeError: no conversion of this rational to integer
     1393        sage: fc(6)
     1394        2
     1395        sage: fc = fast_callable(etb.call(sin, x), domain=ZZ)
     1396        sage: fc(0)
     1397        0
     1398        sage: fc(3)
     1399        Traceback (most recent call last):
     1400        ...
     1401        TypeError: unable to convert x (=sin(3)) to an integer
     1402    """
     1403    cdef ExpressionConstant econst
     1404    cdef ExpressionVariable evar
     1405    cdef ExpressionCall ecall
     1406    cdef ExpressionChoice echoice
     1407
     1408    if isinstance(expr, ExpressionConstant):
     1409        econst = expr
     1410        stream.load_const(econst._value)
     1411    elif isinstance(expr, ExpressionVariable):
     1412        evar = expr
     1413        stream.load_arg(evar._variable_index)
     1414    elif isinstance(expr, ExpressionCall):
     1415        ecall = expr
     1416        fn = ecall._function
     1417        for arg in ecall._arguments:
     1418            generate_code(arg, stream)
     1419        builtins = get_builtin_functions()
     1420        if fn in builtins:
     1421            opname = builtins[fn]
     1422            if stream.has_instr(opname):
     1423                stream.instr(builtin_functions[fn])
     1424                return
     1425        if stream.has_instr('py_call'):
     1426            stream.instr('py_call', fn, len(ecall._arguments))
     1427        else:
     1428            raise ValueError, "Unhandled function %s in generate_code" % fn
     1429    else:
     1430        raise ValueError, "Unhandled expression kind %s in generate_code" % type(expr)
     1431
     1432
     1433cdef class InstructionStream:
     1434    r"""
     1435    An InstructionStream takes a sequence of instructions (passed in by
     1436    a series of method calls) and computes the data structures needed
     1437    by the interpreter.  This is the stage where we switch from operating
     1438    on Expression trees to a linear representation.  If we had a peephole
     1439    optimizer (we don't) it would go here.
     1440
     1441    Currently, this class is not very general; it only works for
     1442    interpreters with a fixed set of memory chunks (with fixed names).
     1443    Basically, it only works for stack-based expression interpreters.
     1444    It should be generalized, so that the interpreter metadata includes
     1445    a description of the memory chunks involved and the instruction stream
     1446    can handle any interpreter.
     1447
     1448    Once you're done adding instructions, you call get_current() to retrieve
     1449    the information needed by the interpreter (as a Python dictionary).
     1450    """
     1451
     1452    cdef object _metadata
     1453    cdef object _instrs
     1454    cdef object _bytecode
     1455    cdef object _constants
     1456    cdef object _constant_locs
     1457    cdef object _py_constants
     1458    cdef object _py_constant_locs
     1459    cdef int _stack_cur_size
     1460    cdef int _stack_max_size
     1461    cdef int _n_args
     1462    cdef object _domain
     1463
     1464    def __init__(self, metadata, n_args, domain=None):
     1465        r"""
     1466        Initialize an InstructionStream.
     1467
     1468        INPUTS:
     1469            metadata - The metadata_by_opname from a wrapper module
     1470            n_args - The number of arguments accessible by the generated code
     1471                     (this is just passed to the wrapper class)
     1472            domain - The domain of interpretation (this is just passed to the
     1473                     wrapper class)
     1474
     1475        EXAMPLES:
     1476            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1477            sage: from sage.ext.fast_callable import InstructionStream
     1478            sage: instr_stream = InstructionStream(metadata, 1)
     1479            sage: instr_stream.get_current()
     1480            {'domain': None, 'code': [], 'py_constants': [], 'args': 1, 'stack': 0, 'constants': []}
     1481            sage: md = instr_stream.get_metadata()
     1482            sage: type(md)
     1483            <class 'sage.ext.fast_callable.InterpreterMetadata'>
     1484            sage: md.by_opname['py_call']
     1485            (CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs']), 3)
     1486            sage: md.by_opcode[3]
     1487            ('py_call', CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs']))
     1488        """
     1489        self._metadata = metadata
     1490        self._instrs = []
     1491        self._bytecode = []
     1492        self._constants = []
     1493        self._constant_locs = {}
     1494        self._py_constants = []
     1495        self._py_constant_locs = {}
     1496        self._stack_cur_size = 0
     1497        self._stack_max_size = 0
     1498        self._domain = domain
     1499        self._n_args = n_args
     1500
     1501    def load_const(self, c):
     1502        r"""
     1503        Add a 'load_const' instruction to this InstructionStream.
     1504
     1505        EXAMPLES:
     1506            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1507            sage: from sage.ext.fast_callable import InstructionStream, op_list
     1508            sage: instr_stream = InstructionStream(metadata, 1)
     1509            sage: instr_stream.load_const(5)
     1510            sage: instr_stream.current_op_list()
     1511            [('load_const', 5)]
     1512            sage: instr_stream.load_const(7)
     1513            sage: instr_stream.load_const(5)
     1514            sage: instr_stream.current_op_list()
     1515            [('load_const', 5), ('load_const', 7), ('load_const', 5)]
     1516
     1517        Note that constants are shared: even though we load 5 twice, it
     1518        only appears once in the constant table.
     1519            sage: instr_stream.get_current()['constants']           
     1520            [5, 7]
     1521        """
     1522        self.instr('load_const', c)
     1523
     1524    def load_arg(self, n):
     1525        r"""
     1526        Add a 'load_arg' instruction to this InstructionStream.
     1527
     1528        EXAMPLES:
     1529            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1530            sage: from sage.ext.fast_callable import InstructionStream
     1531            sage: instr_stream = InstructionStream(metadata, 12)
     1532            sage: instr_stream.load_arg(5)
     1533            sage: instr_stream.current_op_list()
     1534            [('load_arg', 5)]
     1535            sage: instr_stream.load_arg(3)
     1536            sage: instr_stream.current_op_list()
     1537            [('load_arg', 5), ('load_arg', 3)]
     1538        """
     1539        self.instr('load_arg', n)
     1540
     1541    def has_instr(self, opname):
     1542        r"""
     1543        Check whether this InstructionStream knows how to generate code
     1544        for a given instruction.
     1545
     1546        EXAMPLES:
     1547            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1548            sage: from sage.ext.fast_callable import InstructionStream
     1549            sage: instr_stream = InstructionStream(metadata, 1)
     1550            sage: instr_stream.has_instr('return')
     1551            True
     1552            sage: instr_stream.has_instr('factorial')
     1553            False
     1554            sage: instr_stream.has_instr('abs')
     1555            True
     1556        """
     1557        return (opname in self._metadata.by_opname)
     1558
     1559    def instr(self, opname, *args):
     1560        r"""
     1561        Generate code in this InstructionStream for the given instruction
     1562        and arguments.
     1563
     1564        The opname is used to look up a CompilerInstrSpec; the
     1565        CompilerInstrSpec describes how to interpret the arguments.
     1566        (This is documented in the class docstring for CompilerInstrSpec.)
     1567
     1568        EXAMPLES:
     1569            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1570            sage: from sage.ext.fast_callable import InstructionStream
     1571            sage: instr_stream = InstructionStream(metadata, 1)
     1572            sage: instr_stream.instr('load_arg', 0)
     1573            sage: instr_stream.instr('sin')
     1574            sage: instr_stream.instr('py_call', math.sin, 1)
     1575            sage: instr_stream.instr('abs')
     1576            sage: instr_stream.instr('factorial')
     1577            Traceback (most recent call last):
     1578            ...
     1579            KeyError: 'factorial'
     1580            sage: instr_stream.instr('return')
     1581            sage: instr_stream.current_op_list()
     1582            [('load_arg', 0), 'sin', ('py_call', <built-in function sin>, 1), 'abs', 'return']
     1583        """
     1584        cdef int i
     1585
     1586        spec, opcode = self._metadata.by_opname[opname]
     1587        assert len(spec.parameters) == len(args)
     1588
     1589        n_inputs = spec.n_inputs
     1590        n_outputs = spec.n_outputs
     1591
     1592        self._bytecode.append(opcode)
     1593        for i in range(len(args)):
     1594            if spec.parameters[i] == 'constants':
     1595                # XXX bad for strict-mode floating-point constants
     1596                # (doesn't handle signed 0, NaN)
     1597                arg = args[i]
     1598                if arg in self._constant_locs:
     1599                    self._bytecode.append(self._constant_locs[arg])
     1600                else:
     1601                    loc = len(self._constants)
     1602                    self._constants.append(arg)
     1603                    self._constant_locs[arg] = loc
     1604                    self._bytecode.append(loc)
     1605            elif spec.parameters[i] == 'args':
     1606                self._bytecode.append(args[i])
     1607            elif spec.parameters[i] == 'n_inputs':
     1608                self._bytecode.append(args[i])
     1609                n_inputs = args[i]
     1610            elif spec.parameters[i] == 'n_outputs':
     1611                self._bytecode.append(args[i])
     1612                n_outputs = args[i]
     1613            elif spec.parameters[i] == 'py_constants':
     1614                arg = args[i]
     1615                if arg in self._py_constant_locs:
     1616                    self._bytecode.append(self._py_constant_locs[arg])
     1617                else:
     1618                    loc = len(self._py_constants)
     1619                    self._py_constants.append(arg)
     1620                    self._py_constant_locs[arg] = loc
     1621                    self._bytecode.append(loc)
     1622            else:
     1623                raise ValueError
     1624
     1625        self._stack_cur_size -= n_inputs
     1626        self._stack_cur_size += n_outputs
     1627        self._stack_max_size = max(self._stack_max_size, self._stack_cur_size)
     1628
     1629    def get_metadata(self):
     1630        r"""
     1631        Returns the interpreter metadata being used by the current
     1632        InstructionStream.
     1633
     1634        (Probably only useful for writing doctests.)
     1635
     1636        EXAMPLES:
     1637            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1638            sage: from sage.ext.fast_callable import InstructionStream
     1639            sage: instr_stream = InstructionStream(metadata, 1)
     1640            sage: md = instr_stream.get_metadata()
     1641            sage: type(md)
     1642            <class 'sage.ext.fast_callable.InterpreterMetadata'>
     1643        """
     1644        return self._metadata
     1645
     1646    def current_op_list(self):
     1647        r"""
     1648        Returns the list of intructions that have been added to this
     1649        InstructionStream so far.
     1650
     1651        It's OK to call this, then add more instructions.
     1652
     1653        EXAMPLES:
     1654            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1655            sage: from sage.ext.fast_callable import InstructionStream
     1656            sage: instr_stream = InstructionStream(metadata, 1)
     1657            sage: instr_stream.instr('load_arg', 0)
     1658            sage: instr_stream.instr('py_call', math.sin, 1)
     1659            sage: instr_stream.instr('abs')
     1660            sage: instr_stream.instr('return')
     1661            sage: instr_stream.current_op_list()
     1662            [('load_arg', 0), ('py_call', <built-in function sin>, 1), 'abs', 'return']
     1663        """
     1664        return op_list(self.get_current(), self._metadata)
     1665
     1666    def get_current(self):
     1667        r"""
     1668        Return the current state of the InstructionStream, as a dictionary
     1669        suitable for passing to a wrapper class.
     1670
     1671        NOTE: The dictionary includes internal data structures of the
     1672        InstructionStream; you must not modify it.
     1673
     1674        EXAMPLES:
     1675            sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1676            sage: from sage.ext.fast_callable import InstructionStream
     1677            sage: instr_stream = InstructionStream(metadata, 1)
     1678            sage: instr_stream.get_current()
     1679            {'domain': None, 'code': [], 'py_constants': [], 'args': 1, 'stack': 0, 'constants': []}
     1680            sage: instr_stream.instr('load_arg', 0)
     1681            sage: instr_stream.instr('py_call', math.sin, 1)
     1682            sage: instr_stream.instr('abs')
     1683            sage: instr_stream.instr('return')
     1684            sage: instr_stream.current_op_list()
     1685            [('load_arg', 0), ('py_call', <built-in function sin>, 1), 'abs', 'return']
     1686            sage: instr_stream.get_current()
     1687            {'domain': None, 'code': [0, 0, 3, 0, 1, 11, 2], 'py_constants': [<built-in function sin>], 'args': 1, 'stack': 1, 'constants': []}
     1688        """
     1689        d = {'args': self._n_args,
     1690             'constants': self._constants,
     1691             'py_constants': self._py_constants,
     1692             'stack': self._stack_max_size,
     1693             'code': self._bytecode,
     1694             'domain': self._domain}
     1695        return d
     1696
     1697class InterpreterMetadata(object):
     1698    r"""
     1699    The interpreter metadata for a fast_callable interpreter.  Currently
     1700    only consists of a dictionary mapping instruction names to
     1701    (CompilerInstrSpec, opcode) pairs, and a list mapping opcodes to
     1702    (instruction name, CompilerInstrSpec) pairs.
     1703
     1704    See the class docstring for CompilerInstrSpec for more information.
     1705
     1706    NOTE: You must not modify the metadata.
     1707    """
     1708
     1709    def __init__(self, by_opname, by_opcode):
     1710        r"""
     1711        Initialize an InterpreterMetadata object.
     1712
     1713        EXAMPLES:
     1714            sage: from sage.ext.fast_callable import InterpreterMetadata
     1715
     1716        Currently we do no error checking or processing, so we can
     1717        use this simple test:
     1718            sage: metadata = InterpreterMetadata(by_opname='opname dict goes here', by_opcode='opcode list goes here')
     1719            sage: metadata.by_opname
     1720            'opname dict goes here'
     1721            sage: metadata.by_opcode
     1722            'opcode list goes here'
     1723        """
     1724        self.by_opname = by_opname
     1725        self.by_opcode = by_opcode
     1726
     1727class CompilerInstrSpec(object):
     1728    r"""
     1729    Describes a single instruction to the fast_callable code generator.
     1730
     1731    An instruction has a number of stack inputs, a number of stack
     1732    outputs, and a parameter list describing extra arguments that
     1733    must be passed to the InstructionStream.instr method (that end up
     1734    as extra words in the code).
     1735
     1736    The parameter list is a list of strings.  Each string is one of
     1737    the following:
     1738
     1739        - 'args' - The instruction argument refers to an input argument
     1740                   of the wrapper class; it is just appended to the code.
     1741        - 'constants', 'py_constants' - The instruction argument is a value;
     1742                   the value is added to the corresponding list (if it's
     1743                   not already there) and the index is appended to the
     1744                   code.
     1745        - 'n_inputs', 'n_outputs' - The instruction actually takes a variable
     1746                   number of inputs or outputs (the n_inputs and n_outputs
     1747                   attributes of this instruction are ignored).
     1748                   The instruction argument specifies the number of inputs
     1749                   or outputs (respectively); it is just appended to the code.
     1750    """
     1751
     1752    def __init__(self, n_inputs, n_outputs, parameters):
     1753        r"""
     1754        Initialize a CompilerInstrSpec.
     1755
     1756        EXAMPLES:
     1757            sage: from sage.ext.fast_callable import CompilerInstrSpec
     1758            sage: CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])
     1759            CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])
     1760        """
     1761        self.n_inputs = n_inputs
     1762        self.n_outputs = n_outputs
     1763        self.parameters = parameters
     1764
     1765    def __repr__(self):
     1766        r"""
     1767        Give a string representation for this CompilerInstrSpec.
     1768
     1769        EXAMPLES:
     1770            sage: from sage.ext.fast_callable import CompilerInstrSpec
     1771            sage: v = CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])
     1772            sage: v
     1773            CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])
     1774            sage: repr(v)
     1775            "CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])"
     1776            sage: v.__repr__()
     1777            "CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])"
     1778        """
     1779        return "CompilerInstrSpec(%d, %d, %s)" % (self.n_inputs, self.n_outputs, self.parameters)
     1780
     1781def op_list(args, metadata):
     1782    r"""
     1783    Given a dictionary with the result of calling get_current on an
     1784    InstructionStream, and the corresponding interpreter metadata,
     1785    return a list of the instructions, in a simple somewhat
     1786    human-readable format.
     1787
     1788    For debugging only.  (That is, it's probably not a good idea to
     1789    try to programmatically manipulate the result of this function;
     1790    the expected use is just to print the returned list to the
     1791    screen.)
     1792
     1793    There's probably no reason to call this directly; if you
     1794    have a wrapper object, call op_list on it; if you have an
     1795    InstructionStream object, call current_op_list on it.
     1796
     1797    EXAMPLES:
     1798        sage: from sage.ext.interpreters.wrapper_rdf import metadata
     1799        sage: from sage.ext.fast_callable import InstructionStream, op_list
     1800        sage: instr_stream = InstructionStream(metadata, 1)
     1801        sage: instr_stream.instr('load_arg', 0)
     1802        sage: instr_stream.instr('abs')
     1803        sage: instr_stream.instr('return')
     1804        sage: instr_stream.current_op_list()
     1805        [('load_arg', 0), 'abs', 'return']
     1806        sage: op_list(instr_stream.get_current(), metadata)
     1807        [('load_arg', 0), 'abs', 'return']
     1808    """
     1809    ops = []
     1810    code = args['code']
     1811    while len(code):
     1812        opcode = code[0]
     1813        code = code[1:]
     1814        (opname, instr) = metadata.by_opcode[opcode]
     1815        if len(instr.parameters):
     1816            op = [opname]
     1817            for p in instr.parameters:
     1818                p_loc = code[0]
     1819                code = code[1:]
     1820                if p in ('args', 'n_inputs', 'n_outputs'):
     1821                    op.append(p_loc)
     1822                else:
     1823                    op.append(args[p][p_loc])
     1824            ops.append(tuple(op))
     1825        else:
     1826            ops.append(opname)
     1827    return ops
     1828   
     1829
     1830cdef class Wrapper:
     1831    r"""
     1832    The parent class for all fast_callable wrappers.  Implements shared
     1833    behavior (currently only debugging).
     1834    """
     1835
     1836    def __init__(self, args, metadata):
     1837        r"""
     1838        Initialize a Wrapper object.
     1839
     1840        EXAMPLES:
     1841            sage: from sage.ext.fast_callable import ExpressionTreeBuilder, generate_code, InstructionStream
     1842            sage: etb = ExpressionTreeBuilder('x')
     1843            sage: x = etb.var('x')
     1844            sage: expr = ((x+pi)*(x+1))
     1845            sage: from sage.ext.interpreters.wrapper_py import metadata, Wrapper_py
     1846            sage: instr_stream = InstructionStream(metadata, 1)
     1847            sage: generate_code(expr, instr_stream)
     1848            sage: instr_stream.instr('return')
     1849            sage: v = Wrapper_py(instr_stream.get_current())
     1850            sage: v.get_orig_args()
     1851            {'domain': None, 'code': [0, 0, 1, 0, 4, 0, 0, 1, 1, 4, 6, 2], 'py_constants': [], 'args': 1, 'stack': 3, 'constants': [pi, 1]}
     1852            sage: v.op_list()
     1853            [('load_arg', 0), ('load_const', pi), 'add', ('load_arg', 0), ('load_const', 1), 'add', 'mul', 'return']
     1854        """
     1855
     1856        # We only keep the original arguments for debugging (op_list(), etc.);
     1857        # is it worth the memory cost?  (Note that we may be holding on to
     1858        # large objects that could otherwise be garbage collected, for
     1859        # instance.)
     1860        self._orig_args = args
     1861        self._metadata = metadata
     1862
     1863    def get_orig_args(self):
     1864        r"""
     1865        Get the original arguments used when initializing this
     1866        wrapper.
     1867
     1868        (Probably only useful when writing doctests.)
     1869
     1870        EXAMPLES:
     1871            sage: fast_callable(sin(x)/x, domain=RDF).get_orig_args()
     1872            {'domain': None, 'code': [0, 0, 15, 0, 0, 7, 2], 'py_constants': [], 'args': 1, 'stack': 2, 'constants': []}
     1873        """
     1874        return self._orig_args
     1875
     1876    def op_list(self):
     1877        r"""
     1878        Return the list of instructions in this wrapper.
     1879
     1880        EXAMPLES:
     1881            sage: fast_callable(cos(x)*x, domain=RDF).op_list()
     1882            [('load_arg', 0), 'cos', ('load_arg', 0), 'mul', 'return']
     1883        """
     1884        return op_list(self._orig_args, self._metadata)
     1885                       
     1886    def python_calls(self):
     1887        r"""
     1888        List the Python functions that are called in this wrapper.
     1889
     1890        (Python function calls are slow, so ideally this list would
     1891        be empty.  If it is not empty, then perhaps there is an
     1892        optimization opportunity where a Sage developer could speed
     1893        this up by adding a new instruction to the interpreter.)
     1894
     1895        EXAMPLES:
     1896            sage: fast_callable(abs(sin(x)), domain=RDF).python_calls()
     1897            []
     1898            sage: fast_callable(abs(sin(factorial(x)))).python_calls()
     1899            [factorial, sin]
     1900        """
     1901        ops = self.op_list()
     1902        py_calls = []
     1903        for op in ops:
     1904            if isinstance(op, tuple) and op[0] == 'py_call':
     1905                py_calls.append(op[1])
     1906        return py_calls
     1907       
  • sage/ext/fast_eval.pyx

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/ext/fast_eval.pyx
    a b  
    2424See the function \code{fast_float(f, *vars)} to create a fast-callable
    2525version of f.
    2626
     27NOTE: Sage temporarily has two implementations of this functionality;
     28one in this file, which will probably be deprecated soon, and one in
     29fast_callable.pyx.  The following instructions are for the old
     30implementation; you probably want to be looking at fast_callable.pyx
     31instead.
     32
    2733To provide this interface for a class, implement
    2834\code{_fast_float_(self, *vars)}.  The basic building blocks are
    2935provided by the functions \code{fast_float_constant} (returns a
     
    3541
    3642EXAMPLES:
    3743    sage: from sage.ext.fast_eval import fast_float
    38     sage: f = fast_float(sqrt(x^7+1), 'x')
     44    sage: f = fast_float(sqrt(x^7+1), 'x', old=True)
    3945    sage: f(1)
    4046    1.4142135623730951
    4147    sage: f.op_list()
     
    8187#                  http://www.gnu.org/licenses/
    8288#*****************************************************************************
    8389
     90from sage.ext.fast_callable import fast_callable, Wrapper
     91from sage.rings.all import RDF
     92
    8493include "stdsage.pxi"
    8594
    8695cdef extern from "Python.h":
     
    12881297    return FastDoubleFunc('callable', f, *args)
    12891298
    12901299
    1291 def fast_float(f, *vars):
     1300new_fast_float=True
     1301
     1302def fast_float(f, *vars, old=None):
    12921303    """
    1293     Tries to create a fast float function out of the
    1294     input, if possible.
     1304    Tries to create a function that evaluates f quickly using
     1305    floating-point numbers, if possible.  There are two implementations
     1306    of fast_float in Sage; by default we use the newer, which is
     1307    slightly faster on most tests.
    12951308   
    12961309    On failure, returns the input unchanged.
    12971310   
    12981311    INPUT:
    12991312        f    -- an expression
    13001313        vars -- the names of the arguments
     1314        old  -- use the original algorithm for fast_float
    13011315       
    13021316    EXAMPLES:
    13031317        sage: from sage.ext.fast_eval import fast_float
     
    13151329        sage: f(1,2)
    13161330        1.0
    13171331    """
     1332    if old is None:
     1333        old = not new_fast_float
     1334
    13181335    if isinstance(f, (tuple, list)):
    13191336        return tuple([fast_float(x, *vars) for x in f])
    13201337   
     
    13281345            vars = vars[:i] + (v,) + vars[i+1:]
    13291346
    13301347    try:
    1331         return f._fast_float_(*vars)
     1348        if old:
     1349            return f._fast_float_(*vars)
     1350        else:
     1351            return fast_callable(f, vars=vars, domain=RDF)
    13321352    except AttributeError:
    13331353        pass
    13341354
     
    13391359
    13401360    try:
    13411361        from sage.calculus.calculus import SR
    1342         return SR(f)._fast_float_(*vars)
     1362        return fast_float(SR(f), *vars)
    13431363    except TypeError:
    13441364        pass
    13451365
     
    13491369    return f
    13501370
    13511371def is_fast_float(x):
    1352     return PY_TYPE_CHECK(x, FastDoubleFunc)
     1372    return PY_TYPE_CHECK(x, FastDoubleFunc) or PY_TYPE_CHECK(x, Wrapper)
     1373
  • new file sage/ext/gen_interpreters.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/ext/gen_interpreters.py
    - +  
     1r"""
     2Generate interpreters for fast_callable
     3
     4AUTHORS:
     5
     6- Carl Witty
     7
     8This file is part of the Sage support for "planned" computations;
     9that is, computations that are separated into a planning stage and
     10a plan-execution stage.  Here, we generate fast interpreters for plan
     11executions.
     12
     13There are at least two kinds of computations that are often planned in
     14this fashion.  First is arithmetic expression evaluation, where we
     15take an arbitrary user-specified arithmetic expression and compile it
     16into a bytecode form for fast interpretation.  Second is things like
     17FFTs and large multiplications, where large problems are split into
     18multiple smaller problems... we can do the logical "splitting" for a
     19given size only once, producing a plan which can be reused as often as
     20we want for different problems of the same size.  Currently only
     21arithmetic expression evaluation is implemented, but other kinds of
     22planned computations should be easy to add.
     23
     24Typically, for arithmetic expressions, we want the storage of
     25intermediate results to be handled automatically (on a stack); for
     26FFTs/multiplications/etc., the planner will keep track of intermediate
     27results itself.
     28
     29For arithmetic expression evaluation, we want to have lots of
     30interpreters (at least one, and possibly several, per
     31specially-handled type).  Also, for any given type, we have many
     32possible variants of instruction encoding, etc.; some of these could
     33be handled with conditional compilation, but some are more
     34complicated.  So we end up writing an interpreter generator.
     35
     36We want to share as much code as possible across all of these
     37interpreters, while still maintaining the freedom to make drastic
     38changes in the interpretation strategy (which may change the
     39generated code, the calling convention for the interpreter, etc.)
     40
     41To make this work, the interpreter back-end is divided into three
     42parts:
     43
     44  1. The interpreter itself, in C or C++.
     45
     46  2. The wrapper, which is a Cython object holding the
     47     constants, code, etc., and which actually calls the interpreter.
     48
     49  3. The code generator.
     50
     51We generate parts 1 and 2.  The code generator is table-driven,
     52and we generate the tables for the code generator.
     53
     54There are a lot of techniques for fast interpreters that we do not yet
     55use; hopefully at least some of these will eventually be implemented:
     56
     57- using gcc's "labels as values" extension where available
     58
     59- top-of-stack caching
     60
     61- superinstructions and/or superoperators
     62
     63- static stack caching
     64
     65- context threading/subrouting threading
     66
     67- selective inlining/dynamic superinstructions
     68
     69- automatic replication
     70
     71Interpreters may be stack-based or register-based.  Recent research
     72suggests that register-based interpreters are better, but the
     73researchers are investigating interpreters for entire programming
     74languages, rather than interpreters for expressions.  I suspect
     75that stack-based expression interpreters may be better.  However,
     76we'll implement both varieties and see what's best.
     77
     78The relative costs of stack- and register-based interpreters will
     79depend on the costs of moving values.  For complicated types (like
     80mpz_t), a register-based interpreter will quite likely be better,
     81since it will avoid moving values.
     82
     83We will NOT support any sort of storage of bytecode; instead, the
     84code must be re-generated from expression trees in every Sage run.
     85This means that we can trivially experiment with different styles of
     86interpreter, or even use quite different interpreters depending on
     87the architecture, without having to worry about forward and backward
     88compatibility.
     89"""
     90#*****************************************************************************
     91#       Copyright (C) 2009 Carl Witty <Carl.Witty@gmail.com>
     92#
     93#  Distributed under the terms of the GNU General Public License (GPL)
     94#                  http://www.gnu.org/licenses/
     95#*****************************************************************************
     96
     97from __future__ import with_statement
     98
     99import os
     100import re
     101from jinja import Environment
     102from jinja.datastructure import ComplainingUndefined
     103from collections import defaultdict
     104from distutils.extension import Extension
     105
     106##############################
     107# This module is used during the Sage buld process, so it should not
     108# use any other Sage modules.  (In particular, it MUST NOT use any
     109# Cython modules -- they won't be built yet!)
     110# Also, we have some trivial dependency tracking, where we don't
     111# rebuild the interpreters if this file hasn't changed; if
     112# interpreter configuation is split out into a separate file,
     113# that will have to be changed.
     114##############################
     115
     116
     117# We share a single jinja environment among all templating in this file.
     118# We use trim_blocks=True (which means that we ignore white space after
     119# "%}" jinja command endings), and set undefined_singleton to complain
     120# if we use an undefined variable.
     121jinja_env = Environment(trim_blocks=True,
     122                        undefined_singleton=ComplainingUndefined)
     123# Allow 'i' as a shorter alias for the built-in 'indent' filter.
     124jinja_env.filters['i'] = jinja_env.filters['indent']
     125
     126def indent_lines(n, text):
     127    r"""
     128    INPUTS:
     129        n -- indentation amount
     130        text -- text to indent
     131
     132    Indents each line in text by n spaces.
     133
     134    EXAMPLES:
     135        sage: from sage.ext.gen_interpreters import indent_lines
     136        sage: indent_lines(3, "foo")
     137        '   foo'
     138        sage: indent_lines(3, "foo\nbar")
     139        '   foo\n   bar'
     140        sage: indent_lines(3, "foo\nbar\n")
     141        '   foo\n   bar\n'
     142    """
     143    lines = text.splitlines(True)
     144    spaces = ' ' * n
     145    return ''.join(spaces + line for line in lines)
     146
     147def je(template, **kwargs):
     148    r"""
     149    A convenience method for creating strings with Jinja templates.
     150    The name je stands for "Jinja evaluate".
     151
     152    The first argument is the template string; remaining keyword
     153    arguments define Jinja variables.
     154
     155    If the first character in the template string is a newline, it is
     156    removed (this feature is useful when using multi-line templates defined
     157    with triple-quoted strings -- the first line doesn't have to be on
     158    the same line as the quotes, which would screw up the indentation).
     159
     160    (This is very inefficient, because it recompiles the Jinja
     161    template on each call; don't use it in situations where
     162    performance is important.)
     163
     164    EXAMPLES:
     165        sage: from sage.ext.gen_interpreters import je
     166        sage: je("{{ a }} > {{ b }} * {{ c }}", a='"a suffusion of yellow"', b=3, c=7)
     167        u'"a suffusion of yellow" > 3 * 7'
     168    """
     169    if len(template) > 0 and template[0] == '\n':
     170        template = template[1:]
     171
     172    # It looks like Jinja automatically removes one trailing newline?
     173    if len(template) > 0 and template[-1] == '\n':
     174        template = template + '\n'
     175
     176    tmpl = jinja_env.from_string(template)
     177    return tmpl.render(kwargs)
     178
     179class StorageType(object):
     180    r"""
     181    A StorageType specifies the C types used to deal with values of a
     182    given type.
     183
     184    We currently support three categories of types.
     185
     186    First are the "simple" types.  These are types where: the
     187    representation is small, functions expect arguments to be passed
     188    by value, and the C/C++ assignment operator works.  This would
     189    include built-in C types (long, float, etc.) and small structs
     190    (like gsl_complex).
     191
     192    Second is 'PyObject*'.  This is just like a simple type, except
     193    that we have to incref/decref at appropriate places.
     194
     195    Third is "auto-reference" types.  This is how
     196    GMP/MPIR/MPFR/MPFI/FLINT types work.  For these types, functions
     197    expect arguments to be passed by reference, and the C assignment
     198    operator does not do what we want.  In addition, they take
     199    advantage of a quirk in C (where arrays are automatically
     200    converted to pointers) to automatically pass arguments by
     201    reference.
     202
     203    Support for further categories would not be difficult to add (such
     204    as reference-counted types other than PyObject*, or
     205    pass-by-reference types that don't use the GMP auto-reference
     206    trick), if we ever run across a use for them.
     207    """
     208
     209    def __init__(self):
     210        r"""
     211        Initialize an instance of StorageType.
     212
     213        This sets several properties:
     214
     215        class_member_declarations:
     216        A string giving variable declarations that must be members of any
     217        wrapper class using this type.
     218
     219        class_member_initializations:
     220        A string initializing the class_member_declarations; will be
     221        inserted into the __init__ method of any wrapper class using this
     222        type.
     223
     224        local_declarations:
     225        A string giving variable declarations that must be local variables
     226        in Cython methods using this storage type.
     227
     228        EXAMPLES:
     229            sage: from sage.ext.gen_interpreters import *
     230            sage: ty_double.class_member_declarations
     231            ''
     232            sage: ty_double.class_member_initializations
     233            ''
     234            sage: ty_double.local_declarations
     235            ''
     236            sage: ty_mpfr.class_member_declarations
     237            'cdef RealField domain\n'
     238            sage: ty_mpfr.class_member_initializations
     239            "self.domain = args['domain']\n"
     240            sage: ty_mpfr.local_declarations
     241            'cdef RealNumber rn\n'
     242        """
     243        self.class_member_declarations = ''
     244        self.class_member_initializations = ''
     245        self.local_declarations = ''
     246
     247    def cheap_copies(self):
     248        r"""
     249        Returns True or False, depending on whether this StorageType
     250        supports cheap copies -- whether it is cheap to copy values of
     251        this type from one location to another.  This is true for
     252        primitive types, and for types like PyObject* (where you're only
     253        copying a pointer, and possibly changing some reference counts).
     254        It is false for types like mpz_t and mpfr_t, where copying values
     255        can involve arbitrarily much work (including memory allocation).
     256
     257        The practical effect is that if cheap_copies is True,
     258        instructions with outputs of this type write the results into
     259        local variables, and the results are then copied to their
     260        final locations.  If cheap_copies is False, then the addresses
     261        of output locations are passed into the instruction and the
     262        instruction writes outputs directly in the final location.
     263
     264        EXAMPLES:
     265            sage: from sage.ext.gen_interpreters import *
     266            sage: ty_double.cheap_copies()
     267            True
     268            sage: ty_python.cheap_copies()
     269            True
     270            sage: ty_mpfr.cheap_copies()
     271            False
     272        """
     273        return False
     274
     275    def python_refcounted(self):
     276        r"""
     277        Says whether this storage type is a Python type, so we need to
     278        use INCREF/DECREF.
     279
     280        (If we needed to support any non-Python refcounted types, it
     281        might be better to make this object-oriented and have methods
     282        like "generate an incref" and "generate a decref".  But as
     283        long as we only support Python, this way is probably simpler.)
     284
     285        EXAMPLES:
     286            sage: from sage.ext.gen_interpreters import *
     287            sage: ty_double.python_refcounted()
     288            False
     289            sage: ty_python.python_refcounted()
     290            True
     291        """
     292        return False
     293
     294    def cython_decl_type(self):
     295        r"""
     296        Gives the Cython type for a single value of this type (as a string).
     297
     298        EXAMPLES:
     299            sage: from sage.ext.gen_interpreters import *
     300            sage: ty_double.cython_decl_type()
     301            'double'
     302            sage: ty_python.cython_decl_type()
     303            'object'
     304            sage: ty_mpfr.cython_decl_type()
     305            'mpfr_t'
     306        """
     307        return self.c_decl_type()
     308
     309    def cython_array_type(self):
     310        r"""
     311        Gives the Cython type for referring to an array of values of
     312        this type (as a string).
     313
     314        EXAMPLES:
     315            sage: from sage.ext.gen_interpreters import *
     316            sage: ty_double.cython_array_type()
     317            'double*'
     318            sage: ty_python.cython_array_type()
     319            'PyObject**'
     320            sage: ty_mpfr.cython_array_type()
     321            'mpfr_t*'
     322        """
     323        return self.c_ptr_type()
     324
     325    def needs_cython_init_clear(self):
     326        r"""
     327        Says whether values/arrays of this type need to be initialized
     328        before use and cleared before the underlying memory is freed.
     329
     330        (We could remove this method, always call .cython_init() to
     331        generate initialization code, and just let .cython_init()
     332        generate empty code if no initialization is required; that would
     333        generate empty loops, which are ugly and potentially might not
     334        be optimized away.)
     335
     336        EXAMPLES:
     337            sage: from sage.ext.gen_interpreters import *
     338            sage: ty_double.needs_cython_init_clear()
     339            False
     340            sage: ty_mpfr.needs_cython_init_clear()
     341            True
     342            sage: ty_python.needs_cython_init_clear()
     343            True
     344        """
     345        return False
     346
     347    def c_decl_type(self):
     348        r"""
     349        Gives the C type for a single value of this type (as a string).
     350
     351        EXAMPLES:
     352            sage: from sage.ext.gen_interpreters import *
     353            sage: ty_double.c_decl_type()
     354            'double'
     355            sage: ty_python.c_decl_type()
     356            'PyObject*'
     357            sage: ty_mpfr.c_decl_type()
     358            'mpfr_t'
     359        """
     360        raise NotImplementedError
     361
     362    def c_ptr_type(self):
     363        r"""
     364        Gives the C type for a pointer to this type (as a reference to
     365        either a single value or an array) (as a string).
     366
     367        EXAMPLES:
     368            sage: from sage.ext.gen_interpreters import *
     369            sage: ty_double.c_ptr_type()
     370            'double*'
     371            sage: ty_python.c_ptr_type()
     372            'PyObject**'
     373            sage: ty_mpfr.c_ptr_type()
     374            'mpfr_t*'
     375        """
     376        return self.c_decl_type() + '*'
     377
     378    def c_local_type(self):
     379        r"""
     380        Gives the C type used for a value of this type inside an
     381        instruction.  For assignable/cheap_copy types, this is the
     382        same as c_decl_type; for auto-reference types, this is the
     383        pointer type.
     384
     385        EXAMPLES:
     386            sage: from sage.ext.gen_interpreters import *
     387            sage: ty_double.c_local_type()
     388            'double'
     389            sage: ty_python.c_local_type()
     390            'PyObject*'
     391            sage: ty_mpfr.c_local_type()
     392            'mpfr_ptr'
     393        """
     394        raise NotImplementedError
     395
     396    def assign_c_from_py(self, c, py):
     397        r"""
     398        Given a Cython variable/array reference/etc. of this storage type,
     399        and a Python expression, generate code to assign to the Cython
     400        variable from the Python expression.
     401
     402        EXAMPLES:
     403            sage: from sage.ext.gen_interpreters import *
     404            sage: ty_double.assign_c_from_py('foo', 'bar')
     405            u'foo = bar'
     406            sage: ty_python.assign_c_from_py('foo[i]', 'bar[j]')
     407            u'foo[i] = <PyObject *>bar[j]; Py_INCREF(foo[i])'
     408            sage: ty_mpfr.assign_c_from_py('foo', 'bar')
     409            u'rn = self.domain(bar)\nmpfr_set(foo, rn.value, GMP_RNDN)'
     410        """
     411        return je("{{ c }} = {{ py }}", c=c, py=py)
     412
     413    def declare_chunk_class_members(self, name):
     414        r"""
     415        Returns a string giving the declarations of the class members
     416        in a wrapper class for a memory chunk with this storage type
     417        and the given name.
     418
     419        EXAMPLES:
     420            sage: from sage.ext.gen_interpreters import *
     421            sage: ty_mpfr.declare_chunk_class_members('args')
     422            u'    cdef int _n_args\n    cdef mpfr_t* _args\n'
     423        """
     424        return je("""
     425{# XXX Variables here (and everywhere, really) should actually be Py_ssize_t #}
     426    cdef int _n_{{ name }}
     427    cdef {{ self.cython_array_type() }} _{{ name }}
     428""", self=self, name=name)
     429
     430    def alloc_chunk_data(self, name, len):
     431        r"""
     432        Returns a string allocating the memory for the class members for
     433        a memory chunk with this storage type and the given name.
     434
     435        EXAMPLES:
     436            sage: from sage.ext.gen_interpreters import *
     437            sage: print ty_mpfr.alloc_chunk_data('args', 'MY_LENGTH')
     438                    self._n_args = MY_LENGTH
     439                    self._args = <mpfr_t*>sage_malloc(sizeof(mpfr_t) * MY_LENGTH)
     440                    if self._args == NULL: raise MemoryError
     441                    for i in range(MY_LENGTH):
     442                        mpfr_init2(self._args[i], self.domain.prec())
     443            <BLANKLINE>
     444        """
     445        return je("""
     446        self._n_{{ name }} = {{ len }}
     447        self._{{ name }} = <{{ self.c_ptr_type() }}>sage_malloc(sizeof({{ self.c_decl_type() }}) * {{ len }})
     448        if self._{{ name }} == NULL: raise MemoryError
     449{% if self.needs_cython_init_clear() %}
     450        for i in range({{ len }}):
     451            {{ self.cython_init('self._%s[i]' % name) }}
     452{% endif %}
     453""", self=self, name=name, len=len)
     454
     455    def dealloc_chunk_data(self, name):
     456        r"""
     457        Returns a string to be put in the __dealloc__ method of a
     458        wrapper class using a memory chunk with this storage type, to
     459        deallocate the corresponding class members.
     460
     461        EXAMPLES:
     462            sage: from sage.ext.gen_interpreters import *
     463            sage: print ty_double.dealloc_chunk_data('args')
     464                    if self._args:
     465                        sage_free(self._args)
     466            <BLANKLINE>
     467            sage: print ty_mpfr.dealloc_chunk_data('constants')
     468                    if self._constants:
     469                        for i in range(self._n_constants):
     470                            mpfr_clear(self._constants[i])
     471                        sage_free(self._constants)
     472            <BLANKLINE>
     473        """
     474        return je("""
     475        if self._{{ name }}:
     476{%     if self.needs_cython_init_clear() %}
     477            for i in range(self._n_{{ name }}):
     478                {{ self.cython_clear('self._%s[i]' % name) }}
     479{%     endif %}
     480            sage_free(self._{{ name }})
     481""", self=self, name=name)
     482
     483class StorageTypeAssignable(StorageType):
     484    r"""
     485    StorageTypeAssignable is a subtype of StorageType that deals with
     486    types with cheap copies, like primitive types and PyObject*.
     487    """
     488
     489    def __init__(self, ty):
     490        r"""
     491        Initializes the property type (the C/Cython name for this type),
     492        as well as the properties described in the documentation for
     493        StorageType.__init__.
     494
     495        EXAMPLES:
     496            sage: from sage.ext.gen_interpreters import *
     497            sage: ty_double.class_member_declarations
     498            ''
     499            sage: ty_double.class_member_initializations
     500            ''
     501            sage: ty_double.local_declarations
     502            ''
     503            sage: ty_double.type
     504            'double'
     505            sage: ty_python.type
     506            'PyObject*'
     507        """
     508        StorageType.__init__(self)
     509        self.type = ty
     510
     511    def cheap_copies(self):
     512        r"""
     513        Returns True or False, depending on whether this StorageType
     514        supports cheap copies -- whether it is cheap to copy values of
     515        this type from one location to another.  (See StorageType.cheap_copies
     516        for more on this property.)
     517
     518        Since having cheap copies is essentially the definition of
     519        StorageTypeAssignable, this always returns True.
     520
     521        EXAMPLES:
     522            sage: from sage.ext.gen_interpreters import *
     523            sage: ty_double.cheap_copies()
     524            True
     525            sage: ty_python.cheap_copies()
     526            True
     527        """
     528        return True
     529
     530    def c_decl_type(self):
     531        r"""
     532        Gives the C type for a single value of this type (as a string).
     533
     534        EXAMPLES:
     535            sage: from sage.ext.gen_interpreters import *
     536            sage: ty_double.c_decl_type()
     537            'double'
     538            sage: ty_python.c_decl_type()
     539            'PyObject*'
     540        """
     541        return self.type
     542
     543    def c_local_type(self):
     544        r"""
     545        Gives the C type used for a value of this type inside an
     546        instruction.  For assignable/cheap_copy types, this is the
     547        same as c_decl_type; for auto-reference types, this is the
     548        pointer type.
     549
     550        EXAMPLES:
     551            sage: from sage.ext.gen_interpreters import *
     552            sage: ty_double.c_local_type()
     553            'double'
     554            sage: ty_python.c_local_type()
     555            'PyObject*'
     556        """
     557        return self.type
     558
     559class StorageTypeSimple(StorageTypeAssignable):
     560    r"""
     561    StorageTypeSimple is a subtype of StorageTypeAssignable that deals
     562    with non-reference-counted types with cheap copies, like primitive
     563    types.  As of yet, it has no functionality differences from
     564    StorageTypeAssignable.
     565    """
     566    pass
     567
     568ty_int = StorageTypeSimple('int')
     569ty_double = StorageTypeSimple('double')
     570
     571class StorageTypePython(StorageTypeAssignable):
     572    r"""
     573    StorageTypePython is a subtype of StorageTypeAssignable that deals
     574    with Python objects.
     575
     576    Just allocating an array full of PyObject* leads to problems,
     577    because the Python garbage collector must be able to get to every
     578    Python object, and it wouldn't know how to get to these arrays.
     579    So we allocate the array as a Python list, but then we immediately
     580    pull the ob_item out of it and deal only with that from then on.
     581
     582    We often leave these lists with NULL entries.  This is safe for
     583    the garbage collector and the deallocator, which is all we care
     584    about; but it would be unsafe to provide Python-level access to
     585    these lists.
     586
     587    There is one special thing about StorageTypePython: memory that is
     588    used by the interpreter as scratch space (for example, the stack)
     589    must be cleared after each call (so we don't hold on to
     590    potentially-large objects and waste memory).  Since we have to do
     591    this anyway, the interpreter gains a tiny bit of speed by assuming
     592    that the scratch space is cleared on entry; for example, when
     593    pushing a value onto the stack, it doesn't bother to XDECREF the
     594    previous value because it's always NULL.
     595    """
     596
     597    def __init__(self):
     598        r"""
     599        Initializes the properties described in the documentation
     600        for StorageTypeAssignable.__init__.  The type is always
     601        'PyObject*'.
     602
     603        EXAMPLES:
     604            sage: from sage.ext.gen_interpreters import *
     605            sage: ty_python.class_member_declarations
     606            ''
     607            sage: ty_python.class_member_initializations
     608            ''
     609            sage: ty_python.local_declarations
     610            ''
     611            sage: ty_python.type
     612            'PyObject*'
     613        """
     614        StorageTypeAssignable.__init__(self, 'PyObject*')
     615
     616    def python_refcounted(self):
     617        r"""
     618        Says whether this storage type is a Python type, so we need to
     619        use INCREF/DECREF.
     620
     621        Returns True.
     622
     623        EXAMPLES:
     624            sage: from sage.ext.gen_interpreters import *
     625            sage: ty_python.python_refcounted()
     626            True
     627        """
     628        return True
     629
     630    def cython_decl_type(self):
     631        r"""
     632        Gives the Cython type for a single value of this type (as a string).
     633
     634        EXAMPLES:
     635            sage: from sage.ext.gen_interpreters import *
     636            sage: ty_python.cython_decl_type()
     637            'object'
     638        """
     639        return 'object'
     640
     641    def declare_chunk_class_members(self, name):
     642        r"""
     643        Returns a string giving the declarations of the class members
     644        in a wrapper class for a memory chunk with this storage type
     645        and the given name.
     646
     647        EXAMPLES:
     648            sage: from sage.ext.gen_interpreters import *
     649            sage: ty_python.declare_chunk_class_members('args')
     650            u'    cdef object _list_args\n    cdef int _n_args\n    cdef PyObject** _args\n'
     651        """
     652        return je("""
     653    cdef object _list_{{ name }}
     654    cdef int _n_{{ name }}
     655    cdef {{ self.cython_array_type() }} _{{ name }}
     656""", self=self, name=name)
     657
     658    def alloc_chunk_data(self, name, len):
     659        r"""
     660        Returns a string allocating the memory for the class members for
     661        a memory chunk with this storage type and the given name.
     662
     663        EXAMPLES:
     664            sage: from sage.ext.gen_interpreters import *
     665            sage: print ty_python.alloc_chunk_data('args', 'MY_LENGTH')
     666                    self._n_args = MY_LENGTH
     667                    self._list_args = PyList_New(self._n_args)
     668                    self._args = (<PyListObject *>self._list_args).ob_item
     669            <BLANKLINE>
     670        """
     671        return je("""
     672        self._n_{{ name }} = {{ len }}
     673        self._list_{{ name }} = PyList_New(self._n_{{ name }})
     674        self._{{ name }} = (<PyListObject *>self._list_{{ name }}).ob_item
     675""", self=self, name=name, len=len)
     676
     677    def dealloc_chunk_data(self, name):
     678        r"""
     679        Returns a string to be put in the __dealloc__ method of a
     680        wrapper class using a memory chunk with this storage type, to
     681        deallocate the corresponding class members.
     682
     683        Our array was allocated as a Python list; this means we actually
     684        don't need to do anything to deallocate it.
     685
     686        EXAMPLES:
     687            sage: from sage.ext.gen_interpreters import *
     688            sage: ty_python.dealloc_chunk_data('args')
     689            ''
     690        """
     691        return ''
     692
     693    def needs_cython_init_clear(self):
     694        r"""
     695        Says whether values/arrays of this type need to be initialized
     696        before use and cleared before the underlying memory is freed.
     697
     698        Returns True.
     699
     700        EXAMPLES:
     701            sage: from sage.ext.gen_interpreters import *
     702            sage: ty_python.needs_cython_init_clear()
     703            True
     704        """
     705        return True
     706
     707    def assign_c_from_py(self, c, py):
     708        r"""
     709        Given a Cython variable/array reference/etc. of this storage type,
     710        and a Python expression, generate code to assign to the Cython
     711        variable from the Python expression.
     712
     713        EXAMPLES:
     714            sage: from sage.ext.gen_interpreters import *
     715            sage: ty_python.assign_c_from_py('foo[i]', 'bar[j]')
     716            u'foo[i] = <PyObject *>bar[j]; Py_INCREF(foo[i])'
     717        """
     718        return je("""{{ c }} = <PyObject *>{{ py }}; Py_INCREF({{ c }})""",
     719                  c=c, py=py)
     720
     721    def cython_init(self, loc):
     722        r"""
     723        Generates code to initialize a variable (or array reference)
     724        holding a PyObject*.  Sets it to NULL.
     725
     726        EXAMPLES:
     727            sage: from sage.ext.gen_interpreters import *
     728            sage: ty_python.cython_init('foo[i]')
     729            u'foo[i] = NULL'
     730        """
     731        return je("{{ loc }} = NULL", loc=loc)
     732
     733    def cython_clear(self, loc):
     734        r"""
     735        Generates code to clear a variable (or array reference) holding
     736        a PyObject*.
     737
     738        EXAMPLES:
     739            sage: from sage.ext.gen_interpreters import *
     740            sage: ty_python.cython_clear('foo[i]')
     741            u'Py_CLEAR(foo[i])'
     742        """
     743        return je("Py_CLEAR({{ loc }})", loc=loc)
     744
     745ty_python = StorageTypePython()
     746
     747class StorageTypeAutoReference(StorageType):
     748    r"""
     749    StorageTypeAutoReference is a subtype of StorageType that deals with
     750    types in the style of GMP/MPIR/MPFR/MPFI/FLINT, where copies are
     751    not cheap, functions expect arguments to be passed by reference,
     752    and the API takes advantage of the C quirk where arrays are
     753    automatically converted to pointers to automatically pass
     754    arguments by reference.
     755    """
     756
     757    def __init__(self, decl_ty, ref_ty):
     758        r"""
     759        Initializes the properties decl_type and ref_type (the C type
     760        names used when declaring variables and function parameters,
     761        respectively), as well as the properties described in
     762        the documentation for StorageType.__init__.
     763
     764        EXAMPLES:
     765            sage: from sage.ext.gen_interpreters import *
     766            sage: ty_mpfr.class_member_declarations
     767            'cdef RealField domain\n'
     768            sage: ty_mpfr.class_member_initializations
     769            "self.domain = args['domain']\n"
     770            sage: ty_mpfr.local_declarations
     771            'cdef RealNumber rn\n'
     772            sage: ty_mpfr.decl_type
     773            'mpfr_t'
     774            sage: ty_mpfr.ref_type
     775            'mpfr_ptr'
     776        """
     777        StorageType.__init__(self)
     778        self.decl_type = decl_ty
     779        self.ref_type = ref_ty
     780
     781    def c_decl_type(self):
     782        r"""
     783        Gives the C type for a single value of this type (as a string).
     784
     785        EXAMPLES:
     786            sage: from sage.ext.gen_interpreters import *
     787            sage: ty_mpfr.c_decl_type()
     788            'mpfr_t'
     789        """
     790        return self.decl_type
     791
     792    def c_local_type(self):
     793        r"""
     794        Gives the C type used for a value of this type inside an
     795        instruction.  For assignable/cheap_copy types, this is the
     796        same as c_decl_type; for auto-reference types, this is the
     797        pointer type.
     798
     799        EXAMPLES:
     800            sage: from sage.ext.gen_interpreters import *
     801            sage: ty_mpfr.c_local_type()
     802            'mpfr_ptr'
     803        """
     804        return self.ref_type
     805
     806    def needs_cython_init_clear(self):
     807        r"""
     808        Says whether values/arrays of this type need to be initialized
     809        before use and cleared before the underlying memory is freed.
     810
     811        All known examples of auto-reference types do need a special
     812        initialization call, so this always returns True.
     813
     814        EXAMPLES:
     815            sage: from sage.ext.gen_interpreters import *
     816            sage: ty_mpfr.needs_cython_init_clear()
     817            True
     818        """
     819        return True
     820
     821class StorageTypeMPFR(StorageTypeAutoReference):
     822    r"""
     823    StorageTypeMPFR is a subtype of StorageTypeAutoReference that deals
     824    the MPFR's mpfr_t type.
     825
     826    For any given program that we're interpreting, ty_mpfr can only
     827    refer to a single precision.  An interpreter that needs to use
     828    two precisions of mpfr_t in the same program should instantiate two
     829    separate instances of StorageTypeMPFR.  (Interpreters that need
     830    to handle arbitrarily many precisions in the same program are not
     831    handled at all.)
     832    """
     833
     834    def __init__(self, id=''):
     835        r"""
     836        Initializes the id property, as well as the properties described
     837        in the documentation for StorageTypeAutoReference.__init__.
     838
     839        The id property is used if you want to have an interpreter
     840        that handles two instances of StorageTypeMPFR (that is,
     841        handles mpfr_t variables at two different precisions
     842        simultaneously).  It's a string that's used to generate
     843        variable names that don't conflict.  (The id system has
     844        never actually been used, so bugs probably remain.)
     845
     846        EXAMPLES:
     847            sage: from sage.ext.gen_interpreters import *
     848            sage: ty_mpfr.class_member_declarations
     849            'cdef RealField domain\n'
     850            sage: ty_mpfr.class_member_initializations
     851            "self.domain = args['domain']\n"
     852            sage: ty_mpfr.local_declarations
     853            'cdef RealNumber rn\n'
     854            sage: ty_mpfr.decl_type
     855            'mpfr_t'
     856            sage: ty_mpfr.ref_type
     857            'mpfr_ptr'
     858
     859        TESTS:
     860            sage: ty_mpfr2 = StorageTypeMPFR(id='_the_second')
     861            sage: ty_mpfr2.class_member_declarations
     862            'cdef RealField domain_the_second\n'
     863            sage: ty_mpfr2.class_member_initializations
     864            "self.domain_the_second = args['domain_the_second']\n"
     865            sage: ty_mpfr2.local_declarations
     866            'cdef RealNumber rn_the_second\n'
     867        """
     868        StorageTypeAutoReference.__init__(self, 'mpfr_t', 'mpfr_ptr')
     869        self.id = id
     870        self.class_member_declarations = "cdef RealField domain%s\n" % self.id
     871        self.class_member_initializations = \
     872            "self.domain%s = args['domain%s']\n" % (self.id, self.id)
     873        self.local_declarations = "cdef RealNumber rn%s\n" % self.id
     874
     875    def cython_init(self, loc):
     876        r"""
     877        Generates code to initialize an mpfr_t reference (a variable, an
     878        array reference, etc.)
     879
     880        EXAMPLES:
     881            sage: from sage.ext.gen_interpreters import *
     882            sage: ty_mpfr.cython_init('foo[i]')
     883            u'mpfr_init2(foo[i], self.domain.prec())'
     884        """
     885        return je("mpfr_init2({{ loc }}, self.domain{{ self.id }}.prec())",
     886                  self=self, loc=loc)
     887
     888    def cython_clear(self, loc):
     889        r"""
     890        Generates code to clear an mpfr_t reference (a variable, an
     891        array reference, etc.)
     892
     893        EXAMPLES:
     894            sage: from sage.ext.gen_interpreters import *
     895            sage: ty_mpfr.cython_clear('foo[i]')
     896            'mpfr_clear(foo[i])'
     897        """
     898        return 'mpfr_clear(%s)' % loc
     899
     900    def assign_c_from_py(self, c, py):
     901        r"""
     902        Given a Cython variable/array reference/etc. of this storage type,
     903        and a Python expression, generate code to assign to the Cython
     904        variable from the Python expression.
     905
     906        EXAMPLES:
     907            sage: from sage.ext.gen_interpreters import *
     908            sage: ty_mpfr.assign_c_from_py('foo[i]', 'bar[j]')
     909            u'rn = self.domain(bar[j])\nmpfr_set(foo[i], rn.value, GMP_RNDN)'
     910        """
     911        return je("""
     912rn{{ self.id }} = self.domain({{ py }})
     913mpfr_set({{ c }}, rn.value, GMP_RNDN)""", self=self, c=c, py=py)
     914
     915ty_mpfr = StorageTypeMPFR()
     916
     917class MemoryChunk(object):
     918    r"""
     919    Memory chunks control allocation, deallocation, iniialization,
     920    etc.  of the vectors and objects in the interpreter.  Basically,
     921    there is one memory chunk per argument to the C interpreter.
     922
     923    There are three "generic" varieties of memory chunk: "constants",
     924    "arguments", and "scratch".  These are named after their most
     925    common use, but they could be used for other things in some
     926    interpreters.
     927
     928    All three kinds of chunks are allocated in the wrapper class.
     929    Constants are initialized when the wrapper is constructed;
     930    arguments are initialized in the __call__ method, from the
     931    caller's arguments.  "scratch" chunks are not initialized at all;
     932    they are used for scratch storage (often, but not necessarily, for
     933    a stack) in the interpreter.
     934
     935    Interpreters which need memory chunks that don't fit into these
     936    categories can create new subclasses of MemoryChunk.
     937    """
     938
     939    def __init__(self, name, storage_type):
     940        r"""
     941        Initialize an instance of MemoryChunk.
     942
     943        This sets the properties "name" (the name of this memory chunk;
     944        used in generated variable names, etc.) and "storage_type",
     945        which is a StorageType object.
     946
     947        EXAMPLES:
     948            sage: from sage.ext.gen_interpreters import *
     949            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     950            sage: mc.name
     951            'args'
     952            sage: mc.storage_type is ty_mpfr
     953            True
     954        """
     955        self.name = name
     956        self.storage_type = storage_type
     957
     958    def __repr__(self):
     959        r"""
     960        Give a string representation of this memory chunk.
     961
     962        EXAMPLES:
     963            sage: from sage.ext.gen_interpreters import *
     964            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     965            sage: mc
     966            {MC:args}
     967            sage: mc.__repr__()
     968            '{MC:args}'
     969        """
     970        return '{MC:%s}' % self.name
     971
     972    def declare_class_members(self):
     973        r"""
     974        Returns a string giving the declarations of the class members
     975        in a wrapper class for this memory chunk.
     976
     977        EXAMPLES:
     978            sage: from sage.ext.gen_interpreters import *
     979            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     980            sage: mc.declare_class_members()
     981            u'    cdef int _n_args\n    cdef mpfr_t* _args\n'
     982        """
     983        return self.storage_type.declare_chunk_class_members(self.name)
     984
     985    def init_class_members(self):
     986        r"""
     987        Returns a string to be put in the __init__ method of a wrapper
     988        class using this memory chunk, to initialize the corresponding
     989        class members.
     990
     991        EXAMPLES:
     992            sage: from sage.ext.gen_interpreters import *
     993            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     994            sage: print mc.init_class_members()
     995                    count = args['args']
     996                    self._n_args = count
     997                    self._args = <mpfr_t*>sage_malloc(sizeof(mpfr_t) * count)
     998                    if self._args == NULL: raise MemoryError
     999                    for i in range(count):
     1000                        mpfr_init2(self._args[i], self.domain.prec())
     1001            <BLANKLINE>
     1002        """
     1003        return ""
     1004
     1005    def dealloc_class_members(self):
     1006        r"""
     1007        Returns a string to be put in the __dealloc__ method of a wrapper
     1008        class using this memory chunk, to deallocate the corresponding
     1009        class members.
     1010
     1011        EXAMPLES:
     1012            sage: from sage.ext.gen_interpreters import *
     1013            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     1014            sage: print mc.dealloc_class_members()
     1015                    if self._args:
     1016                        for i in range(self._n_args):
     1017                            mpfr_clear(self._args[i])
     1018                        sage_free(self._args)
     1019            <BLANKLINE>
     1020        """
     1021        return ""
     1022
     1023    def declare_parameter(self):
     1024        r"""
     1025        Returns the string to use to declare the interpreter parameter
     1026        corresponding to this memory chunk.
     1027
     1028        EXAMPLES:
     1029            sage: from sage.ext.gen_interpreters import *
     1030            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     1031            sage: mc.declare_parameter()
     1032            'mpfr_t* args'
     1033        """
     1034        return '%s %s' % (self.storage_type.c_ptr_type(), self.name)
     1035
     1036    def declare_call_locals(self):
     1037        r"""
     1038        Returns a string to put in the __call__ method of a wrapper
     1039        class using this memory chunk, to allocate local variables.
     1040       
     1041        EXAMPLES:
     1042            sage: from sage.ext.gen_interpreters import *
     1043            sage: mc = MemoryChunkRRRetval('retval', ty_mpfr)
     1044            sage: mc.declare_call_locals()
     1045            u'        cdef RealNumber retval = (self.domain)()\n'
     1046        """
     1047        return ""
     1048
     1049    def pass_argument(self):
     1050        r"""
     1051        Returns the string to pass the argument corresponding to this
     1052        memory chunk to the interpreter.
     1053
     1054        EXAMPLES:
     1055            sage: from sage.ext.gen_interpreters import *
     1056            sage: mc = MemoryChunkConstants('constants', ty_mpfr)
     1057            sage: mc.pass_argument()
     1058            'self._constants'
     1059        """
     1060        raise NotImplementedError
     1061
     1062    def needs_cleanup_on_error(self):
     1063        r"""
     1064        In an interpreter that can terminate prematurely (due to an
     1065        exception from calling Python code, or divide by zero, or
     1066        whatever) it will just return at the end of the current instruction,
     1067        skipping the rest of the program.  Thus, it may still have
     1068        values pushed on the stack, etc.
     1069
     1070        This method returns True if this memory chunk is modified by the
     1071        interpreter and needs some sort of cleanup when an error happens.
     1072
     1073        EXAMPLES:
     1074            sage: from sage.ext.gen_interpreters import *
     1075            sage: mc = MemoryChunkConstants('constants', ty_mpfr)
     1076            sage: mc.needs_cleanup_on_error()
     1077            False
     1078        """
     1079        return False
     1080
     1081    def is_stack(self):
     1082        r"""
     1083        Says whether this memory chunk is a stack.  This affects code
     1084        generation for instructions using this memory chunk.
     1085
     1086        It would be nicer to make this object-oriented somehow, so
     1087        that the code generator called MemoryChunk methods instead of
     1088        using
     1089            if ch.is_stack():
     1090                ... hardcoded stack code
     1091            else:
     1092                ... hardcoded non-stack code
     1093        but that hasn't been done yet.
     1094
     1095        EXAMPLES:
     1096            sage: from sage.ext.gen_interpreters import *
     1097            sage: mc = MemoryChunkScratch('scratch', ty_mpfr)
     1098            sage: mc.is_stack()
     1099            False
     1100            sage: mc = MemoryChunkScratch('stack', ty_mpfr, is_stack=True)
     1101            sage: mc.is_stack()
     1102            True
     1103        """
     1104        return False
     1105
     1106    def is_python_refcounted_stack(self):
     1107        r"""
     1108        Says whether this memory chunk refers to a stack where the entries
     1109        need to be INCREF/DECREF'ed.
     1110
     1111        It would be nice to make this object-oriented, so that the
     1112        code generator called MemoryChunk methods to do the potential
     1113        INCREF/DECREF and didn't have to explicitly test
     1114        is_python_refcounted_stack.
     1115
     1116        EXAMPLES:
     1117            sage: from sage.ext.gen_interpreters import *
     1118            sage: mc = MemoryChunkScratch('args', ty_python)
     1119            sage: mc.is_python_refcounted_stack()
     1120            False
     1121            sage: mc = MemoryChunkScratch('args', ty_python, is_stack=True)
     1122            sage: mc.is_python_refcounted_stack()
     1123            True
     1124            sage: mc = MemoryChunkScratch('args', ty_mpfr, is_stack=True)
     1125            sage: mc.is_python_refcounted_stack()
     1126            False
     1127        """
     1128        return self.is_stack() and self.storage_type.python_refcounted()
     1129
     1130class MemoryChunkLonglivedArray(MemoryChunk):
     1131    r"""
     1132    MemoryChunkLonglivedArray is a subtype of MemoryChunk that deals
     1133    with memory chunks that are both 1) allocated as class members (rather
     1134    than being allocated in __call__) and 2) are arrays.
     1135    """
     1136
     1137    def init_class_members(self):
     1138        r"""
     1139        Returns a string to be put in the __init__ method of a wrapper
     1140        class using this memory chunk, to initialize the corresponding
     1141        class members.
     1142
     1143        EXAMPLES:
     1144            sage: from sage.ext.gen_interpreters import *
     1145            sage: mc = MemoryChunkArguments('args', ty_double)
     1146            sage: print mc.init_class_members()
     1147                    count = args['args']
     1148                    self._n_args = count
     1149                    self._args = <double*>sage_malloc(sizeof(double) * count)
     1150                    if self._args == NULL: raise MemoryError
     1151            <BLANKLINE>
     1152        """
     1153        return je("""
     1154        count = args['{{ self.name }}']
     1155{% print self.storage_type.alloc_chunk_data(self.name, 'count') %}
     1156""", self=self)
     1157
     1158    def dealloc_class_members(self):
     1159        r"""
     1160        Returns a string to be put in the __dealloc__ method of a wrapper
     1161        class using this memory chunk, to deallocate the corresponding
     1162        class members.
     1163
     1164        EXAMPLES:
     1165            sage: from sage.ext.gen_interpreters import *
     1166            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     1167            sage: print mc.dealloc_class_members()
     1168                    if self._args:
     1169                        for i in range(self._n_args):
     1170                            mpfr_clear(self._args[i])
     1171                        sage_free(self._args)
     1172            <BLANKLINE>
     1173        """
     1174        return self.storage_type.dealloc_chunk_data(self.name)
     1175
     1176    def pass_argument(self):
     1177        r"""
     1178        Returns the string to pass the argument corresponding to this
     1179        memory chunk to the interpreter.
     1180
     1181        EXAMPLES:
     1182            sage: from sage.ext.gen_interpreters import *
     1183            sage: mc = MemoryChunkConstants('constants', ty_mpfr)
     1184            sage: mc.pass_argument()
     1185            'self._constants'
     1186        """
     1187        return 'self._%s' % self.name
     1188
     1189class MemoryChunkConstants(MemoryChunkLonglivedArray):
     1190    r"""
     1191    MemoryChunkConstants is a subtype of MemoryChunkLonglivedArray.
     1192
     1193    MemoryChunkConstants chunks have their contents set in the
     1194    wrapper's __init__ method (and not changed afterward).
     1195    """
     1196
     1197    def init_class_members(self):
     1198        r"""
     1199        Returns a string to be put in the __init__ method of a wrapper
     1200        class using this memory chunk, to initialize the corresponding
     1201        class members.
     1202
     1203        EXAMPLES:
     1204            sage: from sage.ext.gen_interpreters import *
     1205            sage: mc = MemoryChunkConstants('constants', ty_mpfr)
     1206            sage: print mc.init_class_members()
     1207                    val = args['constants']
     1208                    self._n_constants = len(val)
     1209                    self._constants = <mpfr_t*>sage_malloc(sizeof(mpfr_t) * len(val))
     1210                    if self._constants == NULL: raise MemoryError
     1211                    for i in range(len(val)):
     1212                        mpfr_init2(self._constants[i], self.domain.prec())
     1213                    for i in range(len(val)):
     1214                        rn = self.domain(val[i])
     1215                        mpfr_set(self._constants[i], rn.value, GMP_RNDN)
     1216            <BLANKLINE>
     1217        """
     1218        return je("""
     1219        val = args['{{ self.name }}']
     1220{% print self.storage_type.alloc_chunk_data(self.name, 'len(val)') %}
     1221        for i in range(len(val)):
     1222            {{ self.storage_type.assign_c_from_py('self._%s[i]' % self.name, 'val[i]') | i(12) }}
     1223""", self=self)
     1224
     1225class MemoryChunkArguments(MemoryChunkLonglivedArray):
     1226    r"""
     1227    MemoryChunkArguments is a subtype of MemoryChunkLonglivedArray,
     1228    for dealing with arguments to the wrapper's __call__ method.
     1229
     1230    Currently the __call__ method is declared to take a varargs
     1231    *args argument tuple.  We assume that the MemoryChunk named 'args'
     1232    will deal with that tuple.
     1233    """
     1234
     1235    def setup_args(self):
     1236        r"""
     1237        Handle the arguments of __call__ -- copy them into a pre-allocated
     1238        array, ready to pass to the interpreter.
     1239
     1240        EXAMPLES:
     1241            sage: from sage.ext.gen_interpreters import *
     1242            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     1243            sage: print mc.setup_args()
     1244            cdef mpfr_t* c_args = self._args
     1245            cdef int i
     1246            for i from 0 <= i < len(args):
     1247                rn = self.domain(args[i])
     1248                mpfr_set(self._args[i], rn.value, GMP_RNDN)
     1249            <BLANKLINE>
     1250        """
     1251        return je("""
     1252cdef {{ self.storage_type.c_ptr_type() }} c_args = self._args
     1253cdef int i
     1254for i from 0 <= i < len(args):
     1255    {{ self.storage_type.assign_c_from_py('self._args[i]', 'args[i]') | i(4) }}
     1256""", self=self)
     1257
     1258    def pass_argument(self):
     1259        r"""
     1260        Returns the string to pass the argument corresponding to this
     1261        memory chunk to the interpreter.
     1262
     1263        EXAMPLES:
     1264            sage: from sage.ext.gen_interpreters import *
     1265            sage: mc = MemoryChunkArguments('args', ty_mpfr)
     1266            sage: mc.pass_argument()
     1267            'c_args'
     1268        """
     1269        return 'c_args'   
     1270
     1271class MemoryChunkScratch(MemoryChunkLonglivedArray):
     1272    r"""
     1273    MemoryChunkScratch is a subtype of MemoryChunkLonglivedArray
     1274    for dealing with memory chunks that are allocated in the wrapper,
     1275    but only used in the interpreter -- stacks, scratch registers, etc.
     1276
     1277    (Currently these are only used as stacks.)
     1278    """
     1279
     1280    def __init__(self, name, storage_type, is_stack=False):
     1281        r"""
     1282        Initialize an instance of MemoryChunkScratch.
     1283
     1284        Initializes the _is_stack property, as well as
     1285        the properties described in the documentation for
     1286        MemoryChunk.__init__.
     1287
     1288        EXAMPLES:
     1289            sage: from sage.ext.gen_interpreters import *
     1290            sage: mc = MemoryChunkScratch('stack', ty_double, is_stack=True)
     1291            sage: mc.name
     1292            'stack'
     1293            sage: mc.storage_type is ty_double
     1294            True
     1295            sage: mc._is_stack
     1296            True
     1297        """
     1298        MemoryChunkLonglivedArray.__init__(self, name, storage_type)
     1299        self._is_stack = is_stack
     1300
     1301    def is_stack(self):
     1302        r"""
     1303        Says whether this memory chunk is a stack.  This affects code
     1304        generation for instructions using this memory chunk.
     1305
     1306        EXAMPLES:
     1307            sage: from sage.ext.gen_interpreters import *
     1308            sage: mc = MemoryChunkScratch('stack', ty_mpfr, is_stack=True)
     1309            sage: mc.is_stack()
     1310            True
     1311        """
     1312        return self._is_stack
     1313
     1314    def needs_cleanup_on_error(self):
     1315        r"""
     1316        In an interpreter that can terminate prematurely (due to an
     1317        exception from calling Python code, or divide by zero, or
     1318        whatever) it will just return at the end of the current instruction,
     1319        skipping the rest of the program.  Thus, it may still have
     1320        values pushed on the stack, etc.
     1321
     1322        This method returns True if this memory chunk is modified by the
     1323        interpreter and needs some sort of cleanup when an error happens.
     1324
     1325        EXAMPLES:
     1326            sage: from sage.ext.gen_interpreters import *
     1327            sage: mc = MemoryChunkScratch('registers', ty_python)
     1328            sage: mc.needs_cleanup_on_error()
     1329            True
     1330        """
     1331        return self.storage_type.python_refcounted()
     1332
     1333    def handle_cleanup(self):
     1334        r"""
     1335        Handle the cleanup if the interpreter exits with an error.
     1336
     1337        For scratch/stack chunks that hold Python-refcounted values,
     1338        we assume that they are filled with NULL on every entry to the
     1339        interpreter.  If the interpreter exited with an error, it may
     1340        have left values in the chunk, so we need to go through
     1341        the chunk and Py_CLEAR it.
     1342
     1343        EXAMPLES:
     1344            sage: from sage.ext.gen_interpreters import *
     1345            sage: mc = MemoryChunkScratch('registers', ty_python)
     1346            sage: print mc.handle_cleanup()
     1347            for i in range(self._n_registers):
     1348                Py_CLEAR(self._registers[i])
     1349            <BLANKLINE>        """
     1350        # XXX This is a lot slower than it needs to be, because
     1351        # we don't have a "cdef int i" in scope here.
     1352        return je("""
     1353for i in range(self._n_{{ self.name }}):
     1354    Py_CLEAR(self._{{ self.name }}[i])
     1355""", self=self)
     1356
     1357class MemoryChunkRRRetval(MemoryChunk):
     1358    r"""
     1359    A special-purpose memory chunk, for dealing with the return value
     1360    of the RR-based interpreter.
     1361    """
     1362
     1363    def declare_class_members(self):
     1364        r"""
     1365        Returns a string giving the declarations of the class members
     1366        in a wrapper class for this memory chunk.
     1367
     1368        EXAMPLES:
     1369            sage: from sage.ext.gen_interpreters import *
     1370            sage: mc = MemoryChunkRRRetval('retval', ty_mpfr)
     1371            sage: mc.declare_class_members()
     1372            ''
     1373        """
     1374        return ""
     1375
     1376    def declare_call_locals(self):
     1377        r"""
     1378        Returns a string to put in the __call__ method of a wrapper
     1379        class using this memory chunk, to allocate local variables.
     1380       
     1381        EXAMPLES:
     1382            sage: from sage.ext.gen_interpreters import *
     1383            sage: mc = MemoryChunkRRRetval('retval', ty_mpfr)
     1384            sage: mc.declare_call_locals()
     1385            u'        cdef RealNumber retval = (self.domain)()\n'
     1386        """
     1387        return je("""
     1388        cdef RealNumber {{ self.name }} = (self.domain)()
     1389""", self=self)
     1390       
     1391    def pass_argument(self):
     1392        r"""
     1393        Returns the string to pass the argument corresponding to this
     1394        memory chunk to the interpreter.
     1395
     1396        EXAMPLES:
     1397            sage: from sage.ext.gen_interpreters import *
     1398            sage: mc = MemoryChunkRRRetval('retval', ty_mpfr)
     1399            sage: mc.pass_argument()
     1400            u'&retval.value'
     1401        """
     1402        return je("""&{{ self.name }}.value""", self=self)
     1403
     1404class MemoryChunkPythonArguments(MemoryChunk):
     1405    r"""
     1406    A special-purpose memory chunk, for the generic Python-object based
     1407    interpreter.  Rather than copy the arguments into an array allocated
     1408    in the wrapper, we use the PyTupleObject internals and pass the array
     1409    that's inside the argument tuple.
     1410    """
     1411
     1412    def declare_class_members(self):
     1413        r"""
     1414        Returns a string giving the declarations of the class members
     1415        in a wrapper class for this memory chunk.
     1416
     1417        EXAMPLES:
     1418            sage: from sage.ext.gen_interpreters import *
     1419            sage: mc = MemoryChunkPythonArguments('args', ty_python)
     1420        """
     1421        return "    cdef int _n_%s\n" % self.name
     1422
     1423    def init_class_members(self):
     1424        r"""
     1425        Returns a string to be put in the __init__ method of a wrapper
     1426        class using this memory chunk, to initialize the corresponding
     1427        class members.
     1428
     1429        EXAMPLES:
     1430            sage: from sage.ext.gen_interpreters import *
     1431            sage: mc = MemoryChunkPythonArguments('args', ty_python)
     1432            sage: mc.init_class_members()
     1433            u"        count = args['args']\n        self._n_args = count\n"
     1434        """
     1435        return je("""
     1436        count = args['{{ self.name }}']
     1437        self._n_args = count
     1438""", self=self)
     1439
     1440    def setup_args(self):
     1441        r"""
     1442        Handle the arguments of __call__.  Nothing to do.
     1443
     1444        EXAMPLES:
     1445            sage: from sage.ext.gen_interpreters import *
     1446            sage: mc = MemoryChunkPythonArguments('args', ty_python)
     1447            sage: mc.setup_args()
     1448            ''
     1449        """
     1450        return ''
     1451
     1452    def pass_argument(self):
     1453        r"""
     1454        Pass the innards of the argument tuple to the interpreter.
     1455
     1456        EXAMPLES:
     1457            sage: from sage.ext.gen_interpreters import *
     1458            sage: mc = MemoryChunkPythonArguments('args', ty_python)
     1459            sage: mc.pass_argument()
     1460            '(<PyTupleObject*>args).ob_item'
     1461        """
     1462        return "(<PyTupleObject*>args).ob_item"
     1463
     1464class MemoryChunkElementArguments(MemoryChunkPythonArguments):
     1465    r"""
     1466    A special-purpose memory chunk, for the Python-object based
     1467    interpreters that want to process (and perhaps modify) the data.
     1468
     1469    We allocate a new list (via the map function) on every call to
     1470    hold the modified arguments.  That's not strictly necessary --
     1471    we could pre-allocate a list and map into it -- but this lets us
     1472    use simpler code for a very-likely-negligible efficiency cost.
     1473    (The Element interpreter is going to allocate lots of objects
     1474    as it runs, anyway.)
     1475    """
     1476
     1477    def setup_args(self):
     1478        r"""
     1479        Handle the arguments of __call__.  Note: This hardcodes
     1480        "self._domain".
     1481
     1482        EXAMPLES:
     1483            sage: from sage.ext.gen_interpreters import *
     1484            sage: mc = MemoryChunkElementArguments('args', ty_python)
     1485            sage: mc.setup_args()
     1486            'mapped_args = map(self._domain, args)\n'
     1487        """
     1488        return "mapped_args = map(self._domain, args)\n"
     1489
     1490    def pass_argument(self):
     1491        r"""
     1492        Pass the innards of the argument tuple to the interpreter.
     1493
     1494        EXAMPLES:
     1495            sage: from sage.ext.gen_interpreters import *
     1496            sage: mc = MemoryChunkElementArguments('args', ty_python)
     1497            sage: mc.pass_argument()
     1498            '(<PyListObject*>mapped_args).ob_item'
     1499        """
     1500        return "(<PyListObject*>mapped_args).ob_item"
     1501
     1502class MemoryChunkPyConstant(MemoryChunk):
     1503    r"""
     1504    A special-purpose memory chunk, for holding a single Python constant
     1505    and passing it to the interpreter as a PyObject*.
     1506    """
     1507
     1508    def __init__(self, name):
     1509        r"""
     1510        Initialize an instance of MemoryChunkPyConstant.
     1511
     1512        Always uses the type ty_python.
     1513
     1514        EXAMPLES:
     1515            sage: from sage.ext.gen_interpreters import *
     1516            sage: mc = MemoryChunkPyConstant('domain')
     1517            sage: mc.name
     1518            'domain'
     1519            sage: mc.storage_type is ty_python
     1520            True
     1521        """
     1522        MemoryChunk.__init__(self, name, ty_python)
     1523
     1524    def declare_class_members(self):
     1525        r"""
     1526        Returns a string giving the declarations of the class members
     1527        in a wrapper class for this memory chunk.
     1528
     1529        EXAMPLES:
     1530            sage: from sage.ext.gen_interpreters import *
     1531            sage: mc = MemoryChunkPyConstant('domain')
     1532            sage: mc.declare_class_members()
     1533            u'    cdef object _domain\n'
     1534        """
     1535        return je("""
     1536    cdef object _{{ self.name }}
     1537""", self=self)
     1538
     1539    def init_class_members(self):
     1540        r"""
     1541        Returns a string to be put in the __init__ method of a wrapper
     1542        class using this memory chunk, to initialize the corresponding
     1543        class members.
     1544
     1545        EXAMPLES:
     1546            sage: from sage.ext.gen_interpreters import *
     1547            sage: mc = MemoryChunkPyConstant('domain')
     1548            sage: mc.init_class_members()
     1549            u"        self._domain = args['domain']\n"
     1550        """
     1551        return je("""
     1552        self._{{ self.name }} = args['{{ self.name }}']
     1553""", self=self)
     1554
     1555    def declare_parameter(self):
     1556        r"""
     1557        Returns the string to use to declare the interpreter parameter
     1558        corresponding to this memory chunk.
     1559
     1560        EXAMPLES:
     1561            sage: from sage.ext.gen_interpreters import *
     1562            sage: mc = MemoryChunkPyConstant('domain')
     1563            sage: mc.declare_parameter()
     1564            'PyObject* domain'
     1565        """
     1566        return 'PyObject* %s' % self.name
     1567
     1568    def pass_argument(self):
     1569        r"""
     1570        Returns the string to pass the argument corresponding to this
     1571        memory chunk to the interpreter.
     1572
     1573        EXAMPLES:
     1574            sage: from sage.ext.gen_interpreters import *
     1575            sage: mc = MemoryChunkPyConstant('domain')
     1576            sage: mc.pass_argument()
     1577            '<PyObject*>self._domain'
     1578        """
     1579        return '<PyObject*>self._%s' % self.name
     1580
     1581def params_gen(**chunks):
     1582    r"""
     1583    Instructions have a parameter specification that says where they get
     1584    their inputs and where their outputs go.  Each parameter has
     1585    the same form: it is a triple (chunk, addr, len).  The chunk says
     1586    where the parameter is read from/written to.  The addr says which
     1587    value in the chunk is used.  If the chunk is a stack chunk, then
     1588    addr must be null; the value will be read from/written to the top
     1589    of the stack.  Otherwise, addr must be an integer, or another chunk;
     1590    if addr is another chunk, then the next value is read from that chunk
     1591    to be the address.
     1592
     1593    The len says how many values to read/write.  It can be either None
     1594    (meaning to read/write only a single value), an integer, or
     1595    another chunk; if it is a chunk, then the next value is read from that
     1596    chunk to be the len.  Note that specifying len changes the types
     1597    given to the instruction, so len==None is different than len==1 even
     1598    though both mean to use a single value.
     1599
     1600    These parameter specifications are cumbersome to write by hand, so
     1601    there's also a simple string format for them.  This (curried)
     1602    function parses the simple string format and produces parameter
     1603    specifications.  The params_gen function takes keyword arguments
     1604    mapping single-character names to memory chunks.  The string format
     1605    uses these names.  The params_gen function returns another function,
     1606    that takes two strings and returns a pair of lists of parameter
     1607    specifications.
     1608
     1609    Each string is the concatenation of arbitrarily many specifications.
     1610    Each specification consists of an address and a length.  The
     1611    address is either a single character naming a stack chunk,
     1612    or a string of the form 'A[B]' where A names a non-stack chunk
     1613    and B names the code chunk.  The length is either empty, or '@n'
     1614    for a number n (meaning to use that many arguments), or '@C', where
     1615    C is the code chunk.
     1616
     1617    EXAMPLES:
     1618        sage: from sage.ext.gen_interpreters import *
     1619        sage: mc_stack = MemoryChunkScratch('stack', ty_double, is_stack=True)
     1620        sage: mc_args = MemoryChunkArguments('args', ty_double)
     1621        sage: mc_code = MemoryChunkConstants('code', ty_int)
     1622
     1623        sage: pg = params_gen(D=mc_code, A=mc_args, S=mc_stack)
     1624        sage: pg('S', '')
     1625        ([({MC:stack}, None, None)], [])
     1626        sage: pg('A[D]', '')
     1627        ([({MC:args}, {MC:code}, None)], [])
     1628        sage: pg('S@5', '')
     1629        ([({MC:stack}, None, 5)], [])
     1630        sage: pg('S@D', '')
     1631        ([({MC:stack}, None, {MC:code})], [])
     1632        sage: pg('A[D]@D', '')
     1633        ([({MC:args}, {MC:code}, {MC:code})], [])
     1634        sage: pg('SSS@D', 'A[D]S@D')
     1635        ([({MC:stack}, None, None), ({MC:stack}, None, None), ({MC:stack}, None, {MC:code})], [({MC:args}, {MC:code}, None), ({MC:stack}, None, {MC:code})])
     1636    """
     1637   
     1638    def make_params(s):
     1639        p = []
     1640        s = s.strip()
     1641        while s:
     1642            chunk_code = s[0]
     1643            s = s[1:]
     1644            chunk = chunks[chunk_code]
     1645            addr = None
     1646            ch_len = None
     1647            if chunk.is_stack():
     1648                pass
     1649            else:
     1650                m = re.match(r'\[(?:([0-9]+)|([a-zA-Z]))\]', s)
     1651                if m.group(1):
     1652                    addr = int(m.group(1))
     1653                else:
     1654                    ch = chunks[m.group(2)]
     1655                    assert ch.storage_type is ty_int
     1656                    addr = ch
     1657                s = s[m.end():].strip()
     1658            if len(s) and s[0] == '@':
     1659                m = re.match(r'@(?:([0-9]+)|([a-zA-Z]))', s)
     1660                if m.group(1):
     1661                    ch_len = int(m.group(1))
     1662                else:
     1663                    ch = chunks[m.group(2)]
     1664                    assert ch.storage_type is ty_int
     1665                    ch_len = ch
     1666                s = s[m.end():].strip()
     1667            p.append((chunk, addr, ch_len))
     1668        return p
     1669
     1670    def params(s_ins, s_outs):
     1671        ins = make_params(s_ins)
     1672        outs = make_params(s_outs)
     1673        return (ins, outs)
     1674
     1675    return params
     1676
     1677def string_of_addr(a):
     1678    r"""
     1679    An address or a length from a parameter specification may be
     1680    either None, an integer, or a MemoryChunk.  If the address or
     1681    length is an integer or a MemoryChunk, this function will convert
     1682    it to a string giving an expression that will evaluate to the correct
     1683    address or length.  (See the docstring for params_gen for more
     1684    information on parameter specifications.)
     1685
     1686    EXAMPLES:
     1687        sage: from sage.ext.gen_interpreters import *
     1688        sage: mc_code = MemoryChunkConstants('code', ty_int)
     1689        sage: string_of_addr(mc_code)
     1690        '*code++'
     1691        sage: string_of_addr(42r)
     1692        '42'
     1693    """
     1694    if isinstance(a, (int, long)):
     1695        return str(a)
     1696    assert(isinstance(a, MemoryChunk))
     1697    return '*%s++' % a.name
     1698
     1699class InstrSpec(object):
     1700    r"""
     1701    Each instruction in an interpreter is represented as an InstrSpec.
     1702    This contains all the information that we need to generate code
     1703    to interpret the instruction; it also is used to build the tables
     1704    that fast_callable uses, so this is the nexus point between
     1705    users of the interpreter (possibly pure Python) and the
     1706    generated C interpreter.
     1707
     1708    The underlying instructions are matched to the caller by name.
     1709    For instance, fast_callable assumes that if the interpreter has an
     1710    instruction named 'cos', then it will take a single argument,
     1711    return a single result, and implement the cos() function.
     1712   
     1713    The print representation of an instruction (which will probably
     1714    only be used when doctesting this file) consists of the name,
     1715    a simplified stack effect, and the code (truncated if it's long).
     1716    The stack effect has two parts, the input and the output, separated
     1717    by '->'; the input shows what will be popped from the stack,
     1718    the output what will be placed on the stack.  Each consists of
     1719    a sequence of 'S' and '*' characters, where 'S' refers to a single
     1720    argument and '*' refers to a variable number of arguments.
     1721
     1722    The code for an instruction is a small snippet of C code.  It has
     1723    available variables 'i0', 'i1', ..., 'o0', 'o1', ...; one variable
     1724    for each input and output; its job is to assign values to the output
     1725    variables, based on the values of the input variables.
     1726
     1727    Normally, in an interpreter that uses doubles, each of the input
     1728    and output variables will be a double.  If i0 actually represents
     1729    a variable number of arguments, then it will be a pointer to
     1730    double instead, and there will be another variable n_i0 giving
     1731    the actual number of arguments.
     1732
     1733    When instructions refer to auto-reference types, they actually
     1734    get a pointer to the data in its original location; it is
     1735    not copied into a local variable.  Mostly, this makes no difference,
     1736    but there is one potential problem to be aware of.  It is possible
     1737    for an output variable to point to the same object as an input
     1738    variable; in fact, this usually will happen when you're working
     1739    with the stack.  If the instruction maps to a single function call,
     1740    then this is fine; the standard auto-reference implementations
     1741    (GMP, MPFR, etc.) are careful to allow having the input and output
     1742    be the same.  But if the instruction maps to multiple function
     1743    calls, you may need to use a temporary variable.
     1744
     1745    Here's an example of this issue.  Suppose you want to make an
     1746    instruction that does ``out = a+b*c``.  You write code like this:
     1747        out = b*c
     1748        out = a+out
     1749    But out will actually share the same storage as a; so the first line
     1750    modifies a, and you actually end up computing 2*(b+c).  The fix
     1751    is to only write to the output once, at the very end of your
     1752    instruction.
     1753
     1754    Instructions are also allowed to access memory chunks (other than
     1755    the stack and code) directly.  They are available as C variables
     1756    with the same name as the chunk.  This is useful if some type of
     1757    memory chunk doesn't fit well with the params_gen interface.
     1758
     1759    There are additional reference-counting rules that must be
     1760    followed if your interpreter operates on Python objects; these
     1761    rules are described in the docstring of the PythonInterpreter
     1762    class.
     1763
     1764    EXAMPLES:
     1765        sage: from sage.ext.gen_interpreters import *
     1766        sage: pg = RDFInterpreter().pg
     1767        sage: InstrSpec('add', pg('SS','S'), code='o0 = i0+i1;')
     1768        add: SS->S = 'o0 = i0+i1;'
     1769    """
     1770
     1771    def __init__(self, name, io, code=None, uses_error_handler=False, handles_own_decref=False):
     1772        r"""
     1773        Initialize an InstrSpec.
     1774
     1775        INPUTS:
     1776            name -- the name of the instruction
     1777            io -- a pair of lists of parameter specifications for I/O of the
     1778                  instruction
     1779            code -- a string containing a snippet of C code to read
     1780                    from the input variables and write to the output variables
     1781            uses_error_handler -- True if the instruction calls Python
     1782                                  and jumps to error: on a Python error
     1783            handles_own_decref -- True if the instruction handles Python
     1784                                  objects and includes its own
     1785                                  reference-counting
     1786
     1787        EXAMPLES:
     1788            sage: from sage.ext.gen_interpreters import *
     1789           
     1790            sage: pg = RDFInterpreter().pg
     1791            sage: InstrSpec('add', pg('SS','S'), code='o0 = i0+i1;')
     1792            add: SS->S = 'o0 = i0+i1;'
     1793            sage: instr = InstrSpec('py_call', pg('P[D]S@D', 'S'), code=('This is very complicated.  ' + 'blah ' * 30)); instr
     1794            py_call: *->S = 'This is very compli... blah blah blah '
     1795            sage: instr.name
     1796            'py_call'
     1797            sage: instr.inputs
     1798            [({MC:py_constants}, {MC:code}, None), ({MC:stack}, None, {MC:code})]
     1799            sage: instr.outputs
     1800            [({MC:stack}, None, None)]
     1801            sage: instr.code
     1802            'This is very complicated.  blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah blah '
     1803            sage: instr.parameters
     1804            ['py_constants', 'n_inputs']
     1805            sage: instr.n_inputs
     1806            0
     1807            sage: instr.n_outputs
     1808            1
     1809        """
     1810        self.name = name
     1811        self.inputs = io[0]
     1812        self.outputs = io[1]
     1813        self.uses_error_handler = uses_error_handler
     1814        self.handles_own_decref = handles_own_decref
     1815        if code is not None:
     1816            self.code = code
     1817        # XXX We assume that there is only one stack
     1818        n_inputs = 0
     1819        n_outputs = 0
     1820        in_effect = ''
     1821        out_effect = ''
     1822        p = []
     1823        for (ch, addr, len) in self.inputs:
     1824            if ch.is_stack():
     1825                if len is None:
     1826                    n_inputs += 1
     1827                    in_effect += 'S'
     1828                elif isinstance(len, (int, long)):
     1829                    n_inputs += len
     1830                    in_effect += 'S%d' % len
     1831                else:
     1832                    p.append('n_inputs')
     1833                    in_effect += '*'
     1834            else:
     1835                p.append(ch.name)
     1836        for (ch, addr, len) in self.outputs:
     1837            if ch.is_stack():
     1838                if len is None:
     1839                    n_outputs += 1
     1840                    out_effect += 'S'
     1841                elif isinstance(len, (int, long)):
     1842                    n_outputs += len
     1843                    out_effect += 'S%d' % len
     1844                else:
     1845                    p.append('n_outputs')
     1846                    out_effect += '*'
     1847            else:
     1848                p.append(ch.name)
     1849        self.parameters = p
     1850        self.n_inputs = n_inputs
     1851        self.n_outputs = n_outputs
     1852        self.in_effect = in_effect
     1853        self.out_effect = out_effect
     1854
     1855    def __repr__(self):
     1856        r"""
     1857        Produce a string representing a given instruction, consisting
     1858        of its name, a brief stack specification, and its code
     1859        (possibly abbreviated).
     1860
     1861        EXAMPLES:
     1862            sage: from sage.ext.gen_interpreters import *
     1863            sage: pg = RDFInterpreter().pg
     1864            sage: InstrSpec('add', pg('SS','S'), code='o0 = i0+i1;')
     1865            add: SS->S = 'o0 = i0+i1;'
     1866        """
     1867        rcode = repr(self.code)
     1868        if len(rcode) > 40:
     1869            rcode = rcode[:20] + '...' + rcode[-17:]
     1870        return '%s: %s->%s = %s' % \
     1871            (self.name, self.in_effect, self.out_effect, rcode)
     1872
     1873# Now we have a series of helper functions that make it slightly easier
     1874# to create instructions.
     1875
     1876def instr_infix(name, io, op):
     1877    r"""
     1878    A helper function for creating instructions implemented by
     1879    a single infix binary operator.
     1880
     1881    EXAMPLES:
     1882        sage: from sage.ext.gen_interpreters import *
     1883        sage: pg = RDFInterpreter().pg
     1884        sage: instr_infix('mul', pg('SS', 'S'), '*')
     1885        mul: SS->S = 'o0 = i0 * i1;'
     1886    """
     1887    return InstrSpec(name, io, code='o0 = i0 %s i1;' % op)
     1888
     1889def instr_funcall_2args(name, io, op):
     1890    r"""
     1891    A helper function for creating instructions implemented by
     1892    a two-argument function call.
     1893
     1894    EXAMPLES:
     1895        sage: from sage.ext.gen_interpreters import *
     1896        sage: pg = RDFInterpreter().pg
     1897        sage: instr_funcall_2args('atan2', pg('SS', 'S'), 'atan2')
     1898        atan2: SS->S = 'o0 = atan2(i0, i1);'
     1899    """
     1900    return InstrSpec(name, io, code='o0 = %s(i0, i1);' % op)
     1901
     1902def instr_unary(name, io, op):
     1903    r"""
     1904    A helper function for creating instructions with one input
     1905    and one output.
     1906
     1907    EXAMPLES:
     1908        sage: from sage.ext.gen_interpreters import *
     1909        sage: pg = RDFInterpreter().pg
     1910        sage: instr_unary('sin', pg('S','S'), 'sin(i0)')
     1911        sin: S->S = 'o0 = sin(i0);'
     1912        sage: instr_unary('neg', pg('S','S'), '-i0')
     1913        neg: S->S = 'o0 = -i0;'
     1914    """
     1915    return InstrSpec(name, io, code='o0 = ' + op + ';')
     1916
     1917def instr_funcall_2args_mpfr(name, io, op):
     1918    r"""
     1919    A helper function for creating MPFR instructions with two inputs
     1920    and one output.
     1921
     1922    EXAMPLES:
     1923        sage: from sage.ext.gen_interpreters import *
     1924        sage: pg = RRInterpreter().pg
     1925        sage: instr_funcall_2args_mpfr('add', pg('SS','S'), 'mpfr_add')
     1926        add: SS->S = 'mpfr_add(o0, i0, i1, GMP_RNDN);'
     1927    """
     1928    return InstrSpec(name, io, code='%s(o0, i0, i1, GMP_RNDN);' % op)
     1929
     1930def instr_funcall_1arg_mpfr(name, io, op):
     1931    r"""
     1932    A helper function for creating MPFR instructions with one input
     1933    and one output.
     1934
     1935    EXAMPLES:
     1936        sage: from sage.ext.gen_interpreters import *
     1937        sage: pg = RRInterpreter().pg
     1938        sage: instr_funcall_1arg_mpfr('exp', pg('S','S'), 'mpfr_exp')
     1939        exp: S->S = 'mpfr_exp(o0, i0, GMP_RNDN);'
     1940    """
     1941    return InstrSpec(name, io, code='%s(o0, i0, GMP_RNDN);' % op)
     1942
     1943class InterpreterSpec(object):
     1944    r"""
     1945    Each interpreter to be generated by this module is represented
     1946    by an InterpreterSpec.
     1947    """
     1948
     1949    def __init__(self):
     1950        r"""
     1951        Initialize an InterpreterSpec.
     1952
     1953        Initializes the following fields:
     1954
     1955        header -- a code snippet to go at the top of the C interpreter
     1956                  source file
     1957        pyx_header -- a code snippet to go at the top of the wrapper
     1958                      class source file
     1959        err_return -- a string indicating the value to be returned
     1960                      in case of a Python exception
     1961        mc_code -- a memory chunk to use for the interpreted code
     1962                     
     1963        EXAMPLES:
     1964            sage: from sage.ext.gen_interpreters import *
     1965            sage: interp = RDFInterpreter()
     1966            sage: interp.header
     1967            ''
     1968            sage: interp.pyx_header
     1969            ''
     1970            sage: interp.err_return
     1971            '-1094648009105371'
     1972            sage: interp.mc_code
     1973            {MC:code}
     1974        """
     1975        self.header = ''
     1976        self.pyx_header = ''
     1977        self.err_return = 'NULL'
     1978        self.mc_code = MemoryChunkConstants('code', ty_int)
     1979
     1980    def _set_opcodes(self):
     1981        r"""
     1982        Assign opcodes to the instructions in this interpreter.
     1983
     1984        Must be called at the end of __init__ by any subclass of
     1985        InterpreterSpec.
     1986
     1987        EXAMPLES:
     1988            sage: from sage.ext.gen_interpreters import *
     1989            sage: interp = RDFInterpreter()
     1990            sage: interp.instr_descs[5].opcode
     1991            5
     1992        """
     1993        for i in range(len(self.instr_descs)):
     1994            self.instr_descs[i].opcode = i
     1995
     1996
     1997class StackInterpreter(InterpreterSpec):
     1998    r"""
     1999    A subclass of InterpreterSpec, specialized for stack-based
     2000    interpreters.  (Currently all interpreters are stack-based.)
     2001    """
     2002
     2003    def __init__(self, type, mc_retval=None):
     2004        r"""
     2005        Initialize a StackInterpreter.
     2006
     2007        INPUTS:
     2008            type -- A StorageType; the basic type that this interpreter
     2009                    operates on
     2010            mc_retval -- default None; if not None, a special-purpose
     2011                         MemoryChunk to use as a return value
     2012
     2013        Initializes the fields described in the documentation for
     2014        InterpreterSpec.__init__, as well as the following:
     2015
     2016        mc_args, mc_constants, mc_stack -- MemoryChunk values
     2017        return_type -- the type returned by the C interpreter (None for int,
     2018                       where 1 means success and 0 means error)
     2019        mc_retval -- None, or the MemoryChunk to use as a return value
     2020
     2021        EXAMPLES:
     2022            sage: from sage.ext.gen_interpreters import *
     2023            sage: rdf = RDFInterpreter()
     2024            sage: rr = RRInterpreter()
     2025            sage: rdf.mc_args
     2026            {MC:args}
     2027            sage: rdf.mc_constants
     2028            {MC:constants}
     2029            sage: rdf.mc_stack
     2030            {MC:stack}
     2031            sage: rr.mc_retval
     2032            {MC:retval}
     2033            sage: rr.return_type is None
     2034            True
     2035            sage: rdf.return_type.type
     2036            'double'
     2037        """
     2038        InterpreterSpec.__init__(self)
     2039        self.mc_args = MemoryChunkArguments('args', type)
     2040        self.mc_constants = MemoryChunkConstants('constants', type)
     2041        self.mc_stack = MemoryChunkScratch('stack', type, is_stack=True)
     2042        if isinstance(type, StorageTypeAssignable):
     2043            self.return_type = type
     2044        else:
     2045            self.return_type = None
     2046        self.mc_retval = mc_retval
     2047
     2048class RDFInterpreter(StackInterpreter):
     2049    r"""
     2050    A subclass of StackInterpreter, specifying an interpreter over
     2051    machine-floating-point values (C doubles).
     2052    """
     2053
     2054    def __init__(self):
     2055        r"""
     2056        Initialize an RDFInterpreter.
     2057
     2058        EXAMPLES:
     2059            sage: from sage.ext.gen_interpreters import *
     2060            sage: interp = RDFInterpreter()
     2061            sage: interp.name
     2062            'rdf'
     2063            sage: interp.mc_py_constants
     2064            {MC:py_constants}
     2065            sage: interp.chunks
     2066            [{MC:args}, {MC:constants}, {MC:py_constants}, {MC:stack}, {MC:code}]
     2067            sage: interp.pg('A[D]', 'S')
     2068            ([({MC:args}, {MC:code}, None)], [({MC:stack}, None, None)])
     2069            sage: instrs = dict([(ins.name, ins) for ins in interp.instr_descs])
     2070            sage: instrs['add']
     2071            add: SS->S = 'o0 = i0 + i1;'
     2072            sage: instrs['py_call']
     2073            py_call: *->S = ' \nPyObject *py_arg...goto error;\n}\n'
     2074        """
     2075
     2076        StackInterpreter.__init__(self, ty_double)
     2077        self.name = 'rdf'
     2078        self.mc_py_constants = MemoryChunkConstants('py_constants', ty_python)
     2079        # This is a randomly chosen number.  Whenever this number is
     2080        # returned, the wrapper has to check whether an exception actually
     2081        # happened, so if an expression evaluates to this number execution
     2082        # is slightly slower.  Hopefully that won't happen too often :)
     2083        self.err_return = '-1094648009105371'
     2084        self.chunks = [self.mc_args, self.mc_constants, self.mc_py_constants,
     2085                       self.mc_stack,
     2086                       self.mc_code]
     2087        pg = params_gen(A=self.mc_args, C=self.mc_constants, D=self.mc_code,
     2088                        S=self.mc_stack, P=self.mc_py_constants)
     2089        self.pg = pg
     2090        instrs = [
     2091            InstrSpec('load_arg', pg('A[D]', 'S'),
     2092                       code='o0 = i0;'),
     2093            InstrSpec('load_const', pg('C[D]', 'S'),
     2094                       code='o0 = i0;'),
     2095            InstrSpec('return', pg('S', ''),
     2096                       code='return i0;'),
     2097            InstrSpec('py_call', pg('P[D]S@D', 'S'),
     2098                       uses_error_handler=True,
     2099                       code="""
     2100PyObject *py_args = PyTuple_New(n_i1);
     2101if (py_args == NULL) goto error;
     2102int i;
     2103for (i = 0; i < n_i1; i++) {
     2104  PyObject *arg = PyFloat_FromDouble(i1[i]);
     2105  if (arg == NULL) {
     2106    Py_DECREF(py_args);
     2107    goto error;
     2108  }
     2109  PyTuple_SET_ITEM(py_args, i, arg);
     2110}
     2111PyObject *result = PyObject_CallObject(i0, py_args);
     2112Py_DECREF(py_args);
     2113if (result == NULL) goto error;
     2114/* If result is not a float, then this will turn it into a float first. */
     2115o0 = PyFloat_AsDouble(result);
     2116Py_DECREF(result);
     2117if (o0 == -1 && PyErr_Occurred()) {
     2118  goto error;
     2119}
     2120""")
     2121            ]
     2122        for (name, op) in [('add', '+'), ('sub', '-'),
     2123                           ('mul', '*'), ('div', '/')]:
     2124            instrs.append(instr_infix(name, pg('SS', 'S'), op))
     2125        instrs.append(instr_funcall_2args('pow', pg('SS', 'S'), 'pow'))
     2126        for (name, op) in [('neg', '-i0'), ('invert', '1/i0'),
     2127                           ('abs', 'fabs(i0)')]:
     2128            instrs.append(instr_unary(name, pg('S', 'S'), op))
     2129        for name in ['sqrt', 'ceil', 'floor', 'sin', 'cos', 'tan',
     2130                     'asin', 'acos', 'atan', 'sinh', 'cosh', 'tanh',
     2131                     'asinh', 'acosh', 'atanh', 'exp', 'log']:
     2132            instrs.append(instr_unary(name, pg('S',  'S'), "%s(i0)" % name))
     2133        self.instr_descs = instrs
     2134        self._set_opcodes()
     2135       
     2136class RRInterpreter(StackInterpreter):
     2137    r"""
     2138    A subclass of StackInterpreter, specifying an interpreter over
     2139    MPFR arbitrary-precision floating-point numbers.
     2140    """
     2141   
     2142    def __init__(self):
     2143        r"""
     2144        Initialize an RDFInterpreter.
     2145
     2146        EXAMPLES:
     2147            sage: from sage.ext.gen_interpreters import *
     2148            sage: interp = RRInterpreter()
     2149            sage: interp.name
     2150            'rr'
     2151            sage: interp.mc_py_constants
     2152            {MC:py_constants}
     2153            sage: interp.chunks
     2154            [{MC:args}, {MC:retval}, {MC:constants}, {MC:py_constants}, {MC:stack}, {MC:code}, {MC:domain}]
     2155            sage: interp.pg('A[D]', 'S')
     2156            ([({MC:args}, {MC:code}, None)], [({MC:stack}, None, None)])
     2157            sage: instrs = dict([(ins.name, ins) for ins in interp.instr_descs])
     2158            sage: instrs['add']
     2159            add: SS->S = 'mpfr_add(o0, i0, i1, GMP_RNDN);'
     2160            sage: instrs['py_call']
     2161            py_call: *->S = '\nif (!rr_py_call_h...goto error;\n}\n'
     2162
     2163        That py_call instruction is particularly interesting, and
     2164        demonstrates a useful technique to let you use Cython code
     2165        in an interpreter.  Let's look more closely:
     2166
     2167            sage: print instrs['py_call'].code
     2168            if (!rr_py_call_helper(domain, i0, n_i1, i1, o0)) {
     2169              goto error;
     2170            }
     2171
     2172        This instruction makes use of the function rr_py_call_helper,
     2173        which is declared...
     2174
     2175            sage: print interp.header
     2176            #include <mpfr.h>
     2177            #include "wrapper_rr.h"
     2178
     2179        in wrapper_rr.h.  This file is automatically generated by Cython
     2180        from wrapper_rr.pyx, which is automatically generated from
     2181        this spec; in particular, rr_py_call_helper comes from:
     2182           
     2183            sage: print interp.pyx_header
     2184            from ...
     2185            cdef public bint rr_py_call_helper(object domain, object fn,
     2186                                               int n_args,
     2187                                               mpfr_t* args, mpfr_t* retval) except 0:
     2188                py_args = []
     2189                cdef int i
     2190                cdef RealNumber rn
     2191                for i from 0 <= i < n_args:
     2192                    rn = domain()
     2193                    mpfr_set(rn.value, args[i], GMP_RNDN)
     2194                    py_args.append(rn)
     2195                cdef RealNumber result = domain(fn(*py_args))
     2196                mpfr_set(retval[0], result.value, GMP_RNDN)
     2197                return 1
     2198
     2199
     2200        So instructions where you need to interact with Python can
     2201        call back into Cython code fairly easily.
     2202        """
     2203
     2204        StackInterpreter.__init__(self, ty_mpfr, mc_retval= MemoryChunkRRRetval('retval', ty_mpfr))
     2205        self.name = 'rr'
     2206        self.err_return = '0'
     2207        self.mc_py_constants = MemoryChunkConstants('py_constants', ty_python)
     2208        self.mc_domain = MemoryChunkPyConstant('domain')
     2209        self.chunks = [self.mc_args, self.mc_retval, self.mc_constants,
     2210                       self.mc_py_constants,
     2211                       self.mc_stack, self.mc_code, self.mc_domain]
     2212        pg = params_gen(A=self.mc_args, C=self.mc_constants, D=self.mc_code,
     2213                        S=self.mc_stack,
     2214                        P=self.mc_py_constants)
     2215        self.pg = pg
     2216        self.header = """
     2217#include <mpfr.h>
     2218#include "wrapper_rr.h"
     2219"""
     2220        self.pyx_header = """
     2221from sage.rings.real_mpfr cimport RealField, RealNumber
     2222from sage.libs.mpfr cimport *
     2223
     2224cdef public bint rr_py_call_helper(object domain, object fn,
     2225                                   int n_args,
     2226                                   mpfr_t* args, mpfr_t* retval) except 0:
     2227    py_args = []
     2228    cdef int i
     2229    cdef RealNumber rn
     2230    for i from 0 <= i < n_args:
     2231        rn = domain()
     2232        mpfr_set(rn.value, args[i], GMP_RNDN)
     2233        py_args.append(rn)
     2234    cdef RealNumber result = domain(fn(*py_args))
     2235    mpfr_set(retval[0], result.value, GMP_RNDN)
     2236    return 1
     2237
     2238"""[1:]
     2239        instrs = [
     2240            InstrSpec('load_arg', pg('A[D]', 'S'),
     2241                       code='mpfr_set(o0, i0, GMP_RNDN);'),
     2242            InstrSpec('load_const', pg('C[D]', 'S'),
     2243                       code='mpfr_set(o0, i0, GMP_RNDN);'),
     2244            InstrSpec('return', pg('S', ''),
     2245                       code='mpfr_set(retval[0], i0, GMP_RNDN);\nreturn 1;\n'),
     2246            InstrSpec('py_call', pg('P[D]S@D', 'S'),
     2247                       uses_error_handler=True,
     2248                       code="""
     2249if (!rr_py_call_helper(domain, i0, n_i1, i1, o0)) {
     2250  goto error;
     2251}
     2252""")
     2253            ]
     2254        for (name, op) in [('add', 'mpfr_add'), ('sub', 'mpfr_sub'),
     2255                           ('mul', 'mpfr_mul'), ('div', 'mpfr_div'),
     2256                           ('pow', 'mpfr_pow')]:
     2257            instrs.append(instr_funcall_2args_mpfr(name, pg('SS', 'S'), op))
     2258        for name in ['neg', 'abs',
     2259                     'log', 'log2', 'log10',
     2260                     'exp', 'exp2', 'exp10',
     2261                     'cos', 'sin', 'tan',
     2262                     'sec', 'csc', 'cot',
     2263                     'acos', 'asin', 'atan',
     2264                     'cosh', 'sinh', 'tanh',
     2265                     'sech', 'csch', 'coth',
     2266                     'acosh', 'asinh', 'atanh',
     2267                     'log1p', 'expm1', 'eint',
     2268                     'gamma', 'lngamma',
     2269                     'zeta', 'erf', 'erfc',
     2270                     'j0', 'j1', 'y0', 'y1']:
     2271            instrs.append(instr_funcall_1arg_mpfr(name, pg('S', 'S'), 'mpfr_' + name))
     2272        # mpfr_ui_div constructs a temporary mpfr_t and then calls mpfr_div;
     2273        # it would probably be (slightly) faster to use a permanent copy
     2274        # of "one" (on the other hand, the constructed temporary copy is
     2275        # on the stack, so it's very likely to be in the cache).
     2276        instrs.append(InstrSpec('invert', pg('S', 'S'),
     2277                                 code='mpfr_ui_div(o0, 1, i0, GMP_RNDN);'))
     2278        self.instr_descs = instrs
     2279        self._set_opcodes()
     2280
     2281class PythonInterpreter(StackInterpreter):
     2282    r"""
     2283    A subclass of StackInterpreter, specifying an interpreter over
     2284    Python objects.
     2285
     2286    Let's discuss how the reference-counting works in Python-object
     2287    based interpreters.
     2288
     2289    There is a simple rule to remember: when executing the code
     2290    snippets, the input variables contain borrowed references;
     2291    you must fill in the output variables with references you own.
     2292
     2293    As an optimization, an instruction may set .handles_own_decref; in
     2294    that case, it must decref any input variables that came from the
     2295    stack.  (Input variables that came from arguments/constants chunks
     2296    must NOT be decref'ed!)  In addition, with .handles_own_decref, if
     2297    any of your input variables are arbitrary-count, then you must
     2298    NULL out these variables as you decref them.  (Use Py_CLEAR to do
     2299    this, unless you understand the documentation of Py_CLEAR and why
     2300    it's different than Py_XDECREF followed by assigning NULL.)
     2301
     2302    Note that as a tiny optimization, the interpreter always assumes
     2303    (and assures) that empty parts of the stack contain NULL, so
     2304    it doesn't bother to Py_XDECREF before it pushes onto the stack.
     2305    """
     2306
     2307    def __init__(self):
     2308        r"""
     2309        Initialize a PythonInterpreter.
     2310
     2311        EXAMPLES:
     2312            sage: from sage.ext.gen_interpreters import *
     2313            sage: interp = PythonInterpreter()
     2314            sage: interp.name
     2315            'py'
     2316            sage: interp.mc_args
     2317            {MC:args}
     2318            sage: interp.chunks
     2319            [{MC:args}, {MC:constants}, {MC:stack}, {MC:code}]
     2320            sage: instrs = dict([(ins.name, ins) for ins in interp.instr_descs])
     2321            sage: instrs['add']
     2322            add: SS->S = 'o0 = PyNumber_Add(i0, i1);'
     2323            sage: instrs['py_call']
     2324            py_call: *->S = '\nPyObject *py_args...CREF(py_args);\n'
     2325        """
     2326
     2327        StackInterpreter.__init__(self, ty_python)
     2328        self.name = 'py'
     2329        # StackInterpreter.__init__ gave us a MemoryChunkArguments.
     2330        # Override with MemoryChunkPythonArguments.
     2331        self.mc_args = MemoryChunkPythonArguments('args', ty_python)
     2332        self.chunks = [self.mc_args, self.mc_constants, self.mc_stack,
     2333                       self.mc_code]
     2334        pg = params_gen(A=self.mc_args, C=self.mc_constants, D=self.mc_code,
     2335                        S=self.mc_stack)
     2336        self.pg = pg
     2337        self.header = """
     2338#include <Python.h>
     2339#define CHECK(x) (x != NULL)
     2340"""
     2341        instrs = [
     2342            InstrSpec('load_arg', pg('A[D]', 'S'),
     2343                       code='o0 = i0; Py_INCREF(o0);'),
     2344            InstrSpec('load_const', pg('C[D]', 'S'),
     2345                       code='o0 = i0; Py_INCREF(o0);'),
     2346            InstrSpec('return', pg('S', ''),
     2347                       code='return i0;',
     2348                       handles_own_decref=True),
     2349            InstrSpec('py_call', pg('C[D]S@D', 'S'),
     2350                       handles_own_decref=True,
     2351                       code="""
     2352PyObject *py_args = PyTuple_New(n_i1);
     2353if (py_args == NULL) goto error;
     2354int i;
     2355for (i = 0; i < n_i1; i++) {
     2356  PyObject *arg = i1[i];
     2357  PyTuple_SET_ITEM(py_args, i, arg);
     2358  i1[i] = NULL;
     2359}
     2360o0 = PyObject_CallObject(i0, py_args);
     2361Py_DECREF(py_args);
     2362""")
     2363            ]
     2364        for (name, op) in [('add', 'PyNumber_Add'),
     2365                           ('sub', 'PyNumber_Subtract'),
     2366                           ('mul', 'PyNumber_Multiply'),
     2367                           ('div', 'PyNumber_Divide')]:
     2368            instrs.append(instr_funcall_2args(name, pg('SS', 'S'), op))
     2369        instrs.append(InstrSpec('pow', pg('SS', 'S'),
     2370                                code='o0 = PyNumber_Power(i0, i1, Py_None);'))
     2371        for (name, op) in [('neg', 'PyNumber_Negative'),
     2372                           ('invert', 'PyNumber_Invert'),
     2373                           ('abs', 'PyNumber_Absolute')]:
     2374            instrs.append(instr_unary(name, pg('S', 'S'), '%s(i0)'%op))
     2375        self.instr_descs = instrs
     2376        self._set_opcodes()
     2377       
     2378class ElementInterpreter(PythonInterpreter):
     2379    r"""
     2380    A subclass of PythonInterpreter, specifying an interpreter over
     2381    Sage elements with a particular parent.
     2382
     2383    This is very similar to the PythonInterpreter, but after every
     2384    instruction, the result is checked to make sure it actually an
     2385    element with the correct parent; if not, we attempt to convert it.
     2386
     2387    Uses the same instructions (with the same implementation) as
     2388    PythonInterpreter.
     2389    """
     2390
     2391    def __init__(self):
     2392        r"""
     2393        Initialize an ElementInterpreter.
     2394
     2395        EXAMPLES:
     2396            sage: from sage.ext.gen_interpreters import *
     2397            sage: interp = ElementInterpreter()
     2398            sage: interp.name
     2399            'el'
     2400            sage: interp.mc_args
     2401            {MC:args}
     2402            sage: interp.chunks
     2403            [{MC:args}, {MC:constants}, {MC:stack}, {MC:domain}, {MC:code}]
     2404            sage: instrs = dict([(ins.name, ins) for ins in interp.instr_descs])
     2405            sage: instrs['add']
     2406            add: SS->S = 'o0 = PyNumber_Add(i0, i1);'
     2407            sage: instrs['py_call']
     2408            py_call: *->S = '\nPyObject *py_args...CREF(py_args);\n'
     2409        """
     2410
     2411        PythonInterpreter.__init__(self)
     2412        self.name = 'el'
     2413        # PythonInterpreter.__init__ gave us a MemoryChunkPythonArguments.
     2414        # Override with MemoryChunkElementArguments.
     2415        self.mc_args = MemoryChunkElementArguments('args', ty_python)
     2416        self.mc_domain_info = MemoryChunkPyConstant('domain')
     2417        self.chunks = [self.mc_args, self.mc_constants, self.mc_stack,
     2418                       self.mc_domain_info, self.mc_code]
     2419        self.header = """
     2420#include <Python.h>
     2421#include "wrapper_el.h"
     2422
     2423#define CHECK(x) do_check(&(x), domain)
     2424
     2425static inline int do_check(PyObject **x, PyObject *domain) {
     2426  if (*x == NULL) return 0;
     2427  PyObject *new_x = el_check_element(*x, domain);
     2428  Py_DECREF(*x);
     2429  *x = new_x;
     2430  if (*x == NULL) return 0;
     2431  return 1;
     2432}
     2433"""
     2434        self.pyx_header = """
     2435from sage.structure.element cimport Element
     2436
     2437cdef public object el_check_element(object v, parent):
     2438    cdef Element v_el
     2439
     2440    if PY_TYPE_CHECK(v, Element):
     2441        v_el = <Element>v
     2442        if v_el._parent is parent:
     2443            return v_el
     2444
     2445    return parent(v)
     2446
     2447"""[1:]
     2448
     2449class InterpreterGenerator(object):
     2450    r"""
     2451    This class takes an InterpreterSpec and generates the corresponding
     2452    C interpreter and Cython wrapper.
     2453
     2454    See the documentation for methods get_wrapper and get_interpreter
     2455    for more information.
     2456    """
     2457
     2458    def __init__(self, spec):
     2459        r"""
     2460        Initialize an InterpreterGenerator.
     2461
     2462        INPUTS:
     2463            spec -- an InterpreterSpec
     2464
     2465        EXAMPLES:
     2466            sage: from sage.ext.gen_interpreters import *
     2467            sage: interp = RDFInterpreter()
     2468            sage: gen = InterpreterGenerator(interp)
     2469            sage: gen._spec is interp
     2470            True
     2471            sage: gen.uses_error_handler
     2472            False
     2473        """
     2474
     2475        self._spec = spec
     2476        self.uses_error_handler = False
     2477
     2478    def gen_code(self, instr_desc, write):
     2479        r"""
     2480        Generates code for a single instruction.
     2481
     2482        INPUTS:
     2483            instr_desc -- an InstrSpec
     2484            write -- a Python callable
     2485
     2486        This function calls its write parameter successively with
     2487        strings; when these strings are concatenated, the result is
     2488        the code for the given instruction.
     2489
     2490        See the documentation for the get_interpreter method for more
     2491        information.
     2492
     2493        EXAMPLES:
     2494            sage: from sage.ext.gen_interpreters import *
     2495            sage: interp = RDFInterpreter()
     2496            sage: gen = InterpreterGenerator(interp)
     2497            sage: import cStringIO
     2498            sage: buff = cStringIO.StringIO()
     2499            sage: instrs = dict([(ins.name, ins) for ins in interp.instr_descs])
     2500            sage: gen.gen_code(instrs['div'], buff.write)
     2501            sage: print buff.getvalue()
     2502                case 7: /* div */
     2503                  {
     2504                    double i1 = *--stack;
     2505                    double i0 = *--stack;
     2506                    double o0;
     2507                    o0 = i0 / i1;
     2508                    *stack++ = o0;
     2509                  }
     2510                  break;
     2511            <BLANKLINE>
     2512        """
     2513
     2514        d = instr_desc
     2515        w = write
     2516        s = self._spec
     2517
     2518        if d.uses_error_handler:
     2519            self.uses_error_handler = True
     2520
     2521        w(je("""
     2522    case {{ d.opcode }}: /* {{ d.name }} */
     2523      {
     2524""", d=d))
     2525       
     2526        # If the inputs to an instruction come from the stack,
     2527        # then we want to generate code for the inputs in reverse order:
     2528        # for instance, the divide instruction, which takes inputs A and B
     2529        # and generates A/B, needs to pop B off the stack first.
     2530        # On the other hand, if the inputs come from the constant pool,
     2531        # then we want to generate code for the inputs in normal order,
     2532        # because the addresses in the code stream will be in that order.
     2533        # We handle this by running through the inputs in two passes:
     2534        # first a forward pass, where we handle non-stack inputs
     2535        # (and lengths for stack inputs), and then a reverse pass,
     2536        # where we handle stack inputs.
     2537        for i in range(len(d.inputs)):
     2538            (ch, addr, input_len) = d.inputs[i]
     2539            chst = ch.storage_type
     2540            if addr is not None:
     2541                w("        int ai%d = %s;\n" % (i, string_of_addr(addr)))
     2542            if input_len is not None:
     2543                w("        int n_i%d = %s;\n" % (i, string_of_addr(input_len)))
     2544            if not ch.is_stack():
     2545                if input_len is not None:
     2546                    w("        %s i%d = %s + ai%d;\n" %
     2547                      (chst.c_ptr_type(), i, ch.name, i))
     2548                else:
     2549                    w("        %s i%d = %s[ai%d];\n" %
     2550                      (chst.c_local_type(), i, ch.name, i))
     2551
     2552        for i in reversed(range(len(d.inputs))):
     2553            (ch, addr, input_len) = d.inputs[i]
     2554            chst = ch.storage_type
     2555            if ch.is_stack():
     2556                if input_len is not None:
     2557                    w("        %s -= n_i%d;\n" % (ch.name, i))
     2558                    w("        %s i%d = %s;\n" % (chst.c_ptr_type(), i, ch.name))
     2559                else:
     2560                    w("        %s i%d = *--%s;\n" % (chst.c_local_type(), i, ch.name))
     2561                    if ch.is_python_refcounted_stack():
     2562                        w("        *%s = NULL;\n" % ch.name)
     2563
     2564        for i in range(len(d.outputs)):
     2565            (ch, addr, output_len) = d.outputs[i]
     2566            chst = ch.storage_type
     2567            if addr is not None:
     2568                w("        int ao%d = %s;\n" % (i, string_of_addr(addr)))
     2569            if output_len is not None:
     2570                w("        int n_o%d = %s;\n" % (i, string_of_addr(output_len)))
     2571                if ch.is_stack():
     2572                    w("        %s o%d = %s;\n" %
     2573                      (chst.c_ptr_type(), i, ch.name))
     2574                    w("        %s += n_o%d;\n" % (ch.name, i))
     2575                else:
     2576                    w("        %s o%d = %s + ao%d;\n" %
     2577                      (chst.c_ptr_type(), i, ch.name, i))
     2578                     
     2579            else:
     2580                if not chst.cheap_copies():
     2581                    if ch.is_stack():
     2582                        w("        %s o%d = *%s++;\n" %
     2583                          (chst.c_local_type(), i, ch.name))
     2584                    else:
     2585                        w("        %s o%d = %s[ao%d];\n" %
     2586                          (chst.c_local_type(), i, ch.name, i))
     2587                else:
     2588                    w("        %s o%d;\n" % (chst.c_local_type(), i))
     2589               
     2590        w(indent_lines(8, d.code.rstrip('\n') + '\n'))
     2591
     2592        stack_offsets = defaultdict(int)
     2593        for i in range(len(d.inputs)):
     2594            (ch, addr, input_len) = d.inputs[i]
     2595            chst = ch.storage_type
     2596            if ch.is_python_refcounted_stack() and not d.handles_own_decref:
     2597                if input_len is None:
     2598                    w("        Py_DECREF(i%d);\n" % i)
     2599                    stack_offsets[ch] += 1
     2600                else:
     2601                    w(je("""
     2602        int {{ iter }};
     2603        for ({{ iter }} = 0; {{ iter }} < n_i{{ i }}; {{ iter }}++) {
     2604          Py_CLEAR(i{{ i }}[{{ iter }}]);
     2605        }
     2606""", iter='_interp_iter_%d' % i, i=i))
     2607
     2608        for i in range(len(d.outputs)):
     2609            ch = d.outputs[i][0]
     2610            chst = ch.storage_type
     2611            if chst.python_refcounted():
     2612                # We don't yet support code chunks
     2613                # that produce multiple Python values, because of
     2614                # the way it complicates error handling.
     2615                assert i == 0
     2616                w("        if (!CHECK(o%d)) {\n" % i)
     2617                w("          Py_XDECREF(o%d);\n" % i)
     2618                w("          goto error;\n")
     2619                w("        }\n")
     2620                self.uses_error_handler = True
     2621            if chst.cheap_copies():
     2622                if ch.is_stack():
     2623                    w("        *%s++ = o%d;\n" % (ch.name, i))
     2624                else:
     2625                    w("        %s[ao%d] = o%d;\n" % (ch.name, i, i))
     2626
     2627        w(je("""
     2628      }
     2629      break;
     2630"""))
     2631
     2632    def func_header(self, cython=False):
     2633        r"""
     2634        Generates the function header for the declaration (in the Cython
     2635        wrapper) or the definition (in the C interpreter) of the interpreter
     2636        function.
     2637
     2638        EXAMPLES:
     2639            sage: from sage.ext.gen_interpreters import *
     2640            sage: interp = ElementInterpreter()
     2641            sage: gen = InterpreterGenerator(interp)
     2642            sage: print gen.func_header()
     2643            PyObject* interp_el(PyObject** args,
     2644                    PyObject** constants,
     2645                    PyObject** stack,
     2646                    PyObject* domain,
     2647                    int* code)
     2648            sage: print gen.func_header(cython=True)
     2649            object interp_el(PyObject** args,
     2650                    PyObject** constants,
     2651                    PyObject** stack,
     2652                    PyObject* domain,
     2653                    int* code)
     2654        """
     2655        s = self._spec
     2656        ret_ty = 'bint' if cython else 'int'
     2657        if s.return_type:
     2658            ret_ty = s.return_type.c_decl_type()
     2659            if cython:
     2660                ret_ty = s.return_type.cython_decl_type()
     2661        return je("""{{ ret_ty }} interp_{{ s.name }}(
     2662{%- for ch in s.chunks %}
     2663{%    if not loop.first %},
     2664        {% endif %}{{ ch.declare_parameter() }}
     2665{%- endfor %})""", ret_ty=ret_ty, s=s)
     2666
     2667    def write_interpreter(self, write):
     2668        r"""
     2669        Generate the code for the C interpreter.
     2670
     2671        This function calls its write parameter successively with
     2672        strings; when these strings are concatenated, the result is
     2673        the code for the interpreter.
     2674
     2675        See the documentation for the get_interpreter method for more
     2676        information.
     2677
     2678        EXAMPLES:
     2679            sage: from sage.ext.gen_interpreters import *
     2680            sage: interp = RDFInterpreter()
     2681            sage: gen = InterpreterGenerator(interp)
     2682            sage: import cStringIO
     2683            sage: buff = cStringIO.StringIO()
     2684            sage: gen.write_interpreter(buff.write)
     2685            sage: print buff.getvalue()
     2686            /* Automatically generated.  Do not edit! */ ...
     2687        """
     2688        s = self._spec
     2689        w = write
     2690        w(je("""
     2691/* Automatically generated.  Do not edit! */
     2692#include <Python.h>
     2693{% print s.header %}
     2694{{ self.func_header() }} {
     2695  while (1) {
     2696    switch (*code++) {
     2697""", s=s, self=self, i=indent_lines))
     2698        for instr_desc in s.instr_descs:
     2699            self.gen_code(instr_desc, w)
     2700        w(je("""
     2701    }
     2702  }
     2703{% if self.uses_error_handler %}
     2704error:
     2705  return {{ s.err_return }};
     2706{% endif %}
     2707}
     2708
     2709""", s=s, i=indent_lines, self=self))
     2710
     2711    def write_wrapper(self, write):
     2712        r"""
     2713        Generate the code for the Cython wrapper.
     2714        This function calls its write parameter successively with
     2715        strings; when these strings are concatenated, the result is
     2716        the code for the interpreter.
     2717
     2718        See the documentation for the get_wrapper method for more
     2719        information.
     2720
     2721        EXAMPLES:
     2722            sage: from sage.ext.gen_interpreters import *
     2723            sage: interp = RDFInterpreter()
     2724            sage: gen = InterpreterGenerator(interp)
     2725            sage: import cStringIO
     2726            sage: buff = cStringIO.StringIO()
     2727            sage: gen.write_wrapper(buff.write)
     2728            sage: print buff.getvalue()
     2729            # Automatically generated.  Do not edit! ...
     2730        """
     2731        s = self._spec
     2732        w = write
     2733        types = set()
     2734        do_cleanup = False
     2735        for ch in s.chunks:
     2736            if ch.storage_type is not None:
     2737                types.add(ch.storage_type)
     2738            do_cleanup = do_cleanup or ch.needs_cleanup_on_error()
     2739        for ch in s.chunks:
     2740            if ch.name == 'args':
     2741                arg_ch = ch
     2742
     2743        the_call = je("""
     2744        {% if s.return_type %}return {% endif -%}
     2745interp_{{ s.name }}({{ arg_ch.pass_argument() }}
     2746{% for ch in s.chunks[1:] %}
     2747            , {{ ch.pass_argument() }}
     2748{% endfor %}
     2749            )
     2750""", s=s, arg_ch=arg_ch)
     2751
     2752        w(je("""
     2753# Automatically generated.  Do not edit!
     2754
     2755include "../stdsage.pxi"
     2756from python_object cimport PyObject
     2757cdef extern from "Python.h":
     2758    void Py_DECREF(PyObject *o)
     2759    void Py_INCREF(PyObject *o)
     2760    void Py_CLEAR(PyObject *o)
     2761
     2762cdef extern from "listobject.h":
     2763    object PyList_New(Py_ssize_t len)
     2764    ctypedef struct PyListObject:
     2765        PyObject **ob_item
     2766
     2767cdef extern from "tupleobject.h":
     2768    ctypedef struct PyTupleObject:
     2769        PyObject **ob_item
     2770
     2771from sage.ext.fast_callable cimport Wrapper
     2772{% print s.pyx_header %}
     2773
     2774cdef extern {{ self.func_header(cython=true) -}}
     2775{% if s.err_return != 'NULL' %}
     2776 except? {{ s.err_return -}}
     2777{% endif %}
     2778
     2779cdef class Wrapper_{{ s.name }}(Wrapper):
     2780{% for ty in types %}
     2781{% print indent_lines(4, ty.class_member_declarations) %}
     2782{% endfor %}
     2783{% for ch in s.chunks %}
     2784{% print ch.declare_class_members() %}
     2785{% endfor %}
     2786
     2787    def __init__(self, args):
     2788        Wrapper.__init__(self, args, metadata)
     2789        cdef int i
     2790        cdef int count
     2791{% for ty in types %}
     2792{% print indent_lines(8, ty.local_declarations) %}
     2793{% print indent_lines(8, ty.class_member_initializations) %}
     2794{% endfor %}
     2795{% for ch in s.chunks %}
     2796{% print ch.init_class_members() %}
     2797{% endfor %}
     2798
     2799    def __dealloc__(self):
     2800        cdef int i
     2801{% for ch in s.chunks %}
     2802{% print ch.dealloc_class_members() %}
     2803{% endfor %}
     2804
     2805    def __call__(self, *args):
     2806        if self._n_args != len(args): raise ValueError
     2807{% for ty in types %}
     2808{% print indent_lines(8, ty.local_declarations) %}
     2809{% endfor %}
     2810{% print indent_lines(8, arg_ch.setup_args()) %}
     2811{% for ch in s.chunks %}
     2812{% print ch.declare_call_locals() %}
     2813{% endfor %}
     2814{% if do_cleanup %}
     2815        try:
     2816{% print indent_lines(4, the_call) %}
     2817        except:
     2818{%   for ch in s.chunks %}
     2819{%     if ch.needs_cleanup_on_error() %}
     2820{%       print indent_lines(12, ch.handle_cleanup()) %}
     2821{%     endif %}
     2822{%   endfor %}
     2823            raise
     2824{% else %}
     2825{% print the_call %}
     2826{% endif %}
     2827{% if not s.return_type %}
     2828        return retval
     2829{% endif %}
     2830
     2831from sage.ext.fast_callable import CompilerInstrSpec, InterpreterMetadata
     2832metadata = InterpreterMetadata(by_opname={
     2833{% for instr in s.instr_descs %}
     2834  '{{ instr.name }}':
     2835  (CompilerInstrSpec({{ instr.n_inputs }}, {{ instr.n_outputs }}, {{ instr.parameters }}), {{ instr.opcode }}),
     2836{% endfor %}
     2837 },
     2838 by_opcode=[
     2839{% for instr in s.instr_descs %}
     2840  ('{{ instr.name }}',
     2841   CompilerInstrSpec({{ instr.n_inputs }}, {{ instr.n_outputs }}, {{ instr.parameters }})),
     2842{% endfor %}
     2843 ])
     2844""", s=s, self=self, types=types, arg_ch=arg_ch, indent_lines=indent_lines, the_call=the_call, do_cleanup=do_cleanup))
     2845
     2846    def get_interpreter(self):
     2847        r"""
     2848        Returns the code for the C interpreter.
     2849
     2850        EXAMPLES:
     2851
     2852        First we get the InterpreterSpec for several interpreters:
     2853            sage: from sage.ext.gen_interpreters import *
     2854            sage: rdf_spec = RDFInterpreter()
     2855            sage: rr_spec = RRInterpreter()
     2856            sage: el_spec = ElementInterpreter()
     2857
     2858        Then we get the actual interpreter code:
     2859            sage: rdf_interp = InterpreterGenerator(rdf_spec).get_interpreter()
     2860            sage: rr_interp = InterpreterGenerator(rr_spec).get_interpreter()
     2861            sage: el_interp = InterpreterGenerator(el_spec).get_interpreter()
     2862
     2863        Now we can look through these interpreters.
     2864
     2865        Each interpreter starts with a file header; this can be
     2866        customized on a per-interpreter basis:
     2867            sage: print rr_interp
     2868            /* Automatically generated.  Do not edit! */
     2869            #include <Python.h>
     2870            #include <mpfr.h>
     2871            #include "wrapper_rr.h"
     2872            ...
     2873
     2874        Next is the function header, with one argument per memory chunk
     2875        in the interpreter spec.
     2876            sage: print el_interp
     2877            /* ... */ ...
     2878            PyObject* interp_el(PyObject** args,
     2879                    PyObject** constants,
     2880                    PyObject** stack,
     2881                    PyObject* domain,
     2882                    int* code) {
     2883            ...
     2884
     2885        Currently, the interpreters have a very simple structure; just
     2886        grab the next instruction and execute it, in a switch
     2887        statement.
     2888            sage: print rdf_interp
     2889            /* ... */ ...
     2890              while (1) {
     2891                switch (*code++) {
     2892            ...
     2893
     2894        Then comes the code for each instruction.  Here is one of the
     2895        simplest instructions:
     2896            sage: print rdf_interp
     2897            /* ... */ ...
     2898                case 9: /* neg */
     2899                  {
     2900                    double i0 = *--stack;
     2901                    double o0;
     2902                    o0 = -i0;
     2903                    *stack++ = o0;
     2904                  }
     2905                  break;
     2906            ...
     2907
     2908        We simply pull the top of the stack into a variable, negate it,
     2909        and write the result back onto the stack.
     2910
     2911        Let's look at the MPFR-based version of this instruction.
     2912        This is an example of an interpreter with an auto-reference
     2913        type.
     2914            sage: print rr_interp
     2915            /* ... */ ...
     2916                case 9: /* neg */
     2917                  {
     2918                    mpfr_ptr i0 = *--stack;
     2919                    mpfr_ptr o0 = *stack++;
     2920                    mpfr_neg(o0, i0, GMP_RNDN);
     2921                  }
     2922                  break;
     2923            ...
     2924
     2925        Here we see that the input and output variables are actually
     2926        just pointers into the stack.  But due to the auto-reference
     2927        trick, the actual code snippet, ``mpfr_net(o0, i0, GMP_RNDN);``,
     2928        is exactly the same as if i0 and o0 were declared as local
     2929        mpfr_t variables.
     2930
     2931        For completeness, let's look at this instruction in the
     2932        Python-object element interpreter.
     2933            sage: print el_interp
     2934            /* ... */ ...
     2935                case 9: /* neg */
     2936                  {
     2937                    PyObject* i0 = *--stack;
     2938                    *stack = NULL;
     2939                    PyObject* o0;
     2940                    o0 = PyNumber_Negative(i0);
     2941                    Py_DECREF(i0);
     2942                    if (!CHECK(o0)) {
     2943                      Py_XDECREF(o0);
     2944                      goto error;
     2945                    }
     2946                    *stack++ = o0;
     2947                  }
     2948                  break;
     2949            ...
     2950
     2951        The original code snippet was only ``o0 = PyNumber_Negative(i0);``;
     2952        all the rest is automatically generated.  For ElementInterpreter,
     2953        the CHECK macro actually checks for an exception (makes sure that
     2954        o0 is not NULL), tests if the o0 is an element with the correct
     2955        parent, and if not converts it into the correct parent.  (That is,
     2956        it can potentially modify the variable o0.)
     2957        """
     2958        import cStringIO
     2959        buff = cStringIO.StringIO()
     2960        self.write_interpreter(buff.write)
     2961        return buff.getvalue()
     2962
     2963    def get_wrapper(self):
     2964        r"""
     2965        Returns the code for the Cython wrapper.
     2966
     2967        EXAMPLES:
     2968
     2969        First we get the InterpreterSpec for several interpreters:
     2970            sage: from sage.ext.gen_interpreters import *
     2971            sage: rdf_spec = RDFInterpreter()
     2972            sage: rr_spec = RRInterpreter()
     2973            sage: el_spec = ElementInterpreter()
     2974
     2975        Then we get the actual wrapper code:
     2976            sage: rdf_wrapper = InterpreterGenerator(rdf_spec).get_wrapper()
     2977            sage: rr_wrapper = InterpreterGenerator(rr_spec).get_wrapper()
     2978            sage: el_wrapper = InterpreterGenerator(el_spec).get_wrapper()
     2979
     2980        Now we can look through these wrappers.
     2981
     2982        Each wrapper starts with a file header; this can be
     2983        customized on a per-interpreter basis (some blank lines have been
     2984        elided below):
     2985            sage: print rdf_wrapper
     2986            # Automatically generated.  Do not edit!
     2987            include "../stdsage.pxi"
     2988            from python_object cimport PyObject
     2989            cdef extern from "Python.h":
     2990                void Py_DECREF(PyObject *o)
     2991                void Py_INCREF(PyObject *o)
     2992                void Py_CLEAR(PyObject *o)
     2993            cdef extern from "listobject.h":
     2994                object PyList_New(Py_ssize_t len)
     2995                ctypedef struct PyListObject:
     2996                PyObject **ob_item
     2997            cdef extern from "tupleobject.h":
     2998                ctypedef struct PyTupleObject:
     2999                    PyObject **ob_item
     3000            from sage.ext.fast_callable cimport Wrapper
     3001            ...
     3002
     3003        Next is the declaration of the C interpreter function.
     3004            sage: print rdf_wrapper
     3005            # ...
     3006            cdef extern double interp_rdf(double* args,
     3007                    double* constants,
     3008                    PyObject** py_constants,
     3009                    double* stack,
     3010                    int* code) except? -1094648009105371
     3011            ...
     3012
     3013        We need a way to propagate exceptions back to the wrapper,
     3014        even though we only return a double from interp_rdf.  The
     3015        ``except? -1094648009105371`` (that's a randomly chosen
     3016        number) means that we will return that number if there's an
     3017        exception, but the wrapper still has to check whether that's a
     3018        legitimate return or an exception.  (Cython does this
     3019        automatically.)
     3020
     3021        Next comes the actual wrapper class, which starts off with
     3022        a list of member declarations.
     3023            sage: print rdf_wrapper
     3024            # ...
     3025            cdef class Wrapper_rdf(Wrapper):
     3026                cdef int _n_args
     3027                cdef double* _args
     3028                cdef int _n_constants
     3029                cdef double* _constants
     3030                cdef object _list_py_constants
     3031                cdef int _n_py_constants
     3032                cdef PyObject** _py_constants
     3033                cdef int _n_stack
     3034                cdef double* _stack
     3035                cdef int _n_code
     3036                cdef int* _code
     3037            ...
     3038
     3039        Contrast the declaration of ``_stack`` here with the
     3040        ElementInterpreter version.  To simplify our handling of
     3041        reference counting and garbage collection, in a Python-object
     3042        based interpreter, we allocate arrays as Python lists,
     3043        and then pull the array out of the innards of the list.
     3044            sage: print el_wrapper
     3045            # ...
     3046                cdef object _list_stack
     3047                cdef int _n_stack
     3048                cdef PyObject** _stack
     3049            ...
     3050
     3051        Next is the __init__ method, which starts like this:
     3052            sage: print rdf_wrapper
     3053            # ...
     3054                def __init__(self, args):
     3055                    Wrapper.__init__(self, args, metadata)
     3056                    cdef int i
     3057                    cdef int count
     3058            ...
     3059
     3060        To make it possible to generate code for all expression
     3061        interpreters with a single code generator, all wrappers
     3062        have the same API.  The __init__ method takes a single
     3063        argument (here called *args*), which is a dictionary holding
     3064        all the information needed to initialize this wrapper.
     3065
     3066        We call Wrapper.__init__, which saves a copy of this arguments
     3067        object and of the interpreter metadata in the wrapper.  (This is
     3068        only used for debugging.)
     3069
     3070        Now we allocate memory for each memory chunk.  (We allocate
     3071        the memory here, and reuse it on each call of the
     3072        wrapper/interpreter.  This is for speed reasons; in a fast
     3073        interpreter like RDFInterpreter, there are no memory allocations
     3074        involved in a call of the wrapper, except for the ones that
     3075        are required by the Python calling convention.  Eventually
     3076        we will support alternate Cython-only entry points that do
     3077        absolutely no memory allocation.)
     3078
     3079        Basically the same code is repeated, with minor variations, for
     3080        each memory chunk; for brevity, we'll only show the code
     3081        for 'constants'.
     3082
     3083            sage: print rdf_wrapper
     3084            # ...
     3085                    val = args['constants']
     3086                    self._n_constants = len(val)
     3087                    self._constants = <double*>sage_malloc(sizeof(double) * len(val))
     3088                    if self._constants == NULL: raise MemoryError
     3089                    for i in range(len(val)):
     3090                        self._constants[i] = val[i]
     3091            ...
     3092
     3093        Recall that _n_constants is an int, and _constants is a
     3094        double*.
     3095
     3096        The RRInterpreter version is more complicated, because it has to
     3097        call mpfr_init.
     3098            sage: print rr_wrapper
     3099            # ...
     3100                    cdef RealNumber rn
     3101            ...
     3102                    val = args['constants']
     3103                    self._n_constants = len(val)
     3104                    self._constants = <mpfr_t*>sage_malloc(sizeof(mpfr_t) * len(val))
     3105                    if self._constants == NULL: raise MemoryError
     3106                    for i in range(len(val)):
     3107                        mpfr_init2(self._constants[i], self.domain.prec())
     3108                    for i in range(len(val)):
     3109                        rn = self.domain(val[i])
     3110                        mpfr_set(self._constants[i], rn.value, GMP_RNDN)
     3111            ...
     3112
     3113        And as we mentioned above, in Python-object based interpreters
     3114        we actually allocate the memory as a Python list.
     3115            sage: print el_wrapper
     3116            # ...
     3117                    val = args['constants']
     3118                    self._n_constants = len(val)
     3119                    self._list_constants = PyList_New(self._n_constants)
     3120                    self._constants = (<PyListObject *>self._list_constants).ob_item
     3121                    for i in range(len(val)):
     3122                        self._constants[i] = <PyObject *>val[i]; Py_INCREF(self._constants[i])
     3123            ...
     3124
     3125        Of course, once we've allocated the memory, we eventually have
     3126        to free it.  (Again, we'll only look at 'constants'.)
     3127            sage: print rdf_wrapper
     3128            # ...
     3129                def __dealloc__(self):
     3130            ...
     3131                    if self._constants:
     3132                        sage_free(self._constants)
     3133            ...
     3134           
     3135        The RRInterpreter code is more complicated again because it has
     3136        to call mpfr_clear.
     3137            sage: print rr_wrapper
     3138            # ...
     3139                def __dealloc__(self):
     3140                    cdef int i
     3141            ...
     3142                    if self._constants:
     3143                        for i in range(self._n_constants):
     3144                            mpfr_clear(self._constants[i])
     3145                        sage_free(self._constants)
     3146            ...
     3147
     3148        But the ElementInterpreter code is extremely simple --
     3149        it doesn't have to do anything to deallocate constants!
     3150        (Since the memory for constants is actually allocated as a
     3151        Python list, and Cython knows how to deallocate Python lists.)
     3152
     3153        Finally we get to the __call__ method.  We grab the arguments
     3154        passed by the caller, stuff them in our pre-allocated
     3155        argument array, and then call the C interpreter.
     3156            sage: print rdf_wrapper
     3157            # ...
     3158                def __call__(self, *args):
     3159                    if self._n_args != len(args): raise ValueError
     3160                    cdef double* c_args = self._args
     3161                    cdef int i
     3162                    for i from 0 <= i < len(args):
     3163                        self._args[i] = args[i]
     3164                    return interp_rdf(c_args
     3165                        , self._constants
     3166                        , self._py_constants
     3167                        , self._stack
     3168                        , self._code
     3169                        )
     3170            ...
     3171
     3172        In Python-object based interpreters, the call to the C
     3173        interpreter has to be a little more complicated.  We don't
     3174        want to hold on to Python objects from an old computation by
     3175        leaving them referenced from the stack.  In normal operation,
     3176        the C interpreter clears out the stack as it runs, leaving the
     3177        stack totally clear when the interpreter finishes.  However,
     3178        this doesn't happen if the C interpreter raises an exception.
     3179        In that case, we have to clear out any remnants from the stack
     3180        in the wrapper.
     3181            sage: print el_wrapper
     3182            # ...
     3183                    try:
     3184                        return interp_el((<PyListObject*>mapped_args).ob_item
     3185                            , self._constants
     3186                            , self._stack
     3187                            , <PyObject*>self._domain
     3188                            , self._code
     3189                            )
     3190                    except:
     3191                        for i in range(self._n_stack):
     3192                            Py_CLEAR(self._stack[i])
     3193                        raise
     3194            ...
     3195
     3196        That's it for the wrapper class.  The only thing remaining is
     3197        the interpreter metadata.  This is the information necessary
     3198        for the code generator to map instruction names to opcodes; it
     3199        also gives information about stack usage, etc.  This is fully
     3200        documented at InterpreterMetadata; for now, we'll just show
     3201        what it looks like.
     3202
     3203        Currently, there are two parts to the metadata; the first maps
     3204        instruction names to instruction descriptions.  The second one
     3205        maps opcodes to instruction descriptions.  Note that we don't
     3206        use InstrSpec objects here; instead, we use CompilerInstrSpec
     3207        objects, which are much simpler and contain only the information
     3208        we'll need at runtime.
     3209
     3210        First the part that maps instruction names to
     3211        (CompilerInstrSpec, opcode) pairs.
     3212
     3213            sage: print rdf_wrapper
     3214            # ...
     3215            from sage.ext.fast_callable import CompilerInstrSpec, InterpreterMetadata
     3216            metadata = InterpreterMetadata(by_opname={
     3217            ...
     3218              'return':
     3219              (CompilerInstrSpec(1, 0, []), 2),
     3220              'py_call':
     3221              (CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs']), 3),
     3222              'add':
     3223              (CompilerInstrSpec(2, 1, []), 4),
     3224            ...
     3225             }, ...)
     3226
     3227        There's also a table that maps opcodes to (instruction name,
     3228        CompilerInstrSpec) pairs:
     3229            sage: print rdf_wrapper
     3230            # ...
     3231            metadata = InterpreterMetadata(...,  by_opcode=[
     3232            ...
     3233              ('return',
     3234               CompilerInstrSpec(1, 0, [])),
     3235              ('py_call',
     3236               CompilerInstrSpec(0, 1, ['py_constants', 'n_inputs'])),
     3237              ('add',
     3238               CompilerInstrSpec(2, 1, [])),
     3239            ...
     3240             ])
     3241
     3242        And that's it for the wrapper.           
     3243        """
     3244        import cStringIO
     3245        buff = cStringIO.StringIO()
     3246        self.write_wrapper(buff.write)
     3247        return buff.getvalue()
     3248
     3249def write_if_changed(fn, value):
     3250    r"""
     3251    Writes value to the file named fn, if value is different than
     3252    the current contents.
     3253
     3254    EXAMPLES:
     3255        sage: from sage.ext.gen_interpreters import *
     3256        sage: def last_modification(fn): return os.stat(fn).st_mtime
     3257        sage: fn = tmp_filename('gen_interp')
     3258        sage: write_if_changed(fn, 'Hello, world')
     3259        sage: t1 = last_modification(fn)
     3260        sage: open(fn).read()
     3261        'Hello, world'
     3262        sage: sleep(2)            # long time
     3263        sage: write_if_changed(fn, 'Goodbye, world')
     3264        sage: t2 = last_modification(fn)
     3265        sage: open(fn).read()
     3266        'Goodbye, world'
     3267        sage: sleep(2)            # long time
     3268        sage: write_if_changed(fn, 'Goodbye, world')
     3269        sage: t3 = last_modification(fn)
     3270        sage: open(fn).read()
     3271        'Goodbye, world'
     3272        sage: t1 == t2            # long time
     3273        False
     3274        sage: t2 == t3
     3275        True
     3276    """
     3277    old_value = None
     3278    try:
     3279        with open(fn) as file:
     3280            old_value = file.read()
     3281    except IOError:
     3282        pass
     3283
     3284    if value != old_value:
     3285        # We try to remove the file, in case it exists.  This is to
     3286        # automatically break hardlinks... see #5350 for motivation.
     3287        try:
     3288            os.remove(fn)
     3289        except OSError:
     3290            pass
     3291
     3292        with open(fn, 'w') as file:
     3293            file.write(value)
     3294
     3295def build_interp(interp_spec, dir):
     3296    r"""
     3297    Given an InterpreterSpec, writes the C interpreter and the Cython
     3298    wrapper.
     3299
     3300    EXAMPLES:
     3301        sage: from sage.ext.gen_interpreters import *
     3302        sage: testdir = tmp_filename()
     3303        sage: os.mkdir(testdir)
     3304        sage: rdf_interp = RDFInterpreter()
     3305        sage: build_interp(rdf_interp, testdir)
     3306        sage: open(testdir + '/interp_rdf.c').readline()
     3307        '/* Automatically generated.  Do not edit! */\n'
     3308    """
     3309    ig = InterpreterGenerator(interp_spec)
     3310    interp_fn = '%s/interp_%s.c' % (dir, interp_spec.name)
     3311    wrapper_fn = '%s/wrapper_%s.pyx' % (dir, interp_spec.name)
     3312    interp = ig.get_interpreter()
     3313    wrapper = ig.get_wrapper()
     3314    write_if_changed(interp_fn, interp)
     3315    write_if_changed(wrapper_fn, wrapper)
     3316
     3317def rebuild(dir):
     3318    r"""
     3319    Check whether the interpreter and wrapper sources have been written
     3320    since the last time this module was changed.  If not, write them.
     3321
     3322    EXAMPLES:
     3323        sage: from sage.ext.gen_interpreters import *
     3324        sage: testdir = tmp_filename()
     3325        sage: os.mkdir(testdir)
     3326        sage: rebuild(testdir)
     3327        Building interpreters for fast_callable
     3328        sage: rebuild(testdir)
     3329        sage: open(testdir + '/wrapper_el.pyx').readline()
     3330        '# Automatically generated.  Do not edit!\n'
     3331    """
     3332    module_mtime = os.stat(__file__).st_mtime
     3333    try:
     3334        if os.stat(dir + '/timestamp').st_mtime > module_mtime:
     3335            # No need to rebuild.
     3336            return
     3337    except OSError:
     3338        pass
     3339
     3340    # This line will show up "sage -b" (once per upgrade, not every time
     3341    # you run it).
     3342    print "Building interpreters for fast_callable"
     3343
     3344    interp = RDFInterpreter()
     3345    build_interp(interp, dir)
     3346
     3347    interp = RRInterpreter()
     3348    build_interp(interp, dir)
     3349
     3350    interp = PythonInterpreter()
     3351    build_interp(interp, dir)
     3352
     3353    interp = ElementInterpreter()
     3354    build_interp(interp, dir)
     3355
     3356    # Do this last, so we don't do it if there's an error above.
     3357    with open(dir + '/timestamp', 'w'):
     3358        pass
     3359
     3360# This list of modules gets added to the list in module_list.py.
     3361# For now, that's not important -- we could have just put this
     3362# list in module_list.py directly.  But eventually, we'll have
     3363# interpreters that are conditionally built (for example,
     3364# interpreters that rely on SSE so they only work on x86), so
     3365# it makes sense to keep the decisions about which interpreters
     3366# to write and which interpreters to build in the same place.
     3367modules = [
     3368    Extension('sage.ext.interpreters.wrapper_rdf',
     3369              sources = ['sage/ext/interpreters/wrapper_rdf.pyx',
     3370                         'sage/ext/interpreters/interp_rdf.c']),
     3371
     3372    Extension('sage.ext.interpreters.wrapper_rr',
     3373              sources = ['sage/ext/interpreters/wrapper_rr.pyx',
     3374                         'sage/ext/interpreters/interp_rr.c'],
     3375              libraries=['mpfr']),
     3376
     3377    Extension('sage.ext.interpreters.wrapper_py',
     3378              sources = ['sage/ext/interpreters/wrapper_py.pyx',
     3379                         'sage/ext/interpreters/interp_py.c']),
     3380
     3381    Extension('sage.ext.interpreters.wrapper_el',
     3382              sources = ['sage/ext/interpreters/wrapper_el.pyx',
     3383                         'sage/ext/interpreters/interp_el.c']),
     3384
     3385]
  • new file sage/ext/interpreters/__init__.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/ext/interpreters/__init__.py
    - +  
     1
  • sage/functions/constants.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/functions/constants.py
    a b  
    434434        from sage.ext.fast_eval import fast_float_constant
    435435        return fast_float_constant(self)
    436436
     437    def _fast_callable_(self, etb):
     438        r"""
     439        Given an ExpressionTreeBuilder, return an Expression representing
     440        this value.
     441
     442        EXAMPLES::
     443
     444            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     445            sage: etb = ExpressionTreeBuilder(vars=['x'])
     446            sage: pi._fast_callable_(etb)
     447            pi
     448            sage: e._fast_callable_(etb)
     449            e
     450            sage: golden_ratio._fast_callable_(etb)
     451            golden_ratio
     452            sage: log2._fast_callable_(etb)
     453            log2
     454            sage: etb = ExpressionTreeBuilder(vars=['x'], domain=RealField(150))
     455            sage: pi._fast_callable_(etb)
     456            3.1415926535897932384626433832795028841971694
     457            sage: e._fast_callable_(etb)
     458            2.7182818284590452353602874713526624977572471
     459            sage: golden_ratio._fast_callable_(etb)
     460            1.6180339887498948482045868343656381177203092
     461            sage: log2._fast_callable_(etb)
     462            0.69314718055994530941723212145817656807550013
     463        """
     464        return etb.constant(self)
     465
    437466    def floor(self):
    438467        """
    439468        Returns the floor of self.
  • sage/misc/log.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/misc/log.py
    a b  
    6565import latex
    6666import misc
    6767
    68 import sage.plot.all
    69 
    7068from   sage.misc.viewer  import browser, dvi_viewer
    7169
    7270offset = 0
  • sage/misc/mathml.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/misc/mathml.py
    a b  
    2727
    2828import os.path
    2929
    30 import sage.plot.all
    31 
    3230from misc import tmp_dir
    3331
    3432def list_function(x):
  • sage/numerical/optimize.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/numerical/optimize.py
    a b  
    229229       
    230230    """
    231231    from sage.calculus.calculus import SymbolicExpression
     232    from sage.ext.fast_eval import fast_callable
    232233    import scipy
    233234    from scipy import optimize
    234235    if isinstance(func,SymbolicExpression):
    235236        var_list=func.variables()
    236237        var_names=map(str,var_list)
    237         fast_f=func._fast_float_(*var_names)
     238        fast_f=fast_callable(func, vars=var_names, domain=RDF)
    238239        f=lambda p: fast_f(*p)
    239240        gradient_list=func.gradient()
    240         fast_gradient_functions=[gradient_list[i]._fast_float_(*var_names)  for i in xrange(len(gradient_list))]
     241        fast_gradient_functions=[fast_callable(gradient_list[i], vars=var_names, domain=RDF)  for i in xrange(len(gradient_list))]
    241242        gradient=lambda p: scipy.array([ a(*p) for a in fast_gradient_functions])       
    242243    else:
    243244        f=func
     
    259260        elif algorithm=="ncg":
    260261            if isinstance(func,SymbolicExpression):
    261262                hess=func.hessian()
    262                 hess_fast= [ [a._fast_float_(*var_names) for a in row] for row in hess]
     263                hess_fast= [ [fast_callable(a, vars=var_names, domain=RDF) for a in row] for row in hess]
    263264                hessian=lambda p: [[a(*p) for a in row] for row in hess_fast]
    264265                hessian_p=lambda p,v: scipy.dot(scipy.array(hessian(p)),v)
    265266                min= optimize.fmin_ncg(f,map(float,x0),fprime=gradient,fhess=hessian,fhess_p=hessian_p,**args)
  • sage/plot/plot.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/plot/plot.py
    a b  
    27922792   
    27932793        sage: x,y = var('x,y')
    27942794        sage: sage.plot.plot.setup_for_eval_on_grid([x^2 + y^2], (x,0,5), (y,0,pi), 11)
    2795         ([<sage.ext.fast_eval.FastDoubleFunc object at ...>],
     2795        ([<sage.ext... object at ...>],
    27962796         0.5,
    27972797         0.31415926535897931,
    27982798         (0.0, 5.0),
     
    28032803    ::
    28042804   
    28052805        sage: sage.plot.plot.setup_for_eval_on_grid([x^2+y^2], (x,0,1), (y,-1,1), 1)
    2806         ([<sage.ext.fast_eval.FastDoubleFunc object at ...>],
     2806        ([<sage.ext... object at ...>],
    28072807        1.0,
    28082808        2.0,
    28092809        (0.0, 1.0),
  • sage/rings/polynomial/multi_polynomial.pyx

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/rings/polynomial/multi_polynomial.pyx
    a b  
    225225
    226226        TESTS:
    227227            sage: from sage.ext.fast_eval import fast_float
    228             sage: list(fast_float(K(0)))
     228            sage: list(fast_float(K(0), old=True))
    229229            ['push 0.0']
    230             sage: list(fast_float(K(17)))
     230            sage: list(fast_float(K(17), old=True))
    231231            ['push 0.0', 'push 17.0', 'add']
    232             sage: list(fast_float(y))
     232            sage: list(fast_float(y, old=True))
    233233            ['push 0.0', 'push 1.0', 'load 1', 'mul', 'add']
    234234        """
    235235        from sage.ext.fast_eval import fast_float_arg, fast_float_constant
     
    248248            expr = expr + monom
    249249        return expr
    250250
     251    def _fast_callable_(self, etb):
     252        """
     253        Given an ExpressionTreeBuilder, return an Expression representing
     254        this value.
     255
     256        EXAMPLES:
     257            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     258            sage: etb = ExpressionTreeBuilder(vars=['x','y','z'])
     259            sage: K.<x,y,z> = QQ[]
     260            sage: v = K.random_element(degree=3, terms=4); v
     261            -6/5*x*y*z + 2*y*z^2 - x
     262            sage: v._fast_callable_(etb)
     263            add(add(add(0, mul(-6/5, mul(mul(pow(v_0, 1), pow(v_1, 1)), pow(v_2, 1)))), mul(2, mul(pow(v_1, 1), pow(v_2, 2)))), mul(-1, pow(v_0, 1)))
     264       
     265        TESTS:
     266            sage: K.<x,y,z> = QQ[]
     267            sage: from sage.ext.fast_eval import fast_float
     268            sage: fast_float(K(0)).op_list()
     269            [('load_const', 0.0), 'return']
     270            sage: fast_float(K(17)).op_list()
     271            [('load_const', 0.0), ('load_const', 17.0), 'add', 'return']
     272            sage: fast_float(y).op_list()
     273            [('load_const', 0.0), ('load_const', 1.0), ('load_arg', 1), ('load_const', 1.0), 'pow', 'mul', 'add', 'return']
     274        """
     275        my_vars = self.parent().variable_names()
     276        x = [etb.var(v) for v in my_vars]
     277        n = len(x)
     278
     279        expr = etb.constant(0)
     280        for (m, c) in self.dict().iteritems():
     281            monom = misc.mul([ x[i]**m[i] for i in range(n) if m[i] != 0],
     282                             etb.constant(c))
     283            expr = expr + monom
     284        return expr
    251285
    252286    def derivative(self, *args):
    253287        r"""
  • sage/rings/polynomial/polynomial_element.pyx

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c sage/rings/polynomial/polynomial_element.pyx
    a b  
    618618                expr *= x
    619619        return expr
    620620       
     621    def _fast_callable_(self, etb):
     622        r"""
     623        Given an ExpressionTreeBuilder, return an Expression representing
     624        this value.
     625
     626        EXAMPLES::
     627
     628            sage: from sage.ext.fast_callable import ExpressionTreeBuilder
     629            sage: etb = ExpressionTreeBuilder(vars=['t'])
     630            sage: R.<t> = QQ[]
     631            sage: v = R.random_element(6); v
     632            -t^6 - 12*t^5 + 1/2*t^4 - 1/95*t^3 - 1/2*t^2 - 4
     633            sage: v._fast_callable_(etb)
     634            add(mul(mul(add(mul(add(mul(add(mul(add(mul(v_0, -1), -12), v_0), 1/2), v_0), -1/95), v_0), -1/2), v_0), v_0), -4)
     635        """
     636        x = etb.var(self.variable_name())
     637        expr = x
     638        cdef int i, d = self.degree()
     639        coeff = self[d]
     640        # We handle polynomial rings like QQ['x']['y']; that gives us some
     641        # slowdown.  Optimize away some of that:
     642        if len(etb._vars) == 1:
     643            # OK, we're in the (very common) univariate case.
     644            coeff_maker = etb.constant
     645        else:
     646            # There may be variables in our coefficients...
     647            coeff_maker = etb.make
     648        if coeff != 1:
     649            expr *= coeff_maker(coeff)
     650        for i from d > i >= 0:
     651            coeff = self[i]
     652            if coeff:
     653                expr += coeff_maker(coeff)
     654            if i > 0:
     655                expr *= x
     656        return expr
     657
    621658    cdef int _cmp_c_impl(self, Element other) except -2:
    622659        """
    623660        Compare the two polynomials self and other.
  • setup.py

    diff -r 5ba8b66a0ae4 -r 0949e5d7cc4c setup.py
    a b  
    88### List of Extensions
    99###
    1010### Since Sage 3.2 the list of extensions resides in
    11 ### module_list.py in the same directory as this file.
     11### module_list.py in the same directory as this file
     12### (augmented by the list of interpreters
     13### generated by sage/ext/gen_interpreters.py)
    1214#########################################################
    1315
    1416from module_list import ext_modules
     17import sage.ext.gen_interpreters
    1518
    1619#########################################################
    1720### Configuration
     
    7275if DEVEL:
    7376    extra_compile_args.append('-ggdb')
    7477
     78# Generate interpreters
     79
     80sage.ext.gen_interpreters.rebuild(SAGE_DEVEL + 'sage/sage/ext/interpreters')
     81ext_modules = ext_modules + sage.ext.gen_interpreters.modules
    7582
    7683
    7784######################################################################
     
    562569                     'sage.databases',
    563570                     
    564571                     'sage.ext',
     572                     'sage.ext.interpreters',
    565573                     
    566574                     'sage.finance',
    567575