Inverse Harmonic and Inverse Digamma Functions

[Albert Chan]

Compare these expressions, which should be equivalent:

XCas> subst(x>0, x=1)                              → true
XCas> subst(ifte(x>0, true, false), x=1)      → "Ifte: Unable to check test Error: Bad Argument Value"

ifte() tried to evaluate the expression, before substitution. However, "x>0", cannot be evaluated.

XCas> subst(ifte(x-x, true, false), x=1)       → false
XCas> subst(0/x, x=0)                              → 0

Again, expressions were evaluated first, before substitution.

This is branchless invPsi(x), to avoid above ifte() issues.
This version is able to do "plot(invPsi(x), x=-1 .. 1)"

Code:

invPsi(x) := {
 local g, f0, f1;
 g := [x >= -1.8, x < -1.8] * [invPsi_guess(e^x), -1/(x+not(x))];
 g -= (f0 := Psi(g)-x) * (f1 := Psi(g,1)) / (f1*f1 - 0.5*f0*Psi(g,2));
 g -= (f0 := Psi(g)-x) * (f1 := Psi(g,1)) / (f1*f1 - 0.5*f0*Psi(g,2));
};

Both "then" and "else" branch were evaluated, which is not ideal.
Also, -1/x guess replaced by -1/(x+not(x)), to patch against divide-by-zero.

XCas> subst(invPsi(x), x = 0)                          → 1.46163214497
XCas> subst(invPsi(x), x = -euler_gamma)       → 1.0