# Ticket #2516: hyper3.patch

File hyper3.patch, 37.1 KB (added by eviatarbach, 9 years ago)
• ## doc/en/reference/functions/index.rst

# HG changeset patch
# User fredrik@scv
# Date 1280114065 -7200
# Branch hyper
# Parent  42c2de7ba6a5461b04018606bdca93da60991117
working on hypergeometric

diff --git a/doc/en/reference/functions/index.rst b/doc/en/reference/functions/index.rst
 a sage/functions/orthogonal_polys sage/functions/other sage/functions/special sage/functions/hypergeometric sage/functions/exp_integral sage/functions/wigner sage/functions/generalized
• ## sage/calculus/calculus.py

diff --git a/sage/calculus/calculus.py b/sage/calculus/calculus.py
 a maxima_qp = re.compile("\?\%[a-z|A-Z|0-9|_]*")  # e.g., ?%jacobi_cd maxima_var = re.compile("\%[a-z|A-Z|0-9|_]*")  # e.g., ?%jacobi_cd maxima_var = re.compile("\%[a-z|A-Z|0-9|_]*")  # e.g., %jacobi_cd sci_not = re.compile("(-?(?:0|[1-9]\d*))(\.\d+)?([eE][-+]\d+)") maxima_polygamma = re.compile("psi$(\d*)$\(")  # matches psi[n]( where n is a number maxima_hyper = re.compile("\%f$\d+,\d+$")  # matches %f[m,n] def symbolic_expression_from_maxima_string(x, equals_sub=False, maxima=maxima): """ Given a string representation of a Maxima expression, parse it and pass else: syms[X[2:]] = function_factory(X[2:]) s = s.replace("?%","") s = s.replace("?%", "") s = polylog_ex.sub('polylog(\\1,',s) s = maxima_hyper.sub('hypergeometric', s) s = polylog_ex.sub('polylog(\\1,', s) s = multiple_replace(symtable, s) s = s.replace("%","") s = s.replace("#","!=") # a lot of this code should be refactored somewhere... s = maxima_polygamma.sub('psi(\g<1>,',s) # this replaces psi[n](foo) with psi(n,foo), ensuring that derivatives of the digamma function are parsed properly below s = maxima_polygamma.sub('psi(\g<1>,', s) # this replaces psi[n](foo) with psi(n,foo), ensuring that derivatives of the digamma function are parsed properly below if equals_sub: s = s.replace('=','==')
• ## sage/functions/all.py

diff --git a/sage/functions/all.py b/sage/functions/all.py
 a sin_integral, cos_integral, Si, Ci, sinh_integral, cosh_integral, Shi, Chi, exponential_integral_1, Ei) from hypergeometric import hypergeometric
• ## new file sage/functions/hypergeometric.py

diff --git a/sage/functions/hypergeometric.py b/sage/functions/hypergeometric.py
new file mode 100644
 - r""" Hypergeometric Functions This module implements manipulation of infinite hypergeometric series represented in standard parametric form (as $\,_pFq$ functions). AUTHORS: - Fredrik Johansson (2010): initial version - Eviatar Bach (2013): major changes EXAMPLES: Examples from :trac:9908:: sage: maxima('integrate(bessel_j(2, x), x)').sage() 1/24*x^3*hypergeometric((3/2,), (5/2, 3), -1/4*x^2) sage: sum(((2 * I)^x / (x^3 + 1) * (1/4)^x), x, 0, oo) hypergeometric((1, 1, -1/2*I*sqrt(3) - 1/2, 1/2*I*sqrt(3) - 1/2), (2, -1/2*I*sqrt(3) + 1/2, 1/2*I*sqrt(3) + 1/2), 1/2*I) sage: sum((-1)^x / ((2 * x + 1) * factorial(2 * x + 1)), x, 0, oo) hypergeometric((1/2,), (3/2, 3/2), -1/4) Simplification:: sage: hypergeometric((hypergeometric((), (), x),), (), ....:                hypergeometric((), (), x)).simplify_hypergeometric() 1/((-e^x + 1)^e^x) sage: hypergeometric([-2], [], x).simplify_hypergeometric() x^2 - 2*x + 1 sage: hypergeometric([], [], x).simplify_hypergeometric() e^x sage: hypergeometric((hypergeometric((), (), x),), (), ....:                hypergeometric((), (), x)).simplify_hypergeometric(algorithm='sage') (-e^x + 1)^(-e^x) Equality testing:: sage: bool(hypergeometric([], [], x).derivative(x) == ....:      hypergeometric([], [], x))  # diff(e^x, x) == e^x True sage: bool(hypergeometric([], [], x) == hypergeometric([], [1], x)) False Computing terms and series:: sage: z = var('z') sage: hypergeometric([], [], z).series(z, 0) Order(1) sage: hypergeometric([], [], z).series(z, 1) 1 + Order(z) sage: hypergeometric([], [], z).series(z, 2) 1 + 1*z + Order(z^2) sage: hypergeometric([], [], z).series(z, 3) 1 + 1*z + 1/2*z^2 + Order(z^3) sage: hypergeometric([-2], [], z).series(z, 3) 1 + (-2)*z + 1*z^2 sage: hypergeometric([-2], [], z).series(z, 6) 1 + (-2)*z + 1*z^2 sage: hypergeometric([-2], [], z).series(z, 6).is_terminating_series() True sage: hypergeometric([-2], [], z).series(z, 2) 1 + (-2)*z + Order(z^2) sage: hypergeometric([-2], [], z).series(z, 2).is_terminating_series() False sage: hypergeometric([1], [], z).series(z, 6) 1 + 1*z + 1*z^2 + 1*z^3 + 1*z^4 + 1*z^5 + Order(z^6) sage: hypergeometric([], [1/2], -z^2/4).series(z, 11) 1 + (-1/2)*z^2 + 1/24*z^4 + (-1/720)*z^6 + 1/40320*z^8 + (-1/3628800)*z^10 + Order(z^11) sage: hypergeometric([1],[5],x).series(x,5) 1 + 1/5*x + 1/30*x^2 + 1/210*x^3 + 1/1680*x^4 + Order(x^5) sage: sum(hypergeometric([1, 2], [3], 1/3).terms(6)).n() 1.29788359788360 sage: hypergeometric([1, 2], [3], 1/3).n() 1.29837194594696 sage: hypergeometric([], [], x).series(x, 20)(x=1).n() == e.n() True Plotting:: sage: plot(hypergeometric([1, 1], [3, 3, 3], x), x, -30, 30) Numeric evaluation:: sage: hypergeometric([1], [], 1/10).n()  # geometric series 1.11111111111111 sage: hypergeometric([], [], 1).n()  # e 2.71828182845905 sage: hypergeometric([], [], 3., hold=True) hypergeometric((), (), 3.00000000000000) sage: hypergeometric([1,2,3],[4,5,6],1/2).n() 1.02573619590134 sage: hypergeometric([1,2,3],[4,5,6],1/2).n(digits=30) 1.02573619590133865036584139535 sage: hypergeometric([5-3*I],[3/2,2+I,sqrt(2)],4+I).n() #TOLERANCE? 5.52605111678805 - 7.86331357527544*I SymPy conversion:: sage: hypergeometric((5,4), (4,4), 3)._sympy_() hyper((5, 4), (4, 4), 3) """ #***************************************************************************** #       Copyright (C) 2010 Fredrik Johansson #       Copyright (C) 2013 Eviatar Bach # #  Distributed under the terms of the GNU General Public License (GPL) #  as published by the Free Software Foundation; either version 2 of #  the License, or (at your option) any later version. #                  http://www.gnu.org/licenses/ #***************************************************************************** from sage.rings.all import (Integer, ZZ, QQ, Infinity, binomial, rising_factorial, factorial) from sage.functions.other import sqrt, gamma, real_part from sage.functions.log import exp, log from sage.functions.trig import cos, sin from sage.functions.hyperbolic import cosh, sinh from sage.functions.other import erf from sage.symbolic.constants import pi from sage.symbolic.all import I from sage.symbolic.function import BuiltinFunction, is_inexact from sage.symbolic.ring import SR from sage.structure.element import get_coercion_model from sage.misc.latex import latex from sage.misc.flatten import flatten from sage.misc.misc_c import prod from mpmath import hyper from sage.libs.mpmath import utils as mpmath_utils from sage.symbolic.expression import Expression from sage.interfaces.maxima import maxima from sage.calculus.functional import derivative def rational_param_as_tuple(x): """ Utility function for converting rational PFQ parameters to tuples (which mpmath handles more efficiently). sage: from sage.functions.hypergeometric import rational_param_as_tuple sage: rational_param_as_tuple(1/2) (1, 2) sage: rational_param_as_tuple(3) 3 sage: rational_param_as_tuple(pi) pi """ try: x = x.pyobject() except AttributeError: pass try: if x.parent() is QQ: p = int(x.numer()) q = int(x.denom()) return p, q except AttributeError: pass return x class Hypergeometric(BuiltinFunction): r""" Represents a (formal) generalized infinite hypergeometric series. """ def __init__(self): """ EXAMPLES:: sage: hypergeometric([], [], 1) hypergeometric((), (), 1) sage: hypergeometric([], [1], 1) hypergeometric((), (1,), 1) sage: hypergeometric([2, 3], [1], 1) hypergeometric((2, 3), (1,), 1) """ BuiltinFunction.__init__(self, 'hypergeometric', nargs=3, conversions={'mathematica':'HypergeometricPFQ', 'maxima':'hypergeometric', 'sympy':'hyper'}) def __call__(self, a, b, z, **kwargs): return BuiltinFunction.__call__(self, SR._force_pyobject(a), SR._force_pyobject(b), z, **kwargs) def _print_latex_(self, a, b, z): r""" TESTS:: sage: latex(hypergeometric([1, 1], [2], -1)) \,_2F_1\left(\begin{matrix} 1,1 \\ 2 \end{matrix} ; -1 \right) """ aa = ",".join(latex(c) for c in a) bb = ",".join(latex(c) for c in b) z = latex(z) return (r"\,_{}F_{}\left(\begin{{matrix}} {} \\ {} \end{{matrix}} ; " r"{} \right)").format(len(a), len(b), aa, bb, z) def _eval_(self, a, b, z, **kwargs): coercion_model = get_coercion_model() co = reduce(lambda x, y: coercion_model.canonical_coercion(x, y)[0], a + b + (z,)) if is_inexact(co) and not isinstance(co, Expression): from sage.structure.coerce import parent return self._evalf_(a, b, z, parent=parent(co)) if z == 0: return Integer(1) return def _evalf_(self, a, b, z, parent): """ TESTS:: sage: hypergeometric([1, 1], [2], -1).n() 0.693147180559945 sage: hypergeometric([], [], RealField(100)(1)) 2.7182818284590452353602874714 """ aa = [rational_param_as_tuple(c) for c in a] bb = [rational_param_as_tuple(c) for c in b] return mpmath_utils.call(hyper, aa, bb, z, parent=parent) def _tderivative_(self, a, b, z, *args, **kwargs): """ EXAMPLES:: sage: hypergeometric([1/3, 2/3], [5], x^2).diff(x) 4/45*x*hypergeometric((4/3, 5/3), (6,), x^2) sage: hypergeometric([1, 2], [x], 2).diff(x) Traceback (most recent call last): ... NotImplementedError: derivative of hypergeometric function with respect to parameters. Try calling .simplify_hypergeometric() first. """ diff_param = kwargs['diff_param'] if diff_param in hypergeometric(a, b, 1).variables(): raise NotImplementedError("derivative of hypergeometric function " "with respect to parameters. Try calling"                                      " .simplify_hypergeometric() first.") t = (reduce(lambda x, y: x * y, a, 1) * reduce(lambda x, y: x / y, b, Integer(1))) return (t * derivative(z, diff_param) * hypergeometric([c + 1 for c in a], [c + 1 for c in b], z)) '''def _maxima_init_evaled_(self, a, b, z): return "hypergeometric({},{},{})".format(repr(maxima(a)), repr(maxima(b)), repr(maxima(z))) def _maxima_lib_(self, a, b, z): return "hypergeometric({},{},{})".format(repr(maxima(a.operands())), repr(maxima(b.operands())), repr(maxima(z)))''' class EvaluationMethods: '''def operands(self, ex, a, b, z): return [(a), list(b), z]''' def sorted_parameters(self, ex, a, b, z): """ Return with parameters sorted in a canonical order EXAMPLES:: sage: hypergeometric([2,1,3], [5,4], 1/2).sorted_parameters() hypergeometric((1, 2, 3), (4, 5), 1/2) """ return hypergeometric(sorted(a), sorted(b), z) def eliminate_parameters(self, ex, a, b, z): """ Eliminate repeated parameters by pairwise cancellation of identical terms in a and b EXAMPLES:: sage: hypergeometric([1, 1, 2, 5], [5, 1, 4], ....:                1/2).eliminate_parameters() hypergeometric((1, 2), (4,), 1/2) sage: hypergeometric([x], [x], x).eliminate_parameters() hypergeometric((), (), x) sage: hypergeometric((5, 4), (4, 4), 3).eliminate_parameters() hypergeometric((5,), (4,), 3) """ aa = list(a) bb = list(b) p = pp = len(aa) q = qq = len(bb) i = 0 while i < qq and aa: bbb = bb[i] if bbb in aa: aa.remove(bbb) bb.remove(bbb) pp -= 1 qq -= 1 else: i += 1 if (pp, qq) != (p, q): return hypergeometric(aa, bb, z) return ex def is_termwise_finite(self, ex, a, b, z): """ Determine whether all terms of self are finite. Any infinite terms or ambiguous terms beyond the first zero, if one exists, are ignored. Ambiguous cases (where a term is the product of both zero and an infinity) are not considered finite. sage: hypergeometric([2], [3, 4], 5).is_termwise_finite() True sage: hypergeometric([2], [-3, 4], 5).is_termwise_finite() False sage: hypergeometric([-2], [-3, 4], 5).is_termwise_finite() True sage: hypergeometric([-3], [-3, 4], 5).is_termwise_finite()  # ambiguous False sage: hypergeometric([0], [-1], 5).is_termwise_finite() True sage: hypergeometric([0], [0], 5).is_termwise_finite()  # ambiguous False sage: hypergeometric([1], [2], Infinity).is_termwise_finite() False sage: hypergeometric([0], [0], Infinity).is_termwise_finite()  # ambiguous False sage: hypergeometric([0], [], Infinity).is_termwise_finite()  # ambiguous False """ if z == 0: return 0 not in b if abs(z) == Infinity: return False if abs(z) == Infinity: return False for bb in b: if bb in ZZ and bb <= 0: if any((aa in ZZ) and (bb < aa <= 0) for aa in a): continue return False return True def is_terminating(self, ex, a, b, z): """ Determine whether the series represented by self terminates after a finite number of terms, i.e. whether any of the numerator parameters are nonnegative integers (with no preceding nonnegative denominator parameters), or z = 0. If terminating, the series represents a polynomial of z EXAMPLES:: sage: hypergeometric([1, 2], [3, 4], x).is_terminating() False sage: hypergeometric([1, -2], [3, 4], x).is_terminating() True sage: hypergeometric([1, -2], [], x).is_terminating() True """ if z == 0: return True for aa in a: if (aa in ZZ) and (aa <= 0): return ex.is_termwise_finite() return False def is_absolutely_convergent(self, ex, a, b, z): """ Determine whether self converges absolutely as an infinite series. False is returned if not all terms are finite. EXAMPLES: Degree giving infinite radius of convergence:: sage: hypergeometric([2,3],[4,5],6).is_absolutely_convergent() True sage: hypergeometric([2,3],[-4,5],6).is_absolutely_convergent()  # undefined False sage: hypergeometric([2,3],[-4,5],Infinity).is_absolutely_convergent()  # undefined False Ordinary geometric series (unit radius of convergence):: sage: hypergeometric([1],[],1/2).is_absolutely_convergent() True sage: hypergeometric([1],[],2).is_absolutely_convergent() False sage: hypergeometric([1],[],1).is_absolutely_convergent() False sage: hypergeometric([1],[],-1).is_absolutely_convergent() False sage: hypergeometric([1],[],-1).n()    # Sum still exists 0.500000000000000 Degree $p = q+1$ (unit radius of convergence):: sage: hypergeometric([2,3],[4],6).is_absolutely_convergent() False sage: hypergeometric([2,3],[4],1).is_absolutely_convergent() False sage: hypergeometric([2,3],[5],1).is_absolutely_convergent() False sage: hypergeometric([2,3],[6],1).is_absolutely_convergent() True sage: hypergeometric([-2,3],[4],5).is_absolutely_convergent() True sage: hypergeometric([2,-3],[4],5).is_absolutely_convergent() True sage: hypergeometric([2,-3],[-4],5).is_absolutely_convergent() True sage: hypergeometric([2,-3],[-1],5).is_absolutely_convergent() False Degree giving zero radius of convergence:: sage: hypergeometric([1,2,3],[4],2).is_absolutely_convergent() False sage: hypergeometric([1,2,3],[4],1/2).is_absolutely_convergent() False sage: hypergeometric([1,2,-3],[4],1/2).is_absolutely_convergent()  # polynomial True """ p, q = len(a), len(b) if not ex.is_termwise_finite(): return False if p <= q: return True if ex.is_terminating(): return True if p == q + 1: if abs(z) < 1: return True if abs(z) == 1: if real_part(sum(b) - sum(a)) > 0: return True return False def terms(self, ex, a, b, z, n=None): """ Generate the terms of self (optionally only n terms). EXAMPLES:: sage: list(hypergeometric([-2, 1], [3, 4], x).terms()) [1, -1/6*x, 1/120*x^2] sage: list(hypergeometric([-2, 1], [3, 4], x).terms(2)) [1, -1/6*x] sage: list(hypergeometric([-2, 1], [3, 4], x).terms(0)) [] """ if n is None: n = Infinity t = Integer(1) k = 1 while k <= n: yield t for aa in a: t *= (aa + k - 1) for bb in b: t /= (bb + k - 1) t *= z if t == 0: break t /= k k += 1 def deflated(self, ex, a, b, z): """ Rewrite as a linear combination of functions of strictly lower degree by eliminating all parameters a[i] and b[j] such that a[i] = b[i] + m for nonnegative integer m. EXAMPLES:: sage: x = hypergeometric([5], [4], 3) sage: y = x.deflated() sage: y 7/4*hypergeometric((), (), 3) sage: x.n(); y.n() 35.1496896155784 35.1496896155784 sage: x = hypergeometric([6, 1], [3, 4, 5], 10) sage: y = x.deflated() sage: y hypergeometric((1,), (4, 5), 10) + 1/2*hypergeometric((2,), (5, 6), 10) + 1/12*hypergeometric((3,), (6, 7), 10) + 1/252*hypergeometric((4,), (7, 8), 10) sage: x.n(); y.n() 2.87893612686782 2.87893612686782 sage: x = hypergeometric([6, 7], [3, 4, 5], 10) sage: y = x.deflated() sage: y hypergeometric((), (5,), 10) + 5*hypergeometric((), (6,), 10) + 19/3*hypergeometric((), (7,), 10) + 181/63*hypergeometric((), (8,), 10) + 265/504*hypergeometric((), (9,), 10) + 25/648*hypergeometric((), (10,), 10) + 25/27216*hypergeometric((), (11,), 10) sage: x.n(); y.n() 63.0734110716969 63.0734110716969 """ return sum(map(prod, ex._deflated())) def _deflated(self, ex, a, b, z): new = ex.eliminate_parameters() aa = new.operands()[0].operands() bb = new.operands()[1].operands() for i, aaa in enumerate(aa): for j, bbb in enumerate(bb): m = aaa - bbb if m in ZZ and m > 0: aaaa = aa[:i] + aa[i + 1:] bbbb = bb[:j] + bb[j + 1:] terms = [] for k in xrange(m + 1): # TODO: could rewrite prefactors as recurrence term = binomial(m, k) for c in aaaa: term *= rising_factorial(c, k) for c in bbbb: term /= rising_factorial(c, k) term *= z ** k term /= rising_factorial(aaa - m, k) F = hypergeometric([c + k for c in aaaa], [c + k for c in bbbb], z) unique = [] counts = [] for c, f in F._deflated(): if f in unique: counts[unique.index(f)] += c else: unique.append(f) counts.append(c) Fterms = zip(counts, unique) terms += [(term*termG, G) for (termG, G) in Fterms] return terms return ((1, new),) hypergeometric = Hypergeometric() def closed_form(ex): """ Try to evaluate self in closed form using elementary (and other simple) functions. It may be necessary to call :meth:deflated first to find some closed forms. EXAMPLES:: sage: from sage.functions.hypergeometric import closed_form sage: var('a b c z') (a, b, c, z) sage: closed_form(hypergeometric([1], [], 1 + z)) -1/z sage: closed_form(hypergeometric([], [], 1 + z)) e^(z + 1) sage: closed_form(hypergeometric([], [1/2], 4)) cosh(4) sage: closed_form(hypergeometric([], [3/2], 4)) 1/4*sinh(4) sage: closed_form(hypergeometric([], [5/2], 4)) -3/64*sinh(4) + 3/16*cosh(4) sage: closed_form(hypergeometric([], [-3/2], 4)) -4*sinh(4) + 19/3*cosh(4) sage: closed_form(hypergeometric([-3, 1], [var('a')], z)) -6*z^3/((a + 1)*(a + 2)*a) + 6*z^2/((a + 1)*a) - 3*z/a + 1 sage: closed_form(hypergeometric([-3,1/3],[-4],z)) 7/162*z^3 + 1/9*z^2 + 1/4*z + 1 sage: closed_form(hypergeometric([],[],z)) e^z sage: closed_form(hypergeometric([a],[],z)) (-z + 1)^(-a) sage: closed_form(hypergeometric([1,1,2],[1,1],z)) (z - 1)^(-2) sage: closed_form(hypergeometric([2,3],[1],x)) -1/(x - 1)^3 + 3*x/(x - 1)^4 sage: closed_form(hypergeometric([1/2],[3/2],-5)) 1/10*sqrt(pi)*sqrt(5)*erf(sqrt(5)) sage: closed_form(hypergeometric([2],[5],3)) 4 sage: closed_form(hypergeometric([2],[5],5)) 48/625*e^5 + 612/625 sage: closed_form(hypergeometric([1/2,7/2],[3/2],z)) 1/sqrt(-z + 1) + 2/3*z/(-z + 1)^(3/2) + 1/5*z^2/(-z + 1)^(5/2) sage: closed_form(hypergeometric([1/2,1],[2],z)) -2*(sqrt(-z + 1) - 1)/z sage: closed_form(hypergeometric([1,1],[2],z)) -log(-z + 1)/z sage: closed_form(hypergeometric([1,1],[3],z)) -2*((z - 1)*log(-z + 1)/z - 1)/z TODO: Try Python int """ if ex.is_terminating(): return sum(ex.terms()) # necessary? new = ex.eliminate_parameters() #new = ex.deflated() def _closed_form(ex): a, b, z = ex.operands() a, b = a.operands(), b.operands() p, q = len(a), len(b) if z == 0: return Integer(1) if p == q == 0: return exp(z) if p == 1 and q == 0: return (1 - z) ** (-a[0]) if p == 0 and q == 1: # TODO: make this require only linear time def _0f1(b, z): F12 = cosh(2 * sqrt(z)) F32 = sinh(2 * sqrt(z)) / (2 * sqrt(z)) if 2 * b == 1: return F12 if 2 * b == 3: return F32 if 2 * b > 3: return ((b - 2) * (b - 1) / z * (_0f1(b - 2, z) - _0f1(b - 1, z))) if 2 * b < 1: return (_0f1(b + 1, z) + z / (b * (b + 1)) * _0f1(b + 2, z)) raise ValueError # Can evaluate 0F1 in terms of elementary functions when # the parameter is a half-integer if 2*b[0] in ZZ and b[0] not in ZZ: return _0f1(b[0], z) # Confluent hypergeometric function if p == 1 and q == 1: aa, bb = a[0], b[0] if aa * 2 == 1 and bb * 2 == 3: t = sqrt(-z) return sqrt(pi) / 2 * erf(t) / t if a == 1 and b == 2: return (exp(z) - 1) / z n, m = aa, bb if n in ZZ and m in ZZ and m > 0 and n > 0: rf = rising_factorial if m <= n: return (exp(z) * sum(rf(m - n, k) * (-z) ** k / factorial(k) / rf(m, k) for k in xrange(n - m + 1))) else: T = sum(rf(n - m + 1, k) * z ** k / (factorial(k) * rf(2 - m, k)) for k in xrange(m - n)) U = sum(rf(1 - n, k) * (-z) ** k / (factorial(k) * rf(2 - m, k)) for k in xrange(n)) return (factorial(m - 2) * rf(1 - m, n) * z ** (1 - m) / factorial(n - 1) * (T - exp(z) * U)) if p == 2 and q == 1: R12 = QQ('1/2') R32 = QQ('3/2') def _2f1(a, b, c, z): """ Evaluation of 2F1(a, b, c, z), assuming a, b, c positive integers or half-integers """ if b == c: return (1 - z) ** (-a) if a == c: return (1 - z) ** (-b) if a == 0 or b == 0: return Integer(1) if a > b: a, b = b, a if b >= 2: F1 = _2f1(a, b - 1, c, z) F2 = _2f1(a, b - 2, c, z) q = (b - 1) * (z - 1) return (((c - 2 * b + 2 + (b - a - 1) * z) * F1 + (b - c - 1) * F2) / q) if c > 2: # how to handle this case? if a - c + 1 == 0 or b - c + 1 == 0: raise NotImplementedError F1 = _2f1(a, b, c - 1, z) F2 = _2f1(a, b, c - 2, z) r1 = (c - 1) * (2 - c - (a + b - 2 * c + 3) * z) r2 = (c - 1) * (c - 2) * (1 - z) q = (a - c + 1) * (b - c + 1) * z return (r1 * F1 + r2 * F2) / q if (a, b, c) == (R12, 1, 2): return (2 - 2 * sqrt(1 - z)) / z if (a, b, c) == (1, 1, 2): return -log(1 - z) / z if (a, b, c) == (1, R32, R12): return (1 + z) / (1 - z) ** 2 if (a, b, c) == (1, R32, 2): return 2 * (1 / sqrt(1 - z) - 1) / z if (a, b, c) == (R32, 2, R12): return (1 + 3 * z) / (1 - z) ** 3 if (a, b, c) == (R32, 2, 1): return (2 + z) / (2 * (sqrt(1 - z) * (1 - z) ** 2)) if (a, b, c) == (2, 2, 1): return (1 + z) / (1 - z) ** 3 raise NotImplementedError aa, bb = a cc, = b if z == 1: return (gamma(cc) * gamma(cc - aa - bb) / gamma(cc - aa) / gamma(cc - bb)) if ((aa * 2) in ZZ and (bb * 2) in ZZ and (cc * 2) in ZZ and aa > 0 and bb > 0 and cc > 0): try: return _2f1(aa, bb, cc, z) except NotImplementedError: pass return sum([coeff * _closed_form(pfq) for coeff, pfq in new._deflated()]) #    return new
• ## sage/interfaces/maxima_lib.py

diff --git a/sage/interfaces/maxima_lib.py b/sage/interfaces/maxima_lib.py
 a ## Here we build dictionaries for operators needing special conversions. ratdisrep=EclObject("ratdisrep") mrat=EclObject("MRAT") mqapply=EclObject("MQAPPLY") max_li=EclObject("$LI") max_psi=EclObject("$PSI") max_array=EclObject("ARRAY") mdiff=EclObject("%DERIVATIVE") max_lambert_w=sage_op_dict[sage.functions.log.lambert_w] ratdisrep = EclObject("ratdisrep") mrat = EclObject("MRAT") mqapply = EclObject("MQAPPLY") max_li = EclObject("$LI") max_psi = EclObject("$PSI") max_hyper = EclObject("\$%F") max_array = EclObject("ARRAY") mdiff = EclObject("%DERIVATIVE") max_lambert_w = sage_op_dict[sage.functions.log.lambert_w] def mrat_to_sage(expr): r""" """ if caaadr(expr) == max_li: return sage.functions.log.polylog(max_to_sr(cadadr(expr)), max_to_sr(caddr(expr))) max_to_sr(caddr(expr))) if caaadr(expr) == max_psi: return sage.functions.other.psi(max_to_sr(cadadr(expr)), max_to_sr(caddr(expr))) max_to_sr(caddr(expr))) if caaadr(expr) == max_hyper: return sage.functions.hypergeometric.hypergeometric(mlist_to_sage(car(cdr(cdr(expr)))), mlist_to_sage(car(cdr(cdr(cdr(expr))))), max_to_sr(car(cdr(cdr(cdr(cdr(expr))))))) else: op=max_to_sr(cadr(expr)) max_args=cddr(expr) - expr - ECL object; a Maxima MLIST expression (i.e., a list) OUTPUT: a python list of converted expressions. OUTPUT: a Python list of converted expressions. EXAMPLES:: return EclObject(l) elif (op in special_sage_to_max): return EclObject(special_sage_to_max[op](*[sr_to_max(o) for o in expr.operands()])) elif op == tuple: return maxima(expr.operands()).ecl() elif not (op in sage_op_dict): # Maxima does some simplifications automatically by default # so calling maxima(expr) can change the structure of expr
• ## sage/symbolic/expression.pyx

diff --git a/sage/symbolic/expression.pyx b/sage/symbolic/expression.pyx
 a factorial(n) """ x = self x = x.simplify_hypergeometric() x = x.simplify_factorial() x = x.simplify_trig() x = x.simplify_rational() full_simplify = simplify_full def simplify_hypergeometric(self, algorithm='maxima'): """ EXAMPLES:: sage: hypergeometric((5,4), (4,1,2,3),x).simplify_hypergeometric() 1/144*x^2*hypergeometric((), (3, 4), x) + 1/3*x*hypergeometric((), (2, 3), x) + hypergeometric((), (1, 2), x) sage: (2*hypergeometric((), (), x)).simplify_hypergeometric() 2*e^x """ from sage.functions.hypergeometric import hypergeometric, closed_form from sage.calculus.calculus import maxima try: op = self.operator() except RuntimeError: return self ops = self.operands() if op == hypergeometric: if algorithm == 'maxima': return self.parent()(maxima.hgfred(map(lambda o: o.simplify_hypergeometric(algorithm), ops[0].operands()), map(lambda o: o.simplify_hypergeometric(algorithm), ops[1].operands()), ops[2].simplify_hypergeometric(algorithm))) elif algorithm == 'sage': return closed_form(hypergeometric(map(lambda o: o.simplify_hypergeometric(algorithm), ops[0].operands()), map(lambda o: o.simplify_hypergeometric(algorithm), ops[1].operands()), ops[2].simplify_hypergeometric(algorithm))) else: return NotImplementedError('unknown algorithm') if not op: return self return op(*map(lambda o: o.simplify_hypergeometric(algorithm), ops)) hypergeometric_simplify = simplify_hypergeometric def simplify_trig(self,expand=True): r""" Optionally expands and then employs identities such as
• ## sage/symbolic/expression_conversions.py

diff --git a/sage/symbolic/expression_conversions.py b/sage/symbolic/expression_conversions.py
 a return self.relation(ex, operator) elif isinstance(operator, FDerivativeOperator): return self.derivative(ex, operator) elif operator == tuple: return self.tuple(ex, operator) else: return self.composition(ex, operator) return "%s %s %s"%(self(ex.lhs()), self.relation_symbols[operator], self(ex.rhs())) def tuple(self, ex, operator): x = map(self, ex.operands()) X = ','.join(x) return '%s%s%s'%(self.interface._left_list_delim(), X, self.interface._right_list_delim()) def derivative(self, ex, operator): """ EXAMPLES:: 1.41421356237309... """ f = operator g = map(self, ex.operands()) try: g = map(self, ex.operands()) except AttributeError, error: from sage.functions.hypergeometric import hypergeometric if f is hypergeometric: ops = ex.operands() return hypergeometric(*[map(lambda x: self.pyobject(None, x), ops[0]), map(lambda x: self.pyobject(None, x), ops[1]), self(ops[2])]) else: raise AttributeError(error) try: return f(*g) except TypeError: return abs(*g) # special case else: return self.ff.fast_float_func(f, *g) def fast_float(ex, *vars): """ Returns an object which provides fast floating point evaluation of