[Rafael-SOWNet]

Part of #171.

The polygamma functions are nodes. polygamma(n, z) is psi^(n)(z), the (n + 1)-th derivative of ln Gamma(z). digamma(z) is polygamma(0, z), and prints that way.

• Parsing and printing. Both names are parsed, and printed in text, LaTeX (\psi, \psi^{(n)}) and SymPy (digamma, polygamma). Before, digamma(1) was read as the variable digamma times 1, and polygamma(1, 1) did not parse.
• Derivatives. Each order differentiates to the next. An order that depends on the variable declines.
• Domain. The condition is a whole order and an argument other than 0, -1, -2, ....
• The factorial's derivative. It is x! * digamma(x + 1), where it declined (Differentiating a factorial or gamma answers NaN, which claims the derivative does not exist #958). So gamma(x) differentiates to (x - 1)! * digamma(x).

Evaluation is to the working precision, for every complex argument but the poles, where it is NaN.

• Left of Re z = 1/2. Reflection (DLMF 5.15.6), with the derivatives of cot(pi z) as polynomials in cot.
• Right of it. The recurrence (5.15.5) moves the argument to Re w >= R, with R the working digits plus n. Then the asymptotic series (5.11.2, 5.15.8) is summed to the working precision.
• Bernoulli numbers. Exact, from the tangent numbers by Brent and Harvey's integer recurrence, and cached.

Results:

Precision	Accuracy against mpmath	Cost of a fresh evaluation
100 digits	to the last digit, at the test file's points	about 2–9 ms
300 digits	to the last digit, at 2100-digit references	about 0.1 s
2000 digits	to the last digit, at 2100-digit references	2–6 s

The 300- and 2000-digit points are digamma(3/2), polygamma(3, 7/3) and digamma(-7/3 + 2i).

Not in this PR.

• Special values. They stay as written, since digamma(1) = -gamma wants a symbol for Euler's constant, which the library does not have.
• Compiling and the interval evaluator. Neither takes the node yet.
• Integration rules. int psi(x) dx = ln Gamma(x) and the rest are not here.
• The inverse half of the issue. x! = y has no closed-form inverse, and an equation in either function is left as the set of x for which it holds.

PolygammaTest covers:

• values on the real line and off it, against mpmath at 130 digits: reflection, large arguments, orders 0 to 5;
• the poles;
• the factorial's derivative at 3, which is 7.53670601059080…;
• each order differentiating to the next;
• parsing and printing round trips.

Results:

• Unit suite: passes (16,767).
• Corpus: the full Rubi corpus, run against the 2.6.0 run on the same harness (c614b2c8), is 48,305 of 48,305, with no row moved and no answer's allocation grown.
• Parser: regenerated, with its classes internal.
• PublicApi.txt: gains only the node, the kernel and the two factories.

🤖 Generated with Claude Code

https://claude.ai/code/session_012sonx8iAspMiwRwokT1Ura