Nabble has removed Mailing-list integration.
Posts created here DO NOT GET SENT TO THE MAILING LIST.
Mailing-list emails DO NOT GET POSTED TO THE FORUM.
So basically the Forum is now out of date, we are looking into migrating the history.
For now you should send emails, people will see them, discuss@lists.openscad.org.
Recursive Binary Chop
Locked
12
12
 
|
I had a complicated function to calculate an area according to its one argument. The function was monotonic.
There was then a need to find the argument value which gave a specific area - an inverse function was created. Its arguments are the input value to be looked up lowest start value and its function highest start value and its function required accuracy InvFF performs a binary chop until the upper and lower values are within the required accuracy //---------------------------------------------------------------------------------------- // test out with a simple inverse sin // sin(T) goes between 0 and 1 for angles 0 : 90 // if we give the recursive function a number 0 : 1, we want the angle 0 : 90 //------------------------------------------------------------------------------- function FF(A) = sin(A); //------------------------------------------------------------------------------- function InvFF(Value,LowArg,LowValue,HighArg,HighValue,Accuracy) = let(MidArg = (HighArg + LowArg)/2) abs(HighValue - LowValue)/(HighValue + LowValue) < 2*Accuracy ? MidArg // Stop here with best answer : FF(MidArg) > Value ? InvFF(Value,LowArg,LowValue,MidArg,FF(MidArg),Accuracy) // too high : InvFF(Value,MidArg,FF(MidArg),HighArg,HighValue,Accuracy); // too low //------------------------------------------------------------------------------ LowestValue = 0; // for this test HighestValue = 90; // for this test // lookup the angle which gives a sine of 0.5 - ie 30 degrees echo(InvFF(0.5, LowestValue, FF(LowestValue), HighestValue, FF(HighestValue), 0.00001)); returns ECHO: 29.9999 Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
Nice! On Sat, May 15, 2021 at 12:10 PM GregL <[hidden email]> wrote: I had a complicated function to calculate an area according to its one argument. The function was monotonic. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
Yes, nice, I have written one or two binary chops to do things like space two pulleys for a closed belt, which surprisingly has no algebraic formula, like the circumference of an ellipse. But it was before function literals, so not generic. On Sat, 15 May 2021 at 20:36, Leonard Martin Struttmann <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
CONTENTS DELETED
The author has deleted this message.
|
|
|
You start with two extreme values that you assume to be either side of the solution. Then you pick the midpoint and decide if the solution is between the first point and the mid or the mid and the second, so you get a smaller range for the solution and then recurse until the ends points are close enough for your desired accuracy. On Sat, 15 May 2021 at 21:51, craig and heather <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
CONTENTS DELETED
The author has deleted this message.
|
|
|
I think the standard name for this algorithm in the context of root finding is Bisection. I've never heard "binary chop" before and "binary search" to me refers to searching a sorted list, not root finding.
https://en.wikipedia.org/wiki/Bisection_method The spacing for two pulleys seems like a problem that does have a closed form solution. I worked out math to make models like the one below with varying radii at each end and separation and it was all simple trig. 
Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
In reply to this post by GregL
Binary chop can work for this but Newton Raphson should do it in less iterations https://en.wikipedia.org/wiki/Newton%27s_method On Sun, May 16, 2021 at 1:09 AM GregL <[hidden email]> wrote: I had a complicated function to calculate an area according to its one argument. The function was monotonic. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
But you need the derivative of the function for Newton Raphson. On Mon, 17 May 2021 at 14:54, Colin Hercus <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
You don't. You work out the slope from your two initial values. On Mon, May 17, 2021 at 10:27 PM nop head <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
function FF(A) = sin(A); //------------------------------------------------------------------------------- function InvFF(Value,LowArg,LowValue,HighArg,HighValue,Accuracy) = let(MidArg=LowArg + (Value-LowValue)*(HighArg-LowArg)/(HighValue-LowValue)) abs(HighValue - LowValue)/(HighValue + LowValue) < 2*Accuracy ? MidArg // Stop here with best answer: FF(MidArg) > Value ? InvFF(Value,LowArg,LowValue,MidArg,FF(MidArg),Accuracy) // too high : InvFF(Value,MidArg,FF(MidArg),HighArg,HighValue,Accuracy); // too low //------------------------------------------------------------------------------ LowestValue = 0; // for this test HighestValue = 90; // for this test // lookup the angle which gives a sine of 0.5 - ie 30 degrees //echo(LowestValue, FF(LowestValue), HighestValue, FF(HighestValue)); echo(InvFF(0.5, LowestValue, FF(LowestValue), HighestValue, FF(HighestValue), 0.00001)); ECHO: 30 Rendering Polygon Mesh using CGAL... UI-WARNING: No top level geometry to render On Mon, May 17, 2021 at 10:33 PM Colin Hercus <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
And after reworking recursion and termination //---------------------------------------------------------------------------------------- // test out with a simple inverse sin // sin(T) goes between 0 and 1 for angles 0 : 90 // if we give the recursive function a number 0 : 1, we want the angle 0 : 90 //------------------------------------------------------------------------------- function FF(A) = sin(A); //------------------------------------------------------------------------------- function InvFF(Value,LowArg,LowValue,HighArg,HighValue,Accuracy) = let(MidArg=LowArg + (Value-LowValue)*(HighArg-LowArg)/(HighValue-LowValue)) abs(HighValue - LowValue)/(HighValue + LowValue) < 2*Accuracy ? MidArg // Stop here with best answer : InvFF(Value,HighArg,HighValue,MidArg,FF(MidArg),Accuracy); //------------------------------------------------------------------------------ LowestValue = 0; // for this test HighestValue = 90; // for this test // lookup the angle which gives a sine of 0.5 - ie 30 degrees //echo(LowestValue, FF(LowestValue), HighestValue, FF(HighestValue)); echo(InvFF(0.5, LowestValue, FF(LowestValue), HighestValue, FF(HighestValue), 0.00001)); On Mon, May 17, 2021 at 11:00 PM Colin Hercus <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
Sorry, I can't help myself. Even simpler //---------------------------------------------------------------------------------------- // test out with a simple inverse sin // sin(T) goes between 0 and 1 for angles 0 : 90 // if we give the recursive function a number 0 : 1, we want the angle 0 : 90 //------------------------------------------------------------------------------- function FF(A) = sin(A); //------------------------------------------------------------------------------- function InvFF(Value,LowArg,HighArg) = let(MidArg=LowArg + (Value-FF(LowArg))*(HighArg-LowArg)/(FF(HighArg)-FF(LowArg))) MidArg==LowArg || MidArg==HighArg ? MidArg : InvFF(Value,HighArg,MidArg,Accuracy); //------------------------------------------------------------------------------ LowestValue = 0; // for this test HighestValue = 90; // for this test // lookup the angle which gives a sine of 0.5 - ie 30 degrees echo(InvFF(0.5, LowestValue, HighestValue)); On Mon, May 17, 2021 at 11:17 PM Colin Hercus <[hidden email]> wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
Ummmm... "...even simpler.." is broke.
What happened to Accuracy parameter? Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
Accuracy isn't req'd but I forgot to remove it from the recursive call of InvFF. It still worked and I missed the warning message. Sorry, it was very late here. There is now no accuracy setting and it still does less iterations than the Binary Chop. Put and echo(MidArg) after the let() to watch it in action. function InvFF(Value,LowArg,HighArg) = let(MidArg=LowArg + (Value-FF(LowArg))*(HighArg-LowArg)/(FF(HighArg)-FF(LowArg))) MidArg==LowArg || MidArg==HighArg ? MidArg : InvFF(Value,HighArg,MidArg); On Tue, May 18, 2021 at 3:06 AM LenStruttmann <[hidden email]> wrote: Ummmm... "...even simpler.." is broke. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
Much better! Thanks!
Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
In reply to this post by Colin49
Note that this isn't really Newton's method, which requires use of the actual derivative, or even a standard approximate Newton's method, where the derivative would be approximated numerically, e.g. by (f(x+h)-f(x-h))/2h. Anybody planning to use it should be aware that the solution can be any real number and is by no means guaranteed to lie in the interval [LowArg,HighArg]. Is there a convergence proof for this algorithm? Can you prove that the interval shrinks? I have this feeling that it might be possible to construct some kind of infinite loop case where it fails to converge because the derivative approximation is so horrible. It will blow up with an error if the value of the function is equal at both endpoints, and it seems likely to misbehave if the function values are almost equal at the endpoints.
Given these various issues, one might ask, why not just use a regular approximate Newton's method, which will be better behaved, and easier to understand, and wouldn't require the mysterious use of two parameters (the high and low values) to kick off the iteration. If you invoke this with the "high" and "low" values different by only a small amount like 1e-5 then I think might get an actual approximate Newton method. (I'm haven't tried to prove that the interval doesn't grow, which would be necessary to ensure that the derivative approximation stays good.)
Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
On Tue, May 18, 2021 at 10:09 AM adrianv <[hidden email]> wrote: Note that this isn't really Newton's method, which requires use of the actual derivative, or even a standard approximate Newton's method, where the derivative would be approximated numerically, e.g. by (f(x+h)-f(x-h))/2h. Anybody planning to use it should be aware that the solution can be any real number and is by no means guaranteed to lie in the interval [LowArg,HighArg]. Is there a convergence proof for this algorithm? Can you prove that the interval shrinks? I have this feeling that it might be possible to construct some kind of infinite loop case where it fails to converge because the derivative approximation is so horrible. It will blow up with an error if the value of the function is equal at both endpoints, and it seems likely to misbehave if the function values are almost equal at the endpoints. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
In reply to this post by adrianv
Hi, In the original email it was stated that it was used for a monotonic function. The same restriction applies to this. The derivative is (HighArg-LowArg)/(FF(HighArg)-FF(LowArg)) or the inverse of this. It's an approximation for the Hi/Lo range and it's not terrible. It should work fine for any non constant monotonic function where the values of FF at High & Low are not infinite. The binary chop/search will not work if the function is not monotonic and always returns HighArg for a constant function. For a constant function the results are undefined. It should fault with a divide by zero. It would be easy to add a check for this. I can see there is a possible issue if the function has a region of constant values (derivative = 0) and we hit that with our Hi&Lo points. In the original algorithm the user had to choose a High/Low range that was meant to contain the solution. I also think that it works better starting with a large range but my code doesn't require that. With functions like Sine that are not monotonic over their entire range it's important that the hi and low cover a monotonic region and that the target value is inside the region. Much the same as the original binary chop algorithm. I can see there is a possible issue if the function has a region of constant values (derivative = 0) Just adding a test for HiValue == LowValue and returning (HighArg+LowArg)/2 would solve this. But again the assumption was monotonically increasing function. On Tue, May 18, 2021 at 10:09 AM adrianv <[hidden email]> wrote: Note that this isn't really Newton's method, which requires use of the actual derivative, or even a standard approximate Newton's method, where the derivative would be approximated numerically, e.g. by (f(x+h)-f(x-h))/2h. Anybody planning to use it should be aware that the solution can be any real number and is by no means guaranteed to lie in the interval [LowArg,HighArg]. Is there a convergence proof for this algorithm? Can you prove that the interval shrinks? I have this feeling that it might be possible to construct some kind of infinite loop case where it fails to converge because the derivative approximation is so horrible. It will blow up with an error if the value of the function is equal at both endpoints, and it seems likely to misbehave if the function values are almost equal at the endpoints. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
|
|
I see that the original poster said the input should be monotonic, but the bisection method does not require a monotonic input. It requires that the function be continuous and that the function values at the interval endpoints have opposite sign. A properly written implementation would check for this and fail immediately if this condition doesn't hold.
The secant slope is not a good approximation to the derivative where by "good" I mean that we can bound its error. If you want the derivative at zero and use the interval [0,1] and consider a function like ln(x+c) then the secant approximation is ln(1+c)-ln(c) and the actual derivative is 1/c. This approximation is arbitrarily bad, in the limit as c->0. If c = 0.1 then the derivative is 10 and the secant approximation is 2.4. If c = 0.01 the derivative is 100 and the secant approximation is 4.6. What you have implemented is actually the Secant Method: https://en.wikipedia.org/wiki/Secant_method and the referenced page notes that this method (like Newton's method) is not guaranteed to converge. This is in contrast to the bisection method, which *is* guaranteed to converge, albeit with only linear convergence. The rate of convergence of the secant method when it does converge (to a simple root of a C2 function) is superlinear, but not quadratic like Newton's method. Unlike the original method, the solution need not lie in the given interval. You're really initializing this method with two guesses for the root, not an interval that contains the root. So where is the function supposed to be monotonic? Over all the real numbers? If you try to solve for ln(x)=-2.3 (monotonic) on the interval [.001,1] it goes into an infinite loop with nans because the iteration takes it outside the valid bound for the log function. If you ask for exp(x)=0.5 and start it on [0,1] it correctly finds -0.69 as the solution---not in the starting interval. If you give it atan you can get it to diverge by giving it an overly large interval, even though atan is monotonic. The requirement of monotonicity doesn't suffice to force good behavior. If you try to solve sin(x)=0 and give it the interval [89,90] where sin is monotonic (but has no solutions) you get the solution -133920 (which is correct, being 372*360). I suppose you could keep track of the interval and return failure if the final solution is not in the interval. Is it possible that the algorithm finds a root outside the starting interval even though one exists there? (I think the answer is yes for the non-monotonic case.)
Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
«
Return to OpenSCAD
|
1 view|%1 views
| Free forum by Nabble | Edit this page |