Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,9 @@
.vscode
/.cache
/.gdbinit
/.python-version
/.pyodide-xbuildenv-*
/.python-version
/ChangeLog.spell-corrected
/Mathics.egg-info
/Mathics3.egg-info
ChangeLog
Expand Down
101 changes: 62 additions & 39 deletions mathics/builtin/arithmetic.py
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@
IterationFunction,
MPMathFunction,
Predefined,
PrefixOperator,
SympyFunction,
SympyObject,
Test,
Expand Down Expand Up @@ -412,7 +413,7 @@ class Conjugate(MPMathFunction):
rules = {
"Conjugate[Undefined]": "Undefined",
}
summary_text = "complex conjugation"
summary_text = "compute complex conjugation"


class DirectedInfinity(SympyFunction):
Expand Down Expand Up @@ -694,9 +695,11 @@ class Integer_(Builtin):
name = "Integer"


class Product(IterationFunction, SympyFunction):
class Product(IterationFunction, SympyFunction, PrefixOperator):
"""
<url>:WMA link:https://reference.wolfram.com/language/ref/Product.html</url>
<url>:Direct product:https://en.wikipedia.org/wiki/Direct_product</url> (<url>
:SymPy:https://docs.sympy.org/latest/modules/concrete.html#sympy.concrete.products.Product</url>, <url>
:WMA:https://reference.wolfram.com/language/ref/Product.html</url>)

<dl>
<dt>'Product'[$f$, {$i$, $i_{min}$, $i_{max}$}]
Expand All @@ -713,34 +716,49 @@ class Product(IterationFunction, SympyFunction):
outermost-to-innermost order.
</dl>

The product of the first 10 integers is its factorial:
>> Product[k, {k, 1, 10}] == 10!
= True
'Product[k, {k, i, n}]' is defined in terms of <url>
:Factorial:
/doc/reference-of-built-in-symbols/special-functions/gamma-and-related-functions/factorial/</url>:

'Product' can be used symbolically:
>> Product[2 ^ i, {i, 1, n}]
= 2 ^ (n / 2 + n ^ 2 / 2)
>> Product[k, {k, i, n}]
= n! / (-1 + i)!

>> Product[f[i], {i, 1, 7}]
= f[1] f[2] f[3] f[4] f[5] f[6] f[7]
When $i$ is 1, we get the factorial function:
>> Product[k, {k, 1, n}]
= n!

Or more succinctly:
>> Product[k, {k, n}]
= n!

Symbolic products involving the factorial are evaluated:
>> Product[k, {k, 3, n}]
= n! / 2

Evaluate the $n$-th primorial:
>> primorial[0] = 1;
>> primorial[n_Integer] := Product[Prime[k], {k, 1, n}];
>> primorial[12]
= 7420738134810
Examples of numeric evaluation using more complex functions:
>> Product[x^k, {k, 2, 20, 2}]
= x ^ 110

"""
>> Product[2 ^ i, {i, 1, n}]
= 2 ^ (n / 2 + n ^ 2 / 2)

summary_text = "discrete product"
throw_iterb = False
>> Product[f[i], {i, 1, 7}]
= f[1] f[2] f[3] f[4] f[5] f[6] f[7]

sympy_name = "Product"
Evaluate the $n$-th <url>:Primorial:https://en.wikipedia.org/wiki/Primorial</url>:
>> Primorial[0] = 1;
>> Primorial[n_Integer] := Product[Prime[k], {k, 1, n}];
>> Primorial[12]
= 7420738134810

"""

# FIXME Product[k, {k, 3, n}] is rewritten using Factorial via
# Pochhammer rewrite rules. We want this for Product, but WMA
# does not rewrite using Factorial for Pochhammer alone, although it could.
# Nevertheless, if and when our Pochhammer is adjusted to remove
# this transformation to Factorial to match WMA behavior,
# we will need to add a rule that transforms to Factorial here.
rules = IterationFunction.rules.copy()
rules.update(
{
Expand All @@ -752,9 +770,12 @@ class Product(IterationFunction, SympyFunction):
),
}
)
summary_text = "compute the direct product"
sympy_name = "Product"
throw_iterb = False

def get_result(self, items):
return Expression(SymbolTimes, *items)
def get_result(self, elements):
return Expression(SymbolTimes, *elements)

def to_sympy(self, expr, **kwargs):
if expr.has_form("Product", 2) and expr.elements[1].has_form("List", 3):
Expand Down Expand Up @@ -895,22 +916,24 @@ def test(self, expr) -> bool:
)


class Sum(IterationFunction, SympyFunction):
"""
<url>:WMA link:https://reference.wolfram.com/language/ref/Sum.html</url>
class Sum(IterationFunction, SympyFunction, PrefixOperator):
r"""
<url>:Summation:https://en.wikipedia.org/wiki/Summation</url> (<url>
:SymPy:https://docs.sympy.org/latest/modules/concrete.html#sympy.concrete.summations.Sum</url>, <url>
:WMA:https://reference.wolfram.com/language/ref/Sum.html</url>)

<dl>
<dt>'Sum'[$expr$, {$i$, $imin$, $imax$}]
<dd>evaluates the discrete sum of $expr$ with $i$ ranging from $imin$ to $imax$.
<dt>'Sum['$f, \{i, i_{min}, i_{max}\}$']'
<dd>evaluates the discrete sum of $f$ with $i$ ranging from $i_{min}$ to $i_{max}$.

<dt>'Sum'[$expr$, {$i$, $imax$}]
<dd>same as 'Sum[$expr$, {$i$, 1, $imax$}]'.
<dt>'Sum['$f, \{i, i_{max}\}$']'
<dd>same as 'Sum['$f, \{i, 1, i_{max}\}$']'.

<dt>'Sum'[$expr$, {$i$, $imin$, $imax$, $di$}]
<dd>$i$ ranges from $imin$ to $imax$ in steps of $di$.
<dt>'Sum['$f, \{i, i_{min}, i_{max}, di\}$']'
<dd>$i$ ranges from $i_{min}$ to $i_{max}$ in steps of $di$.

<dt>'Sum'[$expr$, {$i$, $imin$, $imax$}, {$j$, $jmin$, $jmax$}, ...]
<dd>evaluates $expr$ as a multiple sum, with {$i$, ...}, {$j$, ...}, ... being \
<dt>'Sum['$f, \{i, i_{min}, i_{max}, \{j, j_{min}, j_{max}, ...$']'
<dd>evaluates $f$ as a multiple sum, with {$i, ...$}, {$j, ...$}, ... being \
in outermost-to-innermost order.
</dl>

Expand Down Expand Up @@ -970,12 +993,6 @@ class Sum(IterationFunction, SympyFunction):
= 1 + 2 I
"""

summary_text = "discrete sum"
# Do not throw warning message for symbolic iteration bounds
throw_iterb = False

sympy_name = "Sum"

rules = IterationFunction.rules.copy()
rules.update(
{
Expand All @@ -988,6 +1005,12 @@ class Sum(IterationFunction, SympyFunction):
}
)

summary_text = "compute a summation"
sympy_name = "Sum"

# Do not throw warning message for symbolic iteration bounds
throw_iterb = False

def get_result(self, elements) -> Expression:
return Expression(SymbolPlus, *elements)

Expand Down
12 changes: 10 additions & 2 deletions mathics/builtin/specialfns/gamma.py
Original file line number Diff line number Diff line change
Expand Up @@ -336,7 +336,6 @@ class Gamma(MPMathMultiFunction):

rules = {
"Gamma[z_, x0_, x1_]": "Gamma[z, x0] - Gamma[z, x1]",
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In WMA, this rule only applies if the three arguments are numeric, and at least one is a Real number

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I imagine then that this kind of check would be better done inside an evaluation method.

If you are up for adding, feel free to do. But if not, let's leave for another PR.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OK, let's postpone for another PR.

"Gamma[1 + z_]": "z!",
Copy link
Contributor

@mmatera mmatera Feb 19, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The rule could be
"Gamma[z_Integer]": "(z-1)!"

Copy link
Member Author

@rocky rocky Feb 19, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this cause any observable difference in output? If not, I'd recommend removing the rule unless there is a good reason to add it.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We have different implementations for Gamma and Factorial, relying on different mpmath functions. Right now, mpmath function gives the same answer, so both outputs are the same.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great - then let's not complicate here when there is no user-observable difference.

"Gamma[Undefined]": "Undefined",
"Gamma[x_, Undefined]": "Undefined",
"Gamma[Undefined, y_]": "Undefined",
Expand Down Expand Up @@ -452,7 +451,16 @@ class Pochhammer(SympyFunction):

rules = {
"Pochhammer[0, 0]": "1",
"Pochhammer[a_, n_]": "Gamma[a + n] / Gamma[a]",
# FIXME: In WMA, if n is an Number with an integer value, it
# is rewritten to expanded terms not using
# Factorial as we do below. For example, Pochhammer[i, 2] is
# i (i + 1) instead of (1 + i)! / (-1 + i)! as the rule below
# gives. Ideally, we should match this behavior. However if
# this is done, we will *also* need to adjust Product, because
# WMA Product is sometimes rewritten as an expression using Factorial.
# In particular, Product[k, {k, 3, n}] == n! / 2 in both Mathics3
# and WMA.
"Pochhammer[a_, n_]": "Factorial[a + n - 1] / Factorial[a - 1]",
Copy link
Contributor

@mmatera mmatera Feb 19, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is what Pochhammer should do:

In[1]:= Pochhammer[a,b]                                                         

Out[1]= Pochhammer[a, b]

In[2]:= Pochhammer[a,12.4]                                                      

Out[2]= Pochhammer[a, 12.4]

In[3]:= Pochhammer[a,12.]                                                       

Out[3]= a (1 + a) (2 + a) (3 + a) (4 + a) (5 + a) (6 + a) (7 + a) (8 + a) 
 
>    (9 + a) (10 + a) (11 + a)

In[4]:= Pochhammer[12,b]                                                        

Out[4]= Pochhammer[12, b]

In[5]:= Pochhammer[12.1,b]                                                      

Out[5]= Pochhammer[12.1, b]

Notice that n must be an Integer to apply this rule.

Suggested change
"Pochhammer[a_, n_]": "Factorial[a + n - 1] / Factorial[a - 1]",
"Pochhammer[a_, n_Integer]": "Factorial[a + n - 1] / Factorial[a - 1]",

or maybe n_?NumerQ

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

According to the examples cited, NumberQ only works when the number can be exactly converted to an Integer.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Indeed. In any case, this is the point when we need to decide if the complexity of a rule that matches the exact behavior worth. The exact behavior would be obtained with something like

        "Pochhammer[a_, /;(IntegerQ[n]|| NumberQ[n] && Round[n]==n)]": "Factorial[a + n - 1] / Factorial[a - 1]",

Copy link
Member Author

@rocky rocky Feb 19, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

IntegerQ[n] can be dropped following the principle: simple as possible.... But even if we change to Pochhammer[a_, n_(NumberQ[n] && Round[n]==n)], Product[k, {k, i, n}] will not transform as WMA does, in particular: n! / (-1+i)!. So WMA may be using some other mechanism to implement Product.

But let's back up a bit. We got into this as a result of a real problem reported by @aravindh-krishnamoorthy, in his day-to-day use of Mathics3. His PR broke something else. If this PR addresses both his problem and does not break what is already there, then let's go with this.

Note that this PR adds a test to ensure we do not rewrite Gamma as Factorial in the kind of way that @aravindh-krishnamoorthy was not expecting.

In the future, if someone has a problem as a result of Pochhammer getting expanded to factorial, then we can address that issue at that time while maintaining satisfaction of all the other checks in place.

In other words: complicate when there is a tangible need to complicate. Again this is also encompassed by the principle: Simple as possible, but no simpler.

"Derivative[1,0][Pochhammer]": "(Pochhammer[#1, #2]*(-PolyGamma[0, #1] + PolyGamma[0, #1 + #2]))&",
"Derivative[0,1][Pochhammer]": "(Pochhammer[#1, #2]*PolyGamma[0, #1 + #2])&",
}
Expand Down
1 change: 1 addition & 0 deletions test/builtin/specialfns/test_gamma.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,7 @@
"Overflow",
),
("Gamma[1., 2.]", None, "Gamma[1., 2.]", "needs mpmath for lowergamma"),
("Gamma[1 + x]", None, "Gamma[1 + x]", "Gamma should not expand to factorial"),
],
)
def test_private_doctests_gamma(str_expr, msgs, str_expected, fail_msg):
Expand Down