Inverse Harmonic and Inverse Digamma Functions
09-01-2020, 09:06 PM (This post was last modified: 10-30-2020 10:09 PM by Albert Chan.)
Post: #1
 Albert Chan Senior Member Posts: 2,351 Joined: Jul 2018
Inverse Harmonic and Inverse Digamma Functions
I don't actually have HP Prime, below code tested with XCas

invPsi_guess(k) based on inverse harmonic function f3, optimized for small k = e^x
For small x, Ψ(x) ≈ -1/x. Flip the order, we have x ≈ -1/Ψ(x)

Guess is improved by Halley's method, 2 times

Code:
invPsi_guess(k) := k - k/(24.*k*k+1.) + 0.5;  // k = e^x

invPsi(x) := {
local g, f0, f1;
g := when(x < -1.8, -1./x, invPsi_guess(e^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));
};

invH(h) := invPsi(h-euler_gamma) - 1;

Test functions for round-trip ...

XCas> map(x -> Psi(invPsi(x)), range(-5,5))
﻿ ﻿ ﻿ ﻿ ﻿ ﻿ → [-5.0, -4.0, -3.0, -2.0, -1.0, 3.19839640883e-17, 1.0, 2.0, 3.0, 4.0]

XCas> map(k -> invPsi(Psi(10.^k)), range(-5,5))
﻿ ﻿ ﻿ ﻿ ﻿ ﻿ → [1e-05, 0.0001, 0.001, 0.01, 0.1, 1.0, 10.0, 100.0, 1000.0, 10000.0]

XCas> map(n -> invH(harmonic(n)), range(1,10))
﻿ ﻿ ﻿ ﻿ ﻿ ﻿ → [1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0]

Update: replaced ifte by when, to avoided the ifte bug.
 « Next Oldest | Next Newest »

 Messages In This Thread Inverse Harmonic and Inverse Digamma Functions - Albert Chan - 09-01-2020 09:06 PM RE: Inverse Harmonic and Inverse Digamma Functions - Albert Chan - 09-02-2020, 11:03 AM RE: Inverse Harmonic and Inverse Digamma Functions - Albert Chan - 10-30-2020, 09:57 PM RE: Inverse Harmonic and Inverse Digamma Functions - Albert Chan - 09-02-2020, 07:48 PM RE: Inverse Harmonic and Inverse Digamma Functions - Albert Chan - 09-03-2020, 02:05 PM

User(s) browsing this thread: 1 Guest(s)