Ticket #7377: trac_7377-fastcalculus-rebased.patch

File trac_7377-fastcalculus-rebased.patch, 14.6 KB (added by jpflori, 9 years ago)

Rebased on Sage 4.6.1, fixed changes in sage.symbolic.integration.external

  • sage/calculus/calculus.py

    # HG changeset patch
    # User Nils Bruin <nbruin@sfu.ca>
    # Date 1265693495 28800
    # Node ID b9db51bba0134fddc27b3af64ca95f987abb4423
    # Parent  d1def36a6fbc238838fca35740574e46fc6f8064
    Fast calculus interface
    
    diff -r d1def36a6fbc -r b9db51bba013 sage/calculus/calculus.py
    a b  
    386386    f1
    387387"""
    388388maxima = sage.interfaces.maxima_lib.maxima
     389#maxima = Maxima(init_code = ['display2d:false', 'domain: complex',
     390#                             'keepfloat: true', 'load(to_poly_solver)', 'load(simplify_sum)'],
     391#                script_subdirectory=None)
    389392
    390393########################################################
    391394def symbolic_sum(expression, v, a, b, algorithm='maxima'):
     
    529532        raise ValueError, "summation limits must not depend on the summation variable"
    530533
    531534    if algorithm == 'maxima':
    532         sum  = "'sum(%s, %s, %s, %s)" % tuple([repr(expr._maxima_()) for expr in (expression, v, a, b)])
    533535        try:
    534             result = maxima.simplify_sum(sum)
    535             result = result.ratsimp()
    536             return expression.parent()(result)
     536            return maxima.sr_sum(expression,v,a,b)
    537537        except TypeError, error:
    538538            s = str(error)
    539539            if "divergent" in s or 'Pole encountered' in s:
     
    11011101
    11021102    if algorithm == 'maxima':
    11031103        if dir is None:
    1104             l = ex._maxima_().limit(v, a)
     1104            l = maxima.sr_limit(ex, v, a)
    11051105        elif dir in ['plus', '+', 'right', 'above']:
    11061106            if dir == 'above':
    11071107                from sage.misc.misc import deprecation
    11081108                deprecation("the keyword 'above' is deprecated. Please use 'right' or '+' instead.", 'Sage version 4.6')
    1109             l = ex._maxima_().limit(v, a, 'plus')
     1109            l = maxima.sr_limit(ex, v, a, 'plus')
    11101110        elif dir in ['minus', '-', 'left', 'below']:
    11111111            if dir == 'below':
    11121112                from sage.misc.misc import deprecation
    11131113                deprecation("the keyword 'below' is deprecated. Please use 'left' or '-' instead.", 'Sage version 4.6')
    1114             l = ex._maxima_().limit(v, a, 'minus')
     1114            l = maxima.sr_limit(ex, v, a, 'minus')
    11151115    elif algorithm == 'maxima_taylor':
    11161116        if dir is None:
    1117             l = ex._maxima_().tlimit(v, a)
     1117            l = maxima.sr_tlimit(ex, v, a)
    11181118        elif dir == 'plus' or dir == 'above' or dir == 'from_right':
    1119             l = ex._maxima_().tlimit(v, a, 'plus')
     1119            l = maxima.sr_tlimit(ex, v, a, 'plus')
    11201120        elif dir == 'minus' or dir == 'below' or dir == 'from_left':
    1121             l = ex._maxima_().tlimit(v, a, 'minus')
     1121            l = maxima.sr_tlimit(ex, v, a, 'minus')
    11221122    elif algorithm == 'sympy':
    11231123        if dir is None:
    11241124            import sympy
  • sage/interfaces/maxima.py

    diff -r d1def36a6fbc -r b9db51bba013 sage/interfaces/maxima.py
    a b  
    10081008        """
    10091009        return maxima_version()
    10101010
     1011##some helper functions to wrap tha calculus use of the maxima interface.
     1012##these routines expect arguments living in the symbolic ring and return something
     1013##that is hopefully coercible into the symbolic ring again.
    10111014
     1015    def sr_integral(self,*args):
     1016        return args[0]._maxima_().integrate(*args[1:])
     1017
     1018    def sr_sum(self,expression,v,a,b):
     1019        sum  = "'sum(%s, %s, %s, %s)" % tuple([repr(expr._maxima_()) for expr in (expression, v, a, b)])
     1020        result = self.simplify_sum(sum)
     1021        result = result.ratsimp()
     1022        return expression.parent()(result)
     1023
     1024    def sr_limit(self,ex,*args):
     1025        return ex._maxima_().limit(*args)
     1026
     1027    def sr_tlimit(self,ex,*args):
     1028        return ex._maxima_().tlimit(*args)
     1029   
    10121030##     def display2d(self, flag=True):
    10131031##         """
    10141032##         Set the flag that determines whether Maxima objects are
     
    10741092def __doctest_cleanup():
    10751093    import sage.interfaces.quit
    10761094    sage.interfaces.quit.expect_quitall()
    1077 
    1078 
  • sage/interfaces/maxima_lib.py

    diff -r d1def36a6fbc -r b9db51bba013 sage/interfaces/maxima_lib.py
    a b  
    470470import sage.server.support
    471471
    472472import maxima_abstract
    473 from maxima_abstract import MaximaFunctionElement, MaximaExpectFunction, MaximaElement, MaximaFunction, maxima_console
     473from maxima_abstract import MaximaFunctionElement, MaximaExpectFunction, MaximaFunction, maxima_console
    474474
    475475# The Maxima "apropos" command, e.g., apropos(det) gives a list
    476476# of all identifiers that begin in a certain way.  This could
     
    485485ecl_eval("(in-package :maxima)")
    486486ecl_eval("(setq $nolabels t))")
    487487ecl_eval("(defvar *MAXIMA-LANG-SUBDIR* NIL)")
     488ecl_eval("(set-locale-subdir)")
    488489ecl_eval("(set-pathnames)")
    489490ecl_eval("(defun add-lineinfo (x) x)")
    490491ecl_eval("""(defun retrieve (msg flag &aux (print? nil))(error (concatenate 'string "Maxima asks:" (meval (list '($string) msg)))))""")
     
    527528
    528529maxprint=EclObject("$STRING")
    529530meval=EclObject("MEVAL")
     531msetq=EclObject("MSETQ")
     532cadadr=EclObject("CADADR")
     533
     534max_integrate=EclObject("$INTEGRATE")
     535max_sum=EclObject("$SUM")
     536max_simplify_sum=EclObject("$SIMPLIFY_SUM")
     537max_ratsimp=EclObject("$RATSIMP")
     538max_limit=EclObject("$LIMIT")
     539max_tlimit=EclObject("$TLIMIT")
     540max_plus=EclObject("$PLUS")
     541max_minus=EclObject("$MINUS")
    530542
    531543def max_to_string(s):
    532544     return meval(EclObject([[maxprint],s])).python()[1:-1]
     
    575587        for l in init_code:
    576588            ecl_eval("#$%s$"%l)
    577589        ecl_eval("(setf *standard-output* original-standard-output)")
     590 
     591    def _coerce_from_special_method(self, x):
     592        if isinstance(x, EclObject):
     593            return MaximaElement(self,self._create(x))
     594        else:
     595            return maxima_abstract.Maxima._coerce_from_special_method(self,x)
     596       
    578597       
    579598    def _function_class(self):
    580599        """
     
    929948        """
    930949        s = self._eval_line('%s;'%var)
    931950        return s
    932    
     951 
     952    def _create(self, value, name=None):
     953        name = self._next_var_name() if name is None else name
     954        if isinstance(value,EclObject):
     955            maxima_eval([[msetq],cadadr("#$%s$#$"%name),value])
     956        else:
     957            self.set(name, value)
     958        return name
     959               
    933960    def version(self):
    934961        """
    935962        Return the version of Maxima that Sage includes.
     
    941968        """
    942969        return maxima_version()
    943970
     971##some helper functions to wrap tha calculus use of the maxima interface.
     972##these routines expect arguments living in the symbolic ring and return something
     973##that is hopefully coercible into the symbolic ring again.
     974
     975    def sr_integral(self,*args):
     976        try:
     977            return max_to_sr(maxima_eval(([max_integrate],[sr_to_max(SR(a)) for a in args])))
     978        except RuntimeError, error:
     979            s = str(error)
     980            if "Divergent" in s or "divergent" in s:
     981                raise ValueError, "Integral is divergent."
     982            else:
     983                raise error
     984
     985    def sr_sum(self,*args):
     986        try:
     987            return max_to_sr(maxima_eval([[max_ratsimp],[[max_simplify_sum],([max_sum],[sr_to_max(SR(a)) for a in args])]]));
     988        except RuntimeError, error:
     989            s = str(error)
     990            if "divergent" in s:
     991                raise ValueError, "Sum is divergent."
     992            else:
     993                raise error
     994
     995    def sr_limit(self,expr,v,a,dir=None):
     996        L=[sr_to_max(SR(a)) for a in [expr,v,a]]
     997        if dir == "plus":
     998            L.append(max_plus)
     999        elif dir == "minus":
     1000            L.append(max_minus)
     1001        return max_to_sr(maxima_eval(([max_limit],L)))
     1002
     1003    def sr_tlimit(self,expr,v,a,dir=None):
     1004        L=[sr_to_max(SR(a)) for a in [expr,v,a]]
     1005        if dir == "plus":
     1006            L.append(max_plus)
     1007        elif dir == "minus":
     1008            L.append(max_minus)
     1009        return max_to_sr(maxima_eval(([max_tlimit],L)))
    9441010
    9451011##     def display2d(self, flag=True):
    9461012##         """
     
    9641030##         """
    9651031##         self._display2d = bool(flag)
    9661032
     1033class MaximaElement(maxima_abstract.MaximaElement):
     1034    """
     1035    """   
     1036    def ecl(self):
     1037        try:
     1038            return self._ecl
     1039        except AttributeError:
     1040            self._ecl=maxima_eval("#$%s$"%self._name)
     1041            return self._ecl
     1042
    9671043def is_MaximaElement(x):
    9681044    """
    9691045    Returns True if x is of type MaximaElement.
     
    10071083    sage.interfaces.quit.expect_quitall()
    10081084
    10091085
     1086import sage.symbolic.expression
     1087from sage.symbolic.ring import SR
     1088
     1089import sage.symbolic.expression
     1090import sage.functions.trig
     1091import sage.functions.log
     1092import sage.functions.other
     1093
     1094car=EclObject("car")
     1095cdr=EclObject("cdr")
     1096caar=EclObject("caar")
     1097cadadr=EclObject("cadadr")
     1098meval=EclObject("meval")
     1099NIL=EclObject("NIL")
     1100ratdisrep=EclObject("ratdisrep")
     1101
     1102sage_op_dict={}
     1103
     1104sage_op_dict = {
     1105    sage.symbolic.expression.operator.abs : "MABS",
     1106    sage.symbolic.expression.operator.add : "MPLUS",
     1107    sage.symbolic.expression.operator.div : "MQUOTIENT",
     1108    sage.symbolic.expression.operator.eq : "MEQUAL",
     1109    sage.symbolic.expression.operator.ge : "MGEQP",
     1110    sage.symbolic.expression.operator.gt : "MGREATERP",
     1111    sage.symbolic.expression.operator.le : "MLEQP",
     1112    sage.symbolic.expression.operator.lt : "MLESSP",
     1113    sage.symbolic.expression.operator.mul : "MTIMES",
     1114    sage.symbolic.expression.operator.ne : "MNOTEQUAL",
     1115    sage.symbolic.expression.operator.neg : "MMINUS",
     1116    sage.symbolic.expression.operator.pow : "MEXPT",
     1117    sage.symbolic.expression.operator.or_ : "MOR",
     1118    sage.symbolic.expression.operator.and_ : "MAND",
     1119    sage.functions.trig.acos : "%ACOS",
     1120    sage.functions.trig.acot : "%ACOT",
     1121    sage.functions.trig.acsc : "%ACSC",
     1122    sage.functions.trig.asec : "%ASEC",
     1123    sage.functions.trig.asin : "%ASIN",
     1124    sage.functions.trig.atan : "%ATAN",
     1125    sage.functions.trig.cos : "%COS",
     1126    sage.functions.trig.cot : "%COT",
     1127    sage.functions.trig.csc : "%CSC",
     1128    sage.functions.trig.sec : "%SEC",
     1129    sage.functions.trig.sin : "%SIN",
     1130    sage.functions.trig.tan : "%TAN",
     1131    sage.functions.log.exp : "%EXP",
     1132    sage.functions.log.ln : "%LOG",
     1133    sage.functions.log.log : "%LOG",
     1134    sage.functions.other.factorial : "MFACTORIAL",
     1135    sage.functions.other.erf : "%ERF",
     1136}
     1137#we compile the dictionary
     1138sage_op_dict = dict([(k,EclObject(sage_op_dict[k])) for k in sage_op_dict])
     1139max_op_dict = dict([(sage_op_dict[k],k) for k in sage_op_dict])
     1140def add_vararg(*args):
     1141    S=0
     1142    for a in args:
     1143        S=S+a
     1144    return S
     1145
     1146def mul_vararg(*args):
     1147    P=1
     1148    for a in args:
     1149        P=P*a
     1150    return P
     1151
     1152def sage_rat(x,y):
     1153    return x/y
     1154
     1155mplus=EclObject("MPLUS")
     1156mtimes=EclObject("MTIMES")
     1157rat=EclObject("RAT")
     1158max_i=EclObject("$%I")
     1159max_op_dict[mplus]=add_vararg
     1160max_op_dict[mtimes]=mul_vararg
     1161max_op_dict[rat]=sage_rat
     1162
     1163def mrat_to_sage(expr):
     1164    r"""
     1165    Convert a maxima MRAT expression to Sage SR
     1166   
     1167    Maxima has an optimised representation for multivariate rational expressions.
     1168    The easiest way to translate those to SR is by first asking maxima to give
     1169    the generic representation of the object. That is what RATDISREP does in
     1170    maxima.
     1171    """
     1172    return max_to_sage(meval(EclObject([[ratdisrep],expr])))
     1173
     1174special_max_to_sage={
     1175    EclObject("MRAT") : mrat_to_sage
     1176}
     1177
     1178sage_sym_dict={}
     1179max_sym_dict={}
     1180
     1181def pyobject_to_max(obj):
     1182    if isinstance(obj,sage.rings.rational.Rational):
     1183        return EclObject(obj) if (obj.denom().is_one()) else EclObject([[rat], obj.numer(),obj.denom()])
     1184    elif isinstance(obj,sage.rings.number_field.number_field_element_quadratic.NumberFieldElement_quadratic) and obj.parent().defining_polynomial().list() == [1,0,1]:
     1185        re, im = obj.list()
     1186        return EclObject([[mplus], pyobject_to_max(re), [[mtimes], pyobject_to_max(im), max_i]])
     1187   
     1188    return EclObject(obj)
     1189
     1190def sr_to_max(expr):
     1191    r"""
     1192    """
     1193    global sage_op_dict, max_op_dict
     1194    global sage_sym_dict, max_sym_dict
     1195    op = expr.operator()
     1196    if op:
     1197        if not (op in sage_op_dict):
     1198            op_max=caar(maxima(expr).ecl())
     1199            sage_op_dict[op]=op_max
     1200            max_op_dict[op_max]=op
     1201        return EclObject(([sage_op_dict[op]],
     1202                     [sr_to_max(o) for o in expr.operands()]))
     1203    elif expr._is_symbol() or expr._is_constant():
     1204        if not expr in sage_sym_dict:
     1205            sym_max=maxima(expr).ecl()
     1206            sage_sym_dict[expr]=sym_max
     1207            max_sym_dict[sym_max]=expr
     1208        return sage_sym_dict[expr]
     1209    else:
     1210        try:
     1211            return pyobject_to_max(expr.pyobject())
     1212        except TypeError:
     1213            return maxima(expr).ecl()
     1214   
     1215def max_to_sr(expr):
     1216    if expr.consp():
     1217        op_max=caar(expr)
     1218        if op_max in special_max_to_sage:
     1219            return special_max_to_sage[op_max](expr)
     1220        if not(op_max in max_op_dict):
     1221            sage_expr=SR(maxima(expr))
     1222            max_op_dict[op_max]=sage_expr.operator()
     1223            sage_op_dict[sage_expr.operator()]=op_max
     1224        op=max_op_dict[op_max]
     1225        max_args=cdr(expr)
     1226        args=[]
     1227        while not(max_args.nullp()):
     1228            args.append(max_to_sr(car(max_args)))
     1229            max_args=cdr(max_args)
     1230        return op(*args)
     1231    elif expr.symbolp():
     1232        if not(expr in max_sym_dict):
     1233            sage_symbol=SR(maxima(expr))
     1234            sage_sym_dict[sage_symbol]=expr
     1235            max_sym_dict[expr]=sage_symbol
     1236        return max_sym_dict[expr]
     1237    else:
     1238        return expr.python()
     1239
  • sage/symbolic/integration/external.py

    diff -r d1def36a6fbc -r b9db51bba013 sage/symbolic/integration/external.py
    a b  
    1212        sage: maxima_integrator(f(x), x)
    1313        integrate(f(x), x)
    1414    """
     15    from sage.calculus.calculus import maxima
    1516    if not isinstance(expression, Expression):
    1617        expression = SR(expression)
    1718    if a is None:
    18         result = expression._maxima_().integrate(v)
     19        result = maxima.sr_integral(expression,v)
    1920    else:
    2021        try:
    21             result = expression._maxima_().integrate(v, a, b)
     22            result = maxima.sr_integral(expression, v, a, b)
    2223        except TypeError, error:
    2324            s = str(error)
    2425            if "divergent" in s or 'Principal Value' in s:
    2526                raise ValueError, "Integral is divergent."
    2627            else:
    2728                raise
    28     return result.sage()
     29    return result
    2930
    3031def sympy_integrator(expression, v, a=None, b=None):
    3132    """