DM42/Free42 - which algorithm used in SOLVE ?

Discussion around the SwissMicros DM42 calculator
User avatar
OlidaBel
Posts: 58
Joined: Thu Mar 11, 2021 8:52 am
Location: Belgium

DM42/Free42 - which algorithm used in SOLVE ?

Post by OlidaBel »

hello,

My DM42 as well as Free42 encounters difficulties with a SOLVE problem (finding roots of a Bessel function).
This Solve is calling a LBL which integrates a function (for Bessel Jn(x)).
On screen I see the process running... one of the boundary reaching the solution, but for the case I'm interested in, it "never" stops,
for a long time (I couldn't wait to the end, too many minutes). The second boundary doesn't converge very well to the root, why ?
The 2 initial guess values are well around the true root, one giving f(x) positive and one giving f(x) negative, so a run without hurdles...

A bit frustrated, and knowing SOLVE is very fast otherwise, I tried to solve myself with a Regula Falsi method, it converges very fast and reaches the zero after 4 iterations.
https://en.wikipedia.org/wiki/Regula_fa ... on)_method

Before giving some pieces of code, I'd like to know how DM42 is searching and also when does it decide to stop.

A buddy tried my code on his HP-42 and could reach a solution without issues...in less than a minute.
thank you !
---
Olivier
48GX, Prime G2, 50G, DM15L, DM42, 28S, HP 15c CE
Thomas Okken
Posts: 1100
Joined: Tue May 02, 2017 5:48 pm
Location: Netherlands
Contact:

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by Thomas Okken »

SOLVE in Free42 uses the secant method, until it has two guesses that evaluate to values with opposite signs. It then proceeds by using Ridders root refinement, falling back on the secant method, or bisection as a last resort, if the function is not well-behaved.

If you could share instructions how to reproduce your problem, I'll be happy to take a look at it.
User avatar
OlidaBel
Posts: 58
Joined: Thu Mar 11, 2021 8:52 am
Location: Belgium

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by OlidaBel »

Thomas Okken wrote:
Sat Sep 04, 2021 7:53 pm
SOLVE in Free42 uses the secant method, until it has two guesses that evaluate to values with opposite signs. It then proceeds by using Ridders root refinement, falling back on the secant method, or bisection as a last resort, if the function is not well-behaved.

If you could share instructions how to reproduce your problem, I'll be happy to take a look at it.
Hi Thomas.
Many thanks for your answer.
I've attached a raw file extracted from the iPhone.
The goal is to reach a zero for the Bessel function of order 0, between 2.4 and 2.5. The solution is 2.404825557...
The Bessel function integration was also described in the HP 42s manual in EN, p.198.
Here the function J is more general, there is a n.theta term in the equation (n is stored in R00).
BesselJnx.jpg
BesselJnx.jpg (11.85 KiB) Viewed 3593 times
BESSEL.raw
(89 Bytes) Downloaded 168 times
The Bessel.raw contains the function to SOLVE (BESSEL) and the function which is integrated (JNX).
Concerning the running :
Must be in Radian.
R00 must=0
LLIM must=0
ULIM must=Pi
Acc = 0.0001 or 0.0001 (both tested)
SOLVE BESSEL on X with [2.4 ; 2.5]

On iPhone Xs, 30 s are necessary to reach the answer.
On DM42 on USB, I can't see the end after more than 5 minutes.

thanks,
Attachments
keisan.jpg
keisan.jpg (230.96 KiB) Viewed 3593 times
---
Olivier
48GX, Prime G2, 50G, DM15L, DM42, 28S, HP 15c CE
Thomas Okken
Posts: 1100
Joined: Tue May 02, 2017 5:48 pm
Location: Netherlands
Contact:

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by Thomas Okken »

Ah, that looks familiar. The function you're solving is an integral, so in the process of approaching a root, INTEG is having to evaluate an integral that is close to zero. And that is something the Free42 INTEG doesn't handle well, because such integrals appear not to converge, because the estimated absolute error is decreasing at the same rate as the function value, which means the estimated relative error doesn't improve, and so INTEG keeps iterating until it reaches the maximum number of iterations. Since the maximum number of iterations is 20 and the number of function evaluations doubles with each iteration, INTEG may evaluate the function up to a million times, and that can take a few seconds even on a fast computer or smartphone... and on the DM42, it is completely impractical.

This is a known problem. INTEG needs to handle integrals of zero, or close to zero, better. Maybe that means the maximum number of iterations simply needs to be lower, but the real HP-42S seems to break out of such integrals much faster... there may be some cleverness there that Free42 should emulate, but I don't know how its termination criterion works for such cases.
User avatar
OlidaBel
Posts: 58
Joined: Thu Mar 11, 2021 8:52 am
Location: Belgium

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by OlidaBel »

Thomas Okken wrote:
Sun Sep 05, 2021 5:07 am
Ah, that looks familiar. The function you're solving is an integral, so in the process of approaching a root, INTEG is having to evaluate an integral that is close to zero. And that is something the Free42 INTEG doesn't handle well, because such integrals appear not to converge, because the estimated absolute error is decreasing at the same rate as the function value, which means the estimated relative error doesn't improve, and so INTEG keeps iterating until it reaches the maximum number of iterations. Since the maximum number of iterations is 20 and the number of function evaluations doubles with each iteration, INTEG may evaluate the function up to a million times, and that can take a few seconds even on a fast computer or smartphone... and on the DM42, it is completely impractical.

This is a known problem. INTEG needs to handle integrals of zero, or close to zero, better. Maybe that means the maximum number of iterations simply needs to be lower, but the real HP-42S seems to break out of such integrals much faster... there may be some cleverness there that Free42 should emulate, but I don't know how its termination criterion works for such cases.
It's a kind of a challenge :o :geek: .
The DM15L also handle this well. Is Mr William Kahan still reachable ? ;)
---
Olivier
48GX, Prime G2, 50G, DM15L, DM42, 28S, HP 15c CE
Thomas Okken
Posts: 1100
Joined: Tue May 02, 2017 5:48 pm
Location: Netherlands
Contact:

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by Thomas Okken »

Actually, this came up in this forum recently, and I think the case under discussion that time didn't involve the integral and absolute error estimates spiraling towards zero, but rather, non-improvement, i.e. the estimate of the integral was just bouncing around and not improving at all.

Either way, work needs to be done to make the termination condition better. It is a challenge because you don't want to give up too soon, in case the first few estimates are bouncing around but then the iteration starts converging a bit later...
User avatar
OlidaBel
Posts: 58
Joined: Thu Mar 11, 2021 8:52 am
Location: Belgium

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by OlidaBel »

Thomas Okken wrote:
Sun Sep 05, 2021 3:06 pm
Actually, this came up in this forum recently, and I think the case under discussion that time didn't involve the integral and absolute error estimates spiraling towards zero, but rather, non-improvement, i.e. the estimate of the integral was just bouncing around and not improving at all.

Either way, work needs to be done to make the termination condition better. It is a challenge because you don't want to give up too soon, in case the first few estimates are bouncing around but then the iteration starts converging a bit later...
yes.
I was also thinking the chosen ACC could have an influence on the integral value quality, when each time SOLVE has to interpret integral values that could "hesitate" to match the real integral value. I'm not sure now. Very small ACC did not improve the speed.
The Regula falsi method did show very good converging results on this case, I see the Ridders method idea is based on it. After 4 or 5 iterations there is a stable solution, each iteration further shows [2.4 ; 2.4048255577]. I'll implement the Ridder's and see how it behaves on some iterations, and think about this (interesting) issue.
---
Olivier
48GX, Prime G2, 50G, DM15L, DM42, 28S, HP 15c CE
User avatar
OlidaBel
Posts: 58
Joined: Thu Mar 11, 2021 8:52 am
Location: Belgium

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by OlidaBel »

Thomas Okken wrote:
Sun Sep 05, 2021 5:07 am
Ah, that looks familiar. The function you're solving is an integral, so in the process of approaching a root, INTEG is having to evaluate an integral that is close to zero. And that is something the Free42 INTEG doesn't handle well, because such integrals appear not to converge, because the estimated absolute error is decreasing at the same rate as the function value, which means the estimated relative error doesn't improve, and so INTEG keeps iterating until it reaches the maximum number of iterations. Since the maximum number of iterations is 20 and the number of function evaluations doubles with each iteration, INTEG may evaluate the function up to a million times, and that can take a few seconds even on a fast computer or smartphone... and on the DM42, it is completely impractical.

This is a known problem. INTEG needs to handle integrals of zero, or close to zero, better. Maybe that means the maximum number of iterations simply needs to be lower, but the real HP-42S seems to break out of such integrals much faster... there may be some cleverness there that Free42 should emulate, but I don't know how its termination criterion works for such cases.
Hello 8-)
Could you define the relative error here, and why those 4 iterations (Ridder's) don't lead to a success ?

I tested again on the same integral equation, with ACC=0.00001 .
I began with [2.4 ; 2.5] and followed the algorithm described here :
https://en.wikipedia.org/wiki/Ridders%27_method
Ridders_4_iterations
Ridders_4_iterations
ridders.JPG (59.8 KiB) Viewed 3421 times
Later observations : successive x3's are converging and the integral goes to something very small. But later also with further iterations, no real progress is done on the x3 and the integral gets smaller and smaller (and takes more and more time), well beyond than the required precision maybe. Is something like 1E-16 not sufficient to stop SOLVE ? Just wondering and trying to understand.
---
Olivier
48GX, Prime G2, 50G, DM15L, DM42, 28S, HP 15c CE
Thomas Okken
Posts: 1100
Joined: Tue May 02, 2017 5:48 pm
Location: Netherlands
Contact:

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by Thomas Okken »

The problem isn't in SOLVE. With ACC=1e-4, the root is found in 22 evaluations of BESSEL. The problem is that INTEG takes a long time: JNX is evaluated 4,218,730 times before the root is found.
User avatar
OlidaBel
Posts: 58
Joined: Thu Mar 11, 2021 8:52 am
Location: Belgium

Re: DM42/Free42 - which algorithm used in SOLVE ?

Post by OlidaBel »

Thomas Okken wrote:
Mon Sep 06, 2021 5:53 pm
The problem isn't in SOLVE. With ACC=1e-4, the root is found in 22 evaluations of BESSEL. The problem is that INTEG takes a long time: JNX is evaluated 4,218,730 times before the root is found.
OK Thomas,
I understand, and confirm this observation.
---
Olivier
48GX, Prime G2, 50G, DM15L, DM42, 28S, HP 15c CE
Post Reply