diff --git a/mathics/builtin/arithmetic.py b/mathics/builtin/arithmetic.py index 9ae9917084..558a1abc5d 100644 --- a/mathics/builtin/arithmetic.py +++ b/mathics/builtin/arithmetic.py @@ -1874,7 +1874,7 @@ class Sum(_IterationFunction, SympyFunction): = 0 >> (-1 + a^n) Sum[a^(k n), {k, 0, m-1}] // Simplify - = Piecewise[{{m (-1 + a ^ n), a ^ n == 1}, {-1 + (a ^ n) ^ m, True}}] + = -1 + (a ^ n) ^ m Infinite sums: >> Sum[1 / 2 ^ i, {i, 1, Infinity}] diff --git a/mathics/builtin/inout.py b/mathics/builtin/inout.py index 06124a4f0c..47c4b8179e 100644 --- a/mathics/builtin/inout.py +++ b/mathics/builtin/inout.py @@ -1974,6 +1974,35 @@ def apply(self, expr, evaluation) -> Expression: return self.apply_python(expr, evaluation) +class SympyForm(Builtin): + """ +
+
'SympyForm[$expr$]' +
returns an Sympy $expr$ in Python. Sympy is used internally + to implement a number of Mathics functions, like Simplify. +
+ + >> SympyForm[Pi^2] + = pi**2 + >> E^2 + 3E // SympyForm + = exp(2) + 3*E + """ + + def apply_sympy(self, expr, evaluation) -> Expression: + 'MakeBoxes[expr_, SympyForm]' + + try: + # from trepan.api import debug; debug() + sympy_equivalent = expr.to_sympy() + except: + return + return StringFromPython(sympy_equivalent) + + def apply(self, expr, evaluation) -> Expression: + "SympyForm[expr_]" + return self.apply_sympy(expr, evaluation) + + class TeXForm(Builtin): r"""
diff --git a/mathics/builtin/linalg.py b/mathics/builtin/linalg.py index 2f80cbdebb..5a1ef63bd5 100644 --- a/mathics/builtin/linalg.py +++ b/mathics/builtin/linalg.py @@ -58,24 +58,24 @@ class Tr(Builtin):
'Tr[$m$]'
computes the trace of the matrix $m$.
- + >> Tr[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}] = 15 - + Symbolic trace: >> Tr[{{a, b, c}, {d, e, f}, {g, h, i}}] = a + e + i """ - + messages = { 'matsq': "The matrix `1` is not square." } - + #TODO: generalize to vectors and higher-rank tensors, and allow function arguments for application - + def apply(self, m, evaluation): 'Tr[m_]' - + matrix = to_sympy_matrix(m) if matrix is None or matrix.cols != matrix.rows or matrix.cols == 0: return evaluation.message('Tr', 'matsq', m) @@ -701,7 +701,7 @@ class Eigenvalues(Builtin): = {-1, 1, 2} >> Eigenvalues[{{Cos[theta],Sin[theta],0},{-Sin[theta],Cos[theta],0},{0,0,1}}] // Sort - = {1, Cos[theta] + Sqrt[-1 + Cos[theta] ^ 2], Cos[theta] - Sqrt[-1 + Cos[theta] ^ 2]} + = {1, Cos[theta] + Sqrt[(-1 + Cos[theta]) (1 + Cos[theta])], Cos[theta] - Sqrt[(-1 + Cos[theta]) (1 + Cos[theta])]} >> Eigenvalues[{{7, 1}, {-4, 3}}] = {5, 5} @@ -731,7 +731,7 @@ def apply(self, m, evaluation): eigenvalues.sort(key=lambda v: (abs(v[0]), - re(v[0]), - im(v[0])), reverse=True) - eigenvalues = [from_sympy(v) for (v, c) in eigenvalues + eigenvalues = [from_sympy(v) for (v, c) in eigenvalues for _ in range(c)] # Sort the eigenvalues in an arbitrary yet deterministic order else: diff --git a/mathics/builtin/lists.py b/mathics/builtin/lists.py index be2342d999..7921004fff 100644 --- a/mathics/builtin/lists.py +++ b/mathics/builtin/lists.py @@ -13,7 +13,7 @@ PartError, PartDepthError, PartRangeError, Predefined, SympyFunction) from mathics.builtin.scoping import dynamic_scoping from mathics.builtin.base import MessageException, NegativeIntegerException, CountableInteger -from mathics.core.expression import Expression, String, Symbol, Integer, Number, Real, strip_context, from_python +from mathics.core.expression import Expression, String, Symbol, SymbolNull, Integer, Number, Real, strip_context, from_python from mathics.core.expression import min_prec, machine_precision from mathics.core.expression import structure from mathics.core.evaluation import BreakInterrupt, ContinueInterrupt, ReturnInterrupt @@ -37,8 +37,8 @@ def deletecases_with_levelspec(expr, pattern, evaluation, levelspec=1, n=-1 ): """ This function walks the expression `expr` and deleting occurrencies of `pattern` - - If levelspec specifies a number, only those positions with `levelspec` "coordinates" are return. By default, it just return occurences in the first level. + + If levelspec specifies a number, only those positions with `levelspec` "coordinates" are return. By default, it just return occurences in the first level. If a tuple (nmin, nmax) is provided, it just return those occurences with a number of "coordinates" between nmin and nmax. n indicates the number of occurrences to return. By default, it returns all the occurences. """ @@ -74,7 +74,7 @@ def deletecases_with_levelspec(expr, pattern, evaluation, levelspec=1, n=-1 ): break idx = curr_index[-1] changed = changed or changed_marks[-1][idx] - changed_marks[-1][idx] = changed + changed_marks[-1][idx] = changed if changed: head = tree[-1][curr_index[-1]].get_head() tree[-1][idx] = Expression(head, *leaves) @@ -83,12 +83,12 @@ def deletecases_with_levelspec(expr, pattern, evaluation, levelspec=1, n=-1 ): curr_index[-1] = curr_index[-1] + 1 continue curr_leave = tree[-1][curr_index[-1]] - if (match(curr_leave, evaluation) and - (len(curr_index) > lsmin) ): + if (match(curr_leave, evaluation) and + (len(curr_index) > lsmin) ): tree[-1][curr_index[-1]] = nothing changed_marks[-1][curr_index[-1]] = True curr_index[-1] = curr_index[-1] + 1 - n = n - 1 + n = n - 1 continue if curr_leave.is_atom() or lsmax == len(curr_index): curr_index[-1] = curr_index[-1] + 1 @@ -104,7 +104,7 @@ def find_matching_indices_with_levelspec(expr, pattern, evaluation,levelspec=1,n """ This function walks the expression `expr` looking for a pattern `pattern` and returns the positions of each occurence. - If levelspec specifies a number, only those positions with `levelspec` "coordinates" are return. By default, it just return occurences in the first level. + If levelspec specifies a number, only those positions with `levelspec` "coordinates" are return. By default, it just return occurences in the first level. If a tuple (nmin, nmax) is provided, it just return those occurences with a number of "coordinates" between nmin and nmax. n indicates the number of occurrences to return. By default, it returns all the occurences. """ @@ -130,11 +130,11 @@ def find_matching_indices_with_levelspec(expr, pattern, evaluation,levelspec=1,n curr_index[-1] = curr_index[-1] + 1 continue curr_leave = tree[-1][curr_index[-1]] - if (match(curr_leave, evaluation) and + if (match(curr_leave, evaluation) and (len(curr_index) >= lsmin) ): found.append([from_python(i) for i in curr_index]) curr_index[-1] = curr_index[-1] + 1 - n = n - 1 + n = n - 1 continue if curr_leave.is_atom() or lsmax == len(curr_index): curr_index[-1] = curr_index[-1] + 1 @@ -1840,11 +1840,11 @@ class DeleteCases(Builtin):
returns the elements of $list$ that do not match $pattern$.
'DeleteCases[$list$, $pattern$, $levelspec$]' -
removes all parts of $list on levels specified by $levelspec$ +
removes all parts of $list on levels specified by $levelspec$ that match pattern (not fully implemented).
'DeleteCases[$list$, $pattern$, $levelspec$, $n$]' -
removes the first $n$ parts of $list$ that match $pattern$. +
removes the first $n$ parts of $list$ that match $pattern$. >> DeleteCases[{a, 1, 2.5, "string"}, _Integer|_Real] @@ -1861,16 +1861,16 @@ class DeleteCases(Builtin): messages = {'level': 'Level specification `1` is not of the form n, {n}, or {m, n}.', 'innf': 'Non-negative integer or Infinity expected at position 4 in `1`', } - + # def apply(self, items, pattern, evaluation): # 'DeleteCases[items_, pattern_]' - # return self.apply_ls_n(items, pattern, Integer(1), Symbol("System`Null"), evaluation) + # return self.apply_ls_n(items, pattern, Integer(1), SymbolNull, evaluation) # def apply_ls(self, items, pattern, levelspec, evaluation): # 'DeleteCases[items_, pattern_, levelspec_]' - # return self.apply_ls_n(items, pattern, levelspec, Symbol("System`Null"), evaluation) - + # return self.apply_ls_n(items, pattern, levelspec, SymbolNull, evaluation) + def apply_ls_n(self, items, pattern, levelspec, n, evaluation): 'DeleteCases[items_, pattern_, levelspec_:1, n_:System`Infinity]' @@ -1893,8 +1893,8 @@ def apply_ls_n(self, items, pattern, levelspec, n, evaluation): evaluation.message('DeleteCases','innf',Expression("DeleteCases", items, pattern, levelspec, n)) else: evaluation.message('DeleteCases','innf',Expression("DeleteCases", items, pattern, levelspec, n)) - return Symbol('System`Null') - + return SymbolNull + if levelspec[0] !=1 or levelspec[1] !=1: return deletecases_with_levelspec(items, pattern, evaluation, levelspec, n) else: diff --git a/mathics/builtin/numbertheory.py b/mathics/builtin/numbertheory.py index 5a54002cd7..baac72b8e4 100644 --- a/mathics/builtin/numbertheory.py +++ b/mathics/builtin/numbertheory.py @@ -294,7 +294,7 @@ class Divisors(Builtin): # TODO: support GaussianIntegers # e.g. Divisors[2, GaussianIntegers -> True] - + attributes = ('Listable',) def apply(self, n, evaluation): @@ -369,102 +369,102 @@ class MantissaExponent(Builtin): >> MantissaExponent[125., 2] = {0.976563, 7} - + >> MantissaExponent[10, b] = MantissaExponent[10, b] - + #> MantissaExponent[E, Pi] = {E / Pi, 1} - + #> MantissaExponent[Pi, Pi] = {1 / Pi, 2} - + #> MantissaExponent[5/2 + 3, Pi] = {11 / (2 Pi ^ 2), 2} - + #> MantissaExponent[b] = MantissaExponent[b] - + #> MantissaExponent[17, E] = {17 / E ^ 3, 3} - + #> MantissaExponent[17., E] = {0.84638, 3} - + #> MantissaExponent[Exp[Pi], 2] = {E ^ Pi / 32, 5} - + #> MantissaExponent[3 + 2 I, 2] : The value 3 + 2 I is not a real number = MantissaExponent[3 + 2 I, 2] - + #> MantissaExponent[25, 0.4] : Base 0.4 is not a real number greater than 1. = MantissaExponent[25, 0.4] - + #> MantissaExponent[0.0000124] = {0.124, -4} - + #> MantissaExponent[0.0000124, 2] = {0.812646, -16} - + #> MantissaExponent[0] = {0, 0} - + #> MantissaExponent[0, 2] = {0, 0} """ - + attributes = ('Listable',) - + rules = { 'MantissaExponent[0]': '{0, 0}', 'MantissaExponent[0, n_]': '{0, 0}', } - + messages = { 'realx': 'The value `1` is not a real number', 'rbase': 'Base `1` is not a real number greater than 1.', } - + def apply(self, n, b, evaluation): 'MantissaExponent[n_, b_]' # Handle Input with special cases such as PI and E n_sympy, b_sympy = n.to_sympy(), b.to_sympy() - + expr = Expression('MantissaExponent', n, b) - + if isinstance(n.to_python(), complex): evaluation.message('MantissaExponent', 'realx', n) return expr - + if n_sympy.is_constant(): temp_n = Expression('N', n).evaluate(evaluation) py_n = temp_n.to_python() else: return expr - + if b_sympy.is_constant(): temp_b = Expression('N', b).evaluate(evaluation) py_b = temp_b.to_python() else: return expr - + if not py_b > 1: evaluation.message('MantissaExponent', 'rbase', b) return expr - + base_exp = int(mpmath.log(py_n, py_b)) - + exp = (base_exp + 1) if base_exp >= 0 else base_exp return Expression('List', Expression('Divide', n , b ** exp), exp) - + def apply_2(self, n, evaluation): 'MantissaExponent[n_]' n_sympy = n.to_sympy() expr = Expression('MantissaExponent', n) - + if isinstance(n.to_python(), complex): evaluation.message('MantissaExponent', 'realx', n) return expr @@ -474,12 +474,12 @@ def apply_2(self, n, evaluation): py_n = temp_n.to_python() else: return expr - - base_exp = int(mpmath.log10(py_n)) + + base_exp = int(mpmath.log10(py_n)) exp = (base_exp + 1) if base_exp >= 0 else base_exp - + return Expression('List', Expression('Divide', n , (10 ** exp)), exp) - + def _fractional_part(self, n, expr, evaluation): n_sympy = n.to_sympy() if n_sympy.is_constant(): @@ -491,9 +491,9 @@ def _fractional_part(self, n, expr, evaluation): result = n - negative_integer_part else: return expr - + return from_python(result) - + class FractionalPart(Builtin): """
@@ -509,34 +509,34 @@ class FractionalPart(Builtin): #> FractionalPart[b] = FractionalPart[b] - + #> FractionalPart[{-2.4, -2.5, -3.0}] = {-0.4, -0.5, 0.} - + #> FractionalPart[14/32] = 7 / 16 - + #> FractionalPart[4/(1 + 3 I)] = 2 / 5 - I / 5 - + #> FractionalPart[Pi^20] = -8769956796 + Pi ^ 20 """ - + attributes = ('Listable', 'NumericFunction', 'ReadProtected') - + def apply(self, n, evaluation): 'FractionalPart[n_]' expr = Expression('FractionalPart', n) return _fractional_part(self.__class__.__name__, n, expr, evaluation) - + def apply_2(self, n, evaluation): 'FractionalPart[n_Complex]' expr = Expression('FractionalPart', n) n_real = Expression("Re", n).evaluate(evaluation) n_image = Expression("Im", n).evaluate(evaluation) - - real_fractional_part = _fractional_part(self.__class__.__name__, n_real, expr, evaluation) + + real_fractional_part = _fractional_part(self.__class__.__name__, n_real, expr, evaluation) image_fractional_part = _fractional_part(self.__class__.__name__, n_image, expr, evaluation) return Expression('Complex', real_fractional_part, image_fractional_part) @@ -634,7 +634,7 @@ class CoprimeQ(Builtin): = True >> CoprimeQ[4+2I, 6+3I] - = False + = True >> CoprimeQ[2, 3, 5] = True @@ -964,4 +964,3 @@ def apply(self, m, n, evaluation): return Expression('List', Integer(py_m // py_n), (py_m % py_n)) else: return Expression('QuotientRemainder', m, n) - diff --git a/mathics/builtin/structure.py b/mathics/builtin/structure.py index ad56b7c98c..cf4a501966 100644 --- a/mathics/builtin/structure.py +++ b/mathics/builtin/structure.py @@ -7,18 +7,26 @@ BinaryOperator, Test, MessageException, + PartRangeError, ) from mathics.core.expression import ( Expression, String, Symbol, + SymbolTrue, + SymbolFalse, Integer, Rational, strip_context, ) -from mathics.core.rules import Pattern +from mathics.core.rules import Pattern, Rule -from mathics.builtin.lists import python_levelspec, walk_levels, InvalidLevelspecError +from mathics.builtin.lists import ( + python_levelspec, + walk_levels, + InvalidLevelspecError, + List, +) from mathics.builtin.functional import Identity import platform @@ -517,6 +525,54 @@ def callback(level): return result +class MapAt(Builtin): + """ +
+
'MapAt[$f$, $expr$, $n$]' +
applies $f$ to the element at position $n$ in $expr$. If $n$ is negative, the position is counted from the end. +
+ + Map $f$ onto the part at position 2: + >> MapAt[f, {a, b, c, d}, 2] + = {a, f[b], c, d} + + Map $f$ onto the at the end: + >> MapAt[f, {a, b, c, d}, -1] + = {a, b, c, f[d]} + + Map $f$ onto an association: + >> MapAt[f, <|"a" -> 1, "b" -> 2, "c" -> 3, "d" -> 4, "e" -> 5|>, 3] + = {a -> 1, b -> 2, c -> f[3], d -> 4, e -> 5} + + Use negative position in an association: + >> MapAt[f, <|"a" -> 1, "b" -> 2, "c" -> 3, "d" -> 4|>, -3] + = {a -> 1, b -> f[2], c -> 3, d -> 4} + """ + + def apply(self, f, expr, n, evaluation, options={}): + "MapAt[f_, expr_, n_Integer]" + i = n.get_int_value() + m = len(expr.leaves) + if 1 <= i <= m: + j = i - 1 + elif -m <= i <= -1: + j = m + i + else: + raise PartRangeError + + new_leaves = list(expr.leaves) + replace_leaf = new_leaves[j] + if hasattr(replace_leaf, "head") and replace_leaf.head == Symbol("System`Rule"): + new_leaves[j] = Expression( + "System`Rule", + replace_leaf.leaves[0], + Expression(f, replace_leaf.leaves[1]), + ) + else: + new_leaves[j] = Expression(f, replace_leaf) + return List(*new_leaves) + + class Scan(Builtin): """
@@ -825,9 +881,9 @@ def apply(self, expr, form, evaluation): form = Pattern.create(form) if expr.is_free(form, evaluation): - return Symbol("True") + return SymbolTrue else: - return Symbol("False") + return SymbolFalse class Flatten(Builtin): @@ -1031,6 +1087,7 @@ def insert_leaf(leaves): new_leaves.append(Expression(h, *insert_leaf(group))) return new_leaves + return Expression(h, *insert_leaf(leaves)) def apply(self, expr, n, h, evaluation): diff --git a/setup.py b/setup.py index 70fda3ec9f..18efc87784 100644 --- a/setup.py +++ b/setup.py @@ -80,7 +80,7 @@ def read(*rnames): # General Requirements INSTALL_REQUIRES += [ - "sympy>=1.6, < 1.7", + "sympy>=1.7, <= 1.8dev", "django >= 3.0, < 3.2", "mpmath>=1.1.0", "numpy", diff --git a/test/test_combinatorica.py b/test/test_combinatorica.py index a25fea1786..6d40c085af 100644 --- a/test/test_combinatorica.py +++ b/test/test_combinatorica.py @@ -13,14 +13,12 @@ """ ) -# A number of examples from Implementing Discrete Mathematics by -# Steven Skiena and -# A number of examples from Computation Discrete Mathematics by -# Sriram Pemmaraju and Steven Skiena. +# A number of examples from: +# * Implementing Discrete Mathematics by Steven Skiena and +# * Computation Discrete Mathematics by Sriram Pemmaraju and Steven Skiena. # Page numbers below come from the first book - def test_combinatorica_permutations_1_1(): for str_expr, str_expected, message in ( @@ -216,8 +214,63 @@ def test_combinatorica_permutations_1_2(): ( "Apply[Or, Map[(# === HideCycles[ToCycles[#]])&, Permutations[Range[5]] ]]", "False", - "None of the permutatoins on five elements is identical to its hidden cycle representation 1.2.4, Page 23", + "None of the permutations on five elements is identical to its hidden cycle representation 1.2.4, Page 23", + ), + ( + "{StirlingFirst[6,3], StirlingS1[6,3]}", + "{225, -225}", + "StirlingFirst 1.2.4, Page 24", + ), + ( + "Select[ Map[ToCycles, Permutations[Range[4]]], (Length[#]==2)&]", + "{{{1}, {3, 4, 2}}, {{1}, {4, 3, 2}}, {{2, 1}, {4, 3}}, " + " {{2, 3, 1}, {4}}, {{2, 4, 1}, {3}}, {{3, 2, 1}, {4}}, " + " {{3, 4, 1}, {2}}, {{3, 1}, {4, 2}}, {{4, 2, 1}, {3}}, " + " {{4, 3, 1}, {2}}, {{4, 1}, {3, 2}}}", + "11 permutations of 4 elements and 2 cycles, Page 24", + ), + ( + "NumberOfPermutationsByCycles[4,2]", + "11", + "NumberOfPermutationsByCycles 1.2.4, Page 24", + ), + ( + "{StirlingSecond[6,3], StirlingS2[6,3]}", + "{90, 90}", + "StirlingSecond 1.2.4, Page 24", + ), + ( + "SignaturePermutation[{1,3,2,4,5,6,7,8}]", + "-1", + "SignaturePermutation 1.2.5, Page 25", + ), + ( + "SignaturePermutation[p] == SignaturePermutation[InversePermutation[p]]", + "True", + "A particular permutation has the same sign as its inverse 1.2.5, Page 25", + ), + ( + "PermutationGroupQ[ Select [ Permutations[Range[4]], (SignaturePermutation[#]==1)&] ]", + "True", + "All permutations have the same sign as their inverse 1.2.5, Page 25", + ), + ( + "Polya[Table[RotateRight[Range[8],i], {i, 8}], m]", + "(4 m + 2 m ^ 2 + m ^ 4 + m ^ 8) / 8", + "Polya counting resulting in polynomial 1.2.6, Page 25", + ), + # MapAt not finished which is probably causing Poly to fail here... + # ( + # "Polya[Automorphisms[Cycle[8]], m]", + # "(4 m + 2 m ^ 2 + m ^ 4 + m ^ 8) / 8", + # "Polya counting resulting in polynomial 1.2.6, Page 26", + # ), + ( + "Factor[(4 m + 2 m^2 + 5m^4 + 4m^5 + m^8)/16]", + "m (1 + m) (4 - 2 m + 2 m ^ 2 + 3 m ^ 3 + m ^ 4 - m ^ 5 + m ^ 6) / 16", + "Factor example in Polya polynomial 1.2.6, Page 26", ), + ): check_evaluation(str_expr, str_expected, message)