When you described it as "horrible", I expected something more along the lines of "a bit awkward", but horrible does not give that solution justice. This goes beyond horrible, and into pure terror... 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 nophead
I took the horror from Wolfram and collected all the common terms and constants. This is my result function fx(a,b,w)= let ( c1=12*2^(1/3), c2=6*2^(2/3), t1 = 1728*w^2*b^10 + 25920*w^4*b^8 + 5184*a^2*w^2*b^8 + 82944*w^6*b^6 - 134784*a^2*w^4*b^6 + 5184*a^4*w^2*b^6 - 110592*w^8*b^4 + 82944*a^2*w^6*b^4 + 25920*a^4*w^4*b^4 + 1728*a^6*w^2*b^4, t2 = 2*a^6 + 6*b^2*a^4 + 48*w^2*a^4 + 6*b^4*a^2 + 384*w^4*a^2 - 120*b^2*w^2*a^2 + 2*b^6 + 1024*w^6 - 1344*b^2*w^4 + 264*b^4*w^2, t3=sqrt3(t2+sqrt(t1)), t4=(a^4+2*b^2*a^2+16*w^2*a^2+b^4+64*w^4-56*b^2*w^2), t5=(a^2+b^2-4*w^2), t6=sqrt( a^2/4 - t5/6 + t3/c1 + t4 / (t3*c2) ), t7=sqrt( a^2/2 - t5/3 - t3/c1 - t4 / (t3*c2) - (a^3-4*w^2*a-t5*a) / (4*t6) ) ) a/4 - t6/2 + t7/2; Still looks odd, but it gives me all I want, The only remaining issues are numeric overflows, when 'w' gets near 'a' ¹). But I can live with this, because in my use case 'w' is almost always something about 0.1*a … 0.3*a. Thank you all for your support, regards. ¹) I assume, this is because OpenSCAD is doing arithmetic in float. Using the same formula in lua with double arithmetic let me come much closer to 'a' without numerical overflows PS: I had to invent function sqrt3(x)=sign(x)*abs(x)^(1/3);because OpenSCAD coult not evaluate cubic root of negative values ( -8^(1/3) := -2) Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
As adrianv already stated, finding and figuring out the algebraic solution is much more difficult than finding the solution in an iterative fashion.
A short analysis of the problem shows that the points [x,0] and [0, y] need to have the same distance from the middle point M=[a/2,b/2] and w from each other. This means we can set up two equations and have two unknowns. However it is not clear, whether there exists some algebraic solution in a closed form, which has already been discussed. Therefore, we better employ some iteration scheme. We choose some interval [x1 ... x2] that surely contains the solution x and test the middle point (x1+x2)/2. Depending on the sign of the difference delta of the distances from M, we recursively continue with the upper respectivley lower half of the interval, until delta < eps holds. For testing purposes, I added some max iteration depth test, which shows that the code needs 26 iteration steps for eps=1E-7 a = 13; b = 15; w = 9.16399; echo(solveit()); function solveit(M=[a/2,b/2],w=w, x1=a, x2=0, eps=1E-7, n = 100) = let(x = (x1+x2)/2) // equation 1 let(y = sqrt(w*w-x*x)) // equation 2 let(delta = norm(M-[x,0])- norm(M-[0,y])) // test condition abs(delta) < eps || n<=0 ? [x, y, n]: delta>0? solveit(M, w, x1, x, n=n-1): solveit(M, w, x, x2, n=n-1); 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 bassklampfe
I did a little googling and found this discussion: https://www.physicsforums.com/threads/rectangle-inscribed-in-another-rectangle-solutions-for-all-cases.508715/
The first attachment is a pretty rigorous write up (https://www.physicsforums.com/attachments/inscribedrectangles-discussion-txt.36642/). Case 3 is the one analogous to the given problem. For his analytical solution, he derives nearly the same quartic as NateTG came up with (except instead of "2 w a x", it's "2 w^2 a x"), which turns into the same result as the mess from Wolfram (and now that I look at it again, the same as adrianv used for the BOSL root finder). He also details an algorithm for determining the inscribed rectangle angle by finding the inflection point of a function of the angle. I thought it was an interesting explanation and covered all of the points raised in this thread. Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
I think the Wolfram mess is only because it produces the overall formula in one step and with some constant folding. Quartics are solved by making various substitutions of variables to simplify the problem and if you code it the same way the expressions are a lot simpler. Quadratics can be solved by completing the square and when you fully expand it you get the familiar formula that is not too complex but as soon as you go to cubics the fully expanded formula starts to get out of hand. A good explanation here: https://www.youtube.com/watch?v=N-KXStupwsc On Sun, 28 Mar 2021 at 22:53, jsc <[hidden email]> wrote: I did a little googling and found this discussion: https://www.physicsforums.com/threads/rectangle-inscribed-in-another-rectangle-solutions-for-all-cases.508715/ _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
Just to be clear, quartics are *actually* solved using an iterative method such as the bisection method posted by Parkinbot or Aberth's method (BOSL2), or perhaps a Newton iteration. Same with cubics. The quartic and cubic formula are not generally used and are probably going to lead to stability problems in general, even if you can implement them in a tidier way as a series of steps with reductions. And the above mentioned methods have the advantage of working on any polynomial (Aberth) or even any equation for the other two, so your labors are better repaid and you have less to debug. The above methods are also all simpler than the cubic and quartic formulas, or at least no more complex.
The only time it makes sense to use the cubic or quartic formulas is if for some reason you want a *symbolic* answer, like you want to know that it's really the cube root of 73 + the cube root of 42 or whatever, instead of having a numerical value. (And I mean that the symbolic answer is the end point...not that you're then going to use it to compute a numerical one.)
Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
Well for this problem the refactored Wolfram code posted by the OP gets the right answer for x such that when back substituted gives -1.81899e-12, so pretty close without any iteration. On Mon, 29 Mar 2021 at 00:46, adrianv <[hidden email]> wrote: Just to be clear, quartics are *actually* solved using an iterative method such as the bisection method posted by Parkinbot or Aberth's method (BOSL2), or perhaps a Newton iteration. Same with cubics. The quartic and cubic formula are not generally used and are probably going to lead to stability problems in general, even if you can implement them in a tidier way as a series of steps with reductions. And the above mentioned methods have the advantage of working on any polynomial (Aberth) or even any equation for the other two, so your labors are better repaid and you have less to debug. The above methods are also all simpler than the cubic and quartic formulas, or at least no more complex. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
An unstable algorithm *sometimes* produces bad answers, not always.
Why are you concerned about avoiding iterations? Also for your approach to be free of iteration you need to be using a square root algorithm, for example, that does not involve iteration. The methods I know use for computing square root are iterative.
Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
I was talking about iteration at the OpenSCAD level, not the hardware level. Perhaps it is slower than 9 inline expressions. But, yes, I can see raising terms to such high powers and then subtracting them might lead to numerical instability. The cube roots and square roots are probably what saves it. On Mon, 29 Mar 2021 at 11:24, adrianv <[hidden email]> wrote: An unstable algorithm *sometimes* produces bad answers, not always. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
In reply to this post by adrianv
Division also uses an iterative algorithm. Intel processors have a SQRT operation, which is roughly as expensive as division. OpenSCAD should be using the SQRT machine code operation if it is linked with glibc on linux. On Mon, Mar 29, 2021, at 6:24 AM, adrianv wrote:
_______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
In reply to this post by Troberg
Newcomer to the forum and beginner in both 3D printing and OpenSCAD, but this thread came early in my browsing.
Probably too late now, but for anyone who has not yet seen the animation: https://www.dropbox.com/s/670k2t3r25xz10w/AnimatedRectangle.mp4?raw=1 Sent from the OpenSCAD mailing list archive at Nabble.com. _______________________________________________ OpenSCAD mailing list To unsubscribe send an email to [hidden email] |
Free forum by Nabble | Edit this page |