Add add FIXMEs to cover of recent discussion

In the future, Pochhamer may be adjusted to match WMA behavior
better. If and when that happens, Product may need to be adjusted to use
Factorial.

mathics/builtin/arithmetic.py

@@ -753,6 +753,12 @@ class Product(IterationFunction, SympyFunction, PrefixOperator):
 
     """
 
+    # 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(
         {

mathics/builtin/specialfns/gamma.py

@@ -451,6 +451,15 @@ class Pochhammer(SympyFunction):
 
     rules = {
         "Pochhammer[0, 0]": "1",
+        # 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]",
         "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])&",