From 7339c34483b1a3296e6aeacb3f66259e2431934e Mon Sep 17 00:00:00 2001 From: Daiki Takahashi Date: Thu, 16 May 2024 22:38:26 +0900 Subject: [PATCH 1/2] Simple fixes on `holonomic/holonomic.py`. - Shallow indentation. - Changed to comprehensions. - Simplified (`if cond: return True; else: return False` to `return cond` etc). --- sympy/holonomic/holonomic.py | 346 ++++++++++++----------------------- 1 file changed, 119 insertions(+), 227 deletions(-) diff --git a/sympy/holonomic/holonomic.py b/sympy/holonomic/holonomic.py index 15e0d6daf2ad..67bf9a85d60c 100644 --- a/sympy/holonomic/holonomic.py +++ b/sympy/holonomic/holonomic.py @@ -151,10 +151,8 @@ def __str__(self): __repr__ = __str__ def __eq__(self, other): - if self.base == other.base and self.gen_symbol == other.gen_symbol: - return True - else: - return False + return self.base == other.base and \ + self.gen_symbol == other.gen_symbol class DifferentialOperator: @@ -233,25 +231,18 @@ def __mul__(self, other): """ listofself = self.listofpoly - - if not isinstance(other, DifferentialOperator): - if not isinstance(other, self.parent.base.dtype): - listofother = [self.parent.base.from_sympy(sympify(other))] - - else: - listofother = [other] - else: + if isinstance(other, DifferentialOperator): listofother = other.listofpoly + elif isinstance(other, self.parent.base.dtype): + listofother = [other] + else: + listofother = [self.parent.base.from_sympy(sympify(other))] # multiplies a polynomial `b` with a list of polynomials def _mul_dmp_diffop(b, listofother): if isinstance(listofother, list): - sol = [] - for i in listofother: - sol.append(i * b) - return sol - else: - return [b * listofother] + return [i * b for i in listofother] + return [b * listofother] sol = _mul_dmp_diffop(listofself[0], listofother) @@ -284,10 +275,7 @@ def __rmul__(self, other): if not isinstance(other, self.parent.base.dtype): other = (self.parent.base).from_sympy(sympify(other)) - sol = [] - for j in self.listofpoly: - sol.append(other * j) - + sol = [other * j for j in self.listofpoly] return DifferentialOperator(sol, self.parent) def __add__(self, other): @@ -296,17 +284,13 @@ def __add__(self, other): sol = _add_lists(self.listofpoly, other.listofpoly) return DifferentialOperator(sol, self.parent) + list_self = self.listofpoly + if not isinstance(other, self.parent.base.dtype): + list_other = [((self.parent).base).from_sympy(sympify(other))] else: - list_self = self.listofpoly - if not isinstance(other, self.parent.base.dtype): - list_other = [((self.parent).base).from_sympy(sympify(other))] - else: - list_other = [other] - sol = [] - sol.append(list_self[0] + list_other[0]) - sol += list_self[1:] - - return DifferentialOperator(sol, self.parent) + list_other = [other] + sol = [list_self[0] + list_other[0]] + list_self[1:] + return DifferentialOperator(sol, self.parent) __radd__ = __add__ @@ -330,18 +314,16 @@ def __pow__(self, n): # if self is `Dx` if self.listofpoly == self.parent.derivative_operator.listofpoly: - sol = [self.parent.base.zero]*n - sol.append(self.parent.base.one) + sol = [self.parent.base.zero]*n + [self.parent.base.one] return DifferentialOperator(sol, self.parent) # the general case - else: - if n % 2 == 1: - powreduce = self**(n - 1) - return powreduce * self - elif n % 2 == 0: - powreduce = self**(n / 2) - return powreduce * powreduce + if n % 2 == 1: + powreduce = self**(n - 1) + return powreduce * self + elif n % 2 == 0: + powreduce = self**(n / 2) + return powreduce * powreduce def __str__(self): listofpoly = self.listofpoly @@ -372,18 +354,10 @@ def __str__(self): def __eq__(self, other): if isinstance(other, DifferentialOperator): - if self.listofpoly == other.listofpoly and self.parent == other.parent: - return True - else: - return False - else: - if self.listofpoly[0] == other: - for i in self.listofpoly[1:]: - if i is not self.parent.base.zero: - return False - return True - else: - return False + return self.listofpoly == other.listofpoly and \ + self.parent == other.parent + return self.listofpoly[0] == other and \ + all(i is self.parent.base.zero for i in self.listofpoly[1:]) def is_singular(self, x0): """ @@ -564,10 +538,8 @@ def _singularics_to_ord(self): b = self.y0[a] if len(self.y0) == 1 and a == int(a) and a > 0: - y0 = [] a = int(a) - for i in range(a): - y0.append(S.Zero) + y0 = [S.Zero] * a y0 += [j * factorial(a + i) for i, j in enumerate(b)] return HolonomicFunction(self.annihilator, self.x, self.x0, y0) @@ -668,26 +640,19 @@ def __add__(self, other): y0 = [a + b for a, b in zip(y1, y2)] return HolonomicFunction(sol, self.x, self.x0, y0) - else: - # change the initial conditions to a same point - selfat0 = self.annihilator.is_singular(0) - otherat0 = other.annihilator.is_singular(0) - - if self.x0 == 0 and not selfat0 and not otherat0: - return self + other.change_ics(0) - - elif other.x0 == 0 and not selfat0 and not otherat0: - return self.change_ics(0) + other + # change the initial conditions to a same point + selfat0 = self.annihilator.is_singular(0) + otherat0 = other.annihilator.is_singular(0) + if self.x0 == 0 and not selfat0 and not otherat0: + return self + other.change_ics(0) + if other.x0 == 0 and not selfat0 and not otherat0: + return self.change_ics(0) + other - else: - selfatx0 = self.annihilator.is_singular(self.x0) - otheratx0 = other.annihilator.is_singular(self.x0) - - if not selfatx0 and not otheratx0: - return self + other.change_ics(self.x0) - - else: - return self.change_ics(other.x0) + other + selfatx0 = self.annihilator.is_singular(self.x0) + otheratx0 = other.annihilator.is_singular(self.x0) + if not selfatx0 and not otheratx0: + return self + other.change_ics(self.x0) + return self.change_ics(other.x0) + other if self.x0 != other.x0: return HolonomicFunction(sol, self.x) @@ -929,19 +894,11 @@ def diff(self, *args, **kwargs): return HolonomicFunction(sol, self.x, self.x0, y0) def __eq__(self, other): - if self.annihilator == other.annihilator: - if self.x == other.x: - if self._have_init_cond() and other._have_init_cond(): - if self.x0 == other.x0 and self.y0 == other.y0: - return True - else: - return False - else: - return True - else: - return False - else: + if self.annihilator != other.annihilator or self.x != other.x: return False + if self._have_init_cond() and other._have_init_cond(): + return self.x0 == other.x0 and self.y0 == other.y0 + return True def __mul__(self, other): ann_self = self.annihilator @@ -954,14 +911,9 @@ def __mul__(self, other): if not self._have_init_cond(): return self - else: - y0 = _extend_y0(self, ann_self.order) - y1 = [] - - for j in y0: - y1.append((Poly.new(j, self.x) * other).rep) - - return HolonomicFunction(ann_self, self.x, self.x0, y1) + y0 = _extend_y0(self, ann_self.order) + y1 = [(Poly.new(j, self.x) * other).rep for j in y0] + return HolonomicFunction(ann_self, self.x, self.x0, y1) if self.annihilator.parent.base != other.annihilator.parent.base: a, b = self.unify(other) @@ -969,20 +921,14 @@ def __mul__(self, other): ann_other = other.annihilator - list_self = [] - list_other = [] - a = ann_self.order b = ann_other.order R = ann_self.parent.base K = R.get_field() - for j in ann_self.listofpoly: - list_self.append(K.new(j.to_list())) - - for j in ann_other.listofpoly: - list_other.append(K.new(j.to_list())) + list_self = [K.new(j.to_list()) for j in ann_self.listofpoly] + list_other = [K.new(j.to_list()) for j in ann_other.listofpoly] # will be used to reduce the degree self_red = [-list_self[i] / list_self[a] for i in range(a)] @@ -1016,19 +962,19 @@ def __mul__(self, other): # reduce the terms to lower power using annihilators of f, g for i in range(a + 1): - if not coeff_mul[i][b].is_zero: - for j in range(b): - coeff_mul[i][j] += other_red[j] * \ - coeff_mul[i][b] - coeff_mul[i][b] = K.zero + if coeff_mul[i][b].is_zero: + continue + for j in range(b): + coeff_mul[i][j] += other_red[j] * coeff_mul[i][b] + coeff_mul[i][b] = K.zero # not d2 + 1, as that is already covered in previous loop for j in range(b): - if not coeff_mul[a][j] == 0: - for i in range(a): - coeff_mul[i][j] += self_red[i] * \ - coeff_mul[a][j] - coeff_mul[a][j] = K.zero + if coeff_mul[a][j] == 0: + continue + for i in range(a): + coeff_mul[i][j] += self_red[i] * coeff_mul[a][j] + coeff_mul[a][j] = K.zero lin_sys_elements.append([coeff_mul[i][j] for i in range(a) for j in range(b)]) lin_sys = DomainMatrix(lin_sys_elements, (len(lin_sys_elements), a*b), K).transpose() @@ -1069,26 +1015,19 @@ def __mul__(self, other): return HolonomicFunction(sol_ann, self.x, self.x0, y0) # if the points are different, consider one - else: + selfat0 = self.annihilator.is_singular(0) + otherat0 = other.annihilator.is_singular(0) - selfat0 = self.annihilator.is_singular(0) - otherat0 = other.annihilator.is_singular(0) + if self.x0 == 0 and not selfat0 and not otherat0: + return self * other.change_ics(0) + if other.x0 == 0 and not selfat0 and not otherat0: + return self.change_ics(0) * other - if self.x0 == 0 and not selfat0 and not otherat0: - return self * other.change_ics(0) - - elif other.x0 == 0 and not selfat0 and not otherat0: - return self.change_ics(0) * other - - else: - selfatx0 = self.annihilator.is_singular(self.x0) - otheratx0 = other.annihilator.is_singular(self.x0) - - if not selfatx0 and not otheratx0: - return self * other.change_ics(self.x0) - - else: - return self.change_ics(other.x0) * other + selfatx0 = self.annihilator.is_singular(self.x0) + otheratx0 = other.annihilator.is_singular(self.x0) + if not selfatx0 and not otheratx0: + return self * other.change_ics(self.x0) + return self.change_ics(other.x0) * other if self.x0 != other.x0: return HolonomicFunction(sol_ann, self.x) @@ -1114,12 +1053,8 @@ def __mul__(self, other): for i in y1: for j in y2: k = min(len(y1[i]), len(y2[j])) - c = [] - for a in range(k): - s = S.Zero - for b in range(a + 1): - s += y1[i][b] * y2[j][a - b] - c.append(s) + c = [sum((y1[i][b] * y2[j][a - b] for b in range(a + 1)), + start=S.Zero) for a in range(k)] if not i + j in y0: y0[i + j] = c else: @@ -1165,13 +1100,12 @@ def __pow__(self, n): return HolonomicFunction(Dx, self.x, S.Zero, [S.One]) if n == 1: return self - else: - if n % 2 == 1: - powreduce = self**(n - 1) - return powreduce * self - elif n % 2 == 0: - powreduce = self**(n / 2) - return powreduce * powreduce + if n % 2 == 1: + powreduce = self**(n - 1) + return powreduce * self + if n % 2 == 0: + powreduce = self**(n / 2) + return powreduce * powreduce def degree(self): """ @@ -1344,10 +1278,9 @@ def to_sequence(self, lb=True): # an appropriate shift of the recurrence for j in range(lower, upper + 1): if j in keylist: - temp = S.Zero - for k in dict1.keys(): - if k[0] == j: - temp += dict1[k].subs(n, n - lower) + temp = sum((v.subs(n, n - lower) + for k, v in dict1.items() if k[0] == j), + start=S.Zero) sol.append(temp) else: sol.append(S.Zero) @@ -1366,11 +1299,8 @@ def to_sequence(self, lb=True): order += smallest_n y0 = _extend_y0(self, order) - u0 = [] - # u(n) = y^n(0)/factorial(n) - for i, j in enumerate(y0): - u0.append(j / factorial(i)) + u0 = [j / factorial(i) for i, j in enumerate(y0)] # if sufficient conditions can't be computed then # try to use the series method i.e. @@ -1458,21 +1388,19 @@ def _frobenius(self, lb=True): grp = [] for i in reals: - intdiff = False if len(grp) == 0: grp.append([i]) continue for j in grp: if int_valued(j[0] - i): j.append(i) - intdiff = True break - if not intdiff: + else: grp.append([i]) # True if none of the roots differ by an integer i.e. # each element in group have only one member - independent = True if all(len(i) == 1 for i in grp) else False + independent = all(len(i) == 1 for i in grp) allpos = all(i >= 0 for i in reals) allint = all(int_valued(i) for i in reals) @@ -1493,10 +1421,7 @@ def _frobenius(self, lb=True): rootstoconsider = [i[0] for i in grp] + [j[0] for j in compl] elif not allint: - rootstoconsider = [] - for i in reals: - if not int(i) == i: - rootstoconsider.append(i) + rootstoconsider = [i for i in reals if not int(i) == i] elif not allpos: @@ -1504,10 +1429,7 @@ def _frobenius(self, lb=True): rootstoconsider = [min(reals)] else: - posroots = [] - for i in reals: - if i >= 0: - posroots.append(i) + posroots = [i for i in reals if i >= 0] rootstoconsider = [min(posroots)] n = Symbol('n', integer=True) @@ -1550,10 +1472,9 @@ def _frobenius(self, lb=True): for j in range(lower, upper + 1): if j in keylist: - temp = S.Zero - for k in dict1.keys(): - if k[0] == j: - temp += dict1[k].subs(n, n - lower) + temp = sum((v.subs(n, n - lower) + for k, v in dict1.items() if k[0] == j), + start=S.Zero) sol.append(temp) else: sol.append(S.Zero) @@ -1580,8 +1501,7 @@ def _frobenius(self, lb=True): y0 = _extend_y0(self, order + int(p)) # u(n) = y^n(0)/factorial(n) if len(y0) > int(p): - for i in range(int(p), len(y0)): - u0.append(y0[i] / factorial(i)) + u0 = [y0[i] / factorial(i) for i in range(int(p), len(y0))] if len(u0) < order: @@ -1698,10 +1618,7 @@ def series(self, n=6, coefficient=False, order=True, _recur=None): constantpower = recurrence[0][1] recurrence = recurrence[0][0] else: - sol = [] - for i in recurrence: - sol.append(self.series(_recur=i)) - return sol + return [self.series(_recur=i) for i in recurrence] n = n - int(constantpower) l = len(recurrence.u0) - 1 @@ -1711,31 +1628,22 @@ def series(self, n=6, coefficient=False, order=True, _recur=None): seq_dmp = recurrence.recurrence.listofpoly R = recurrence.recurrence.parent.base K = R.get_field() - seq = [] - - for i, j in enumerate(seq_dmp): - seq.append(K.new(j.to_list())) - + seq = [K.new(j.to_list()) for j in seq_dmp] sub = [-seq[i] / seq[k] for i in range(k)] sol = list(recurrence.u0) - if l + 1 >= n: - pass - else: + if l + 1 < n: # use the initial conditions to find the next term for i in range(l + 1 - k, n - k): - coeff = S.Zero - for j in range(k): - if i + j >= 0: - coeff += DMFsubs(sub[j], i) * sol[i + j] + coeff = sum((DMFsubs(sub[j], i) * sol[i + j] + for j in range(k) if i + j >= 0), start=S.Zero) sol.append(coeff) if coefficient: return sol - ser = S.Zero - for i, j in enumerate(sol): - ser += x**(i + constantpower) * j + ser = sum((x**(i + constantpower) * j for i, j in enumerate(sol)), + start=S.Zero) if order: ser += Order(x**(n + int(constantpower)), x) if x0 != 0: @@ -1772,7 +1680,7 @@ def _pole_degree(poly): for i, j in enumerate(list_coeff): listofdmp = j.all_coeffs() degree = len(listofdmp) - 1 - if - i - b <= 0 and degree - i - b >= 0: + if 0 <= i + b <= degree: s = s + listofdmp[degree - i - b] * y y *= R.from_sympy(x - i) @@ -1862,9 +1770,7 @@ def change_x(self, z): dom = self.annihilator.parent.base.dom R = dom.old_poly_ring(z) parent, _ = DifferentialOperators(R, 'Dx') - sol = [] - for j in self.annihilator.listofpoly: - sol.append(R(j.to_list())) + sol = [R(j.to_list()) for j in self.annihilator.listofpoly] sol = DifferentialOperator(sol, parent) return HolonomicFunction(sol, z, self.x0, self.y0) @@ -1989,14 +1895,7 @@ def to_hyper(self, as_list=False, _recur=None): raise NotImplementedError("Can't compute sufficient Initial Conditions") # check if the recurrence represents a hypergeometric series - is_hyper = True - - for i in range(1, len(r.listofpoly)-1): - if r.listofpoly[i] != r.parent.base.zero: - is_hyper = False - break - - if not is_hyper: + if any(i != r.parent.base.zero for i in r.listofpoly[1:-1]): raise NotHyperSeriesError(self, self.x0) a = r.listofpoly[0] @@ -2569,8 +2468,7 @@ def _hyper_to_meijerg(func): ap = func.ap bq = func.bq - ispoly = any(i <= 0 and int(i) == i for i in ap) - if ispoly: + if any(i <= 0 and int(i) == i for i in ap): return hyperexpand(func) z = func.args[2] @@ -2627,26 +2525,21 @@ def _extend_y0(Holonomic, n): if len(y0) < a or n <= len(y0): return y0 - else: - list_red = [-listofpoly[i] / listofpoly[a] - for i in range(a)] - if len(y0) > a: - y1 = [y0[i] for i in range(a)] - else: - y1 = list(y0) - for i in range(n - a): - sol = 0 - for a, b in zip(y1, list_red): - r = DMFsubs(b, Holonomic.x0) - if not getattr(r, 'is_finite', True): - return y0 - if isinstance(r, (PolyElement, FracElement)): - r = r.as_expr() - sol += a * r - y1.append(sol) - list_red = _derivate_diff_eq(list_red, K) - - return y0 + y1[len(y0):] + list_red = [-listofpoly[i] / listofpoly[a] + for i in range(a)] + y1 = y0[:min(len(y0), a)] + for _ in range(n - a): + sol = 0 + for a, b in zip(y1, list_red): + r = DMFsubs(b, Holonomic.x0) + if not getattr(r, 'is_finite', True): + return y0 + if isinstance(r, (PolyElement, FracElement)): + r = r.as_expr() + sol += a * r + y1.append(sol) + list_red = _derivate_diff_eq(list_red, K) + return y0 + y1[len(y0):] def DMFdiff(frac, K): @@ -2737,10 +2630,9 @@ def _convert_poly_rat_alg(func, x, x0=0, y0=None, lenics=None, domain=QQ, initco for i, j in enumerate(reversed(rep)): if j == 0: continue - else: - coeff = list(reversed(rep))[i:] - indicial = i - break + coeff = list(reversed(rep))[i:] + indicial = i + break for i, j in enumerate(coeff): if isinstance(j, (PolyElement, FracElement)): coeff[i] = j.as_expr() From 77cc5ec16edca6a6fe77f50ae3bb529f4d1de49e Mon Sep 17 00:00:00 2001 From: haru-44 <36563693+haru-44@users.noreply.github.com> Date: Fri, 17 May 2024 23:48:51 +0900 Subject: [PATCH 2/2] Update sympy/holonomic/holonomic.py Co-authored-by: S.Y. Lee --- sympy/holonomic/holonomic.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/sympy/holonomic/holonomic.py b/sympy/holonomic/holonomic.py index 67bf9a85d60c..9ac4b6515e05 100644 --- a/sympy/holonomic/holonomic.py +++ b/sympy/holonomic/holonomic.py @@ -321,9 +321,8 @@ def __pow__(self, n): if n % 2 == 1: powreduce = self**(n - 1) return powreduce * self - elif n % 2 == 0: - powreduce = self**(n / 2) - return powreduce * powreduce + powreduce = self**(n // 2) + return powreduce * powreduce def __str__(self): listofpoly = self.listofpoly