Small Solver Program - Printable Version +- HP Forums (https://www.hpmuseum.org/forum) +-- Forum: HP Calculators (and very old HP Computers) (/forum-3.html) +--- Forum: General Forum (/forum-4.html) +--- Thread: Small Solver Program (/thread-12423.html) Pages: 1 2 Small Solver Program - Gamo - 02-14-2019 05:25 AM Recently I try to make a small but versatile solver program specifically for HP Voyager Series calculator. This solver solve for f(x) = 0 and use some type of counter for iteration. Due to the iteration process this program will run slow on Original Calculator. With modern 12C, 15C, and all other emulator this will run extremely fast. This program is intended for HP-11C but this example will use HP-15C emulator so that this solver program result can be compare with the 15C build-in SOLVE function. --------------------------------------------------------- To run program. Enter 2 guesses closest to the root. X1 [ENTER] X2 [A] display answer // Use this Solver X1 [ENTER] X2 f [SOLVE] [E] display answer // Use 15C build-in SOLVE --------------------------------------------------------- Example for X^X = Y // X to the power of X equal to Y X on R1 // Register Storage 1 // For this Sover Y on R5 // Register Storage 5 X on R2 // Register Storage 2 // For 15C SOLVE Y on R5 // Register Storage 5 -------------------------------------------------------- Program: Code:  LBL A   // This Solver STO 0 Rv STO 1 90    STO I CLx LBL 1 RCL 0 / CHS RCL 1 + STO 1 DSE I GTO B RCL 1 RTN ------------------------------ LBL B  // f(x) program for X^X = Y RCL 1 LN RCL 1 x RCL 5 LN - GTO 1 ------------------------------- LBL E  // 15C SOLVE STO 2 // f(x) program for X^X = Y LN RCL 2 x RCL 5 LN - RTN ------------------------------------------- Example: X^X = 1000 STORE 1000 to R5 // [STO] 5 This Solver 4 [ENTER] 6 f [A] display 4.555535704 ------------------------------------------------------- 15C SOLVE 4 [ENTER] 6 f [SOLVE] [E] display 4.555535705 ------------------------------------------------------- Remark: For this solver answer is in R1 For 15C SOLVE answer is in R2 Gamo RE: Small Solver Program - Thomas Klemm - 02-14-2019 07:06 AM (02-14-2019 05:25 AM)Gamo Wrote:  Example: X^X = 1000 You can rewrite this equation such that $$x$$ is a fixed point of $$f(x)$$, that is $$x=f(x)$$: $$x=\frac{3}{LOG(x)}$$ Then you can iterate the following program with a starting value, say $$x=4$$: Code: 001-  43  13 :    LOG 002-       3 :    3 003-      34 :    x<>Y 004-      10 :    ÷ Example: 4 R/S 4.982892143 R/S 4.301189432 R/S 4.734933901 R/S 4.442378438 R/S 4.632377942 R/S 4.505830646 R/S 4.588735607 R/S 4.533824357 R/S 4.569933524 R/S 4.546075273 R/S 4.561789745 R/S 4.551417837 R/S 4.558254212 R/S 4.553744141 R/S 4.556717747 R/S 4.554756405 R/S 4.556049741 R/S 4.555196752 R/S 4.555759258 R/S 4.555388285 R/S 4.555632930 R/S 4.555471588 R/S 4.555577990 R/S 4.555507820 R/S 4.555554095 R/S 4.555523578 R/S 4.555543703 R/S 4.555530431 R/S 4.555539183 R/S 4.555533412 R/S 4.555537217 R/S 4.555534708 … 4.555535704 R/S 4.555535706 R/S 4.555535704 R/S 4.555535706 … From then on the result flips between these last two values. Cheers Thomas Addendum: Small Solver Program - Thomas Klemm - 02-14-2019 07:15 AM If you don't want to press R/S repeatedly: Code: 001-  43  13 :    LOG 002-       3 :    3 003-      34 :    x<>Y 004-      10 :    ÷ 005-  42  31 :    PSE 006-  43  40 :    x=0 RE: Small Solver Program - Albert Chan - 02-15-2019 12:07 AM (02-14-2019 07:06 AM)Thomas Klemm Wrote:  $$x=\frac{3}{LOG(x)}$$ 4 R/S 4.982892143 R/S 4.301189432... It would be nice if we can temper the oscillation, or slow convergence. Let x0 = 4, x1, x2 = the first two iterated values. Rate = (x2-x1)/(x1-x0) ≈ (4.30 - 4.98) / (4.98 - 4) = -0.694 If the same trend continued, we expect final % = 1/(1-r) ~ 60% x∞ ≈ x0/(1-r) = (x1 - (x1-x0))/(1-r) = (x1 - r x0) / (1-r) Use weighted fixed-point equation x = 0.6 * 3/log10(x) + 0.4 x 4.6 4.555924149 4.555537395 4.555535712 4.555535705 (converged) Edit: compare with Newton's method, x = (ln(1000) + x) / (ln(x) + 1) 4 4.571001573 4.555546101 4.555535705 (converged) RE: Small Solver Program - Thomas Klemm - 02-15-2019 06:26 PM (02-15-2019 12:07 AM)Albert Chan Wrote:  It would be nice if we can temper the oscillation, or slow convergence. We can also use Aitken's delta-squared process to accelerate the speed of convergence: $$a_{n}=x_{n+2}-\frac {(x_{n+2}-x_{n+1})^{2}}{(x_{n+2}-x_{n+1})-(x_{n+1}-x_{n})}$$ This leads to the following program for the HP-11C: Code: 001-42,21, 1 :    ▸LBL 01 002-      36 :     ENTER 003-    32 0 :     GSB 00 004-      30 :     - 005-  43  36 :     LSTx 006-      36 :     ENTER 007-    32 0 :     GSB 00 008-      30 :     - 009-  43  36 :     LSTx 010-      34 :     x<>y 011-  43  33 :     R↑ 012-      34 :     x<>y 013-      30 :     - 014-  43  36 :     LSTx 015-  43  11 :     x² 016-      34 :     x<>y 017-      10 :     ÷ 018-      30 :     - 019-  43  32 :     RTN 020-42,21, 0 :    ▸LBL 00 021-  43  12 :     LN 022-    45 0 :     RCL 00 023-      34 :     x<>y 024-      10 :     ÷ 025-  43  32 :     RTN Example: Solve: $$x^x=1000$$ 1000 LN STO 0 4.6 GSB 1 4.555665779 R/S 4.555535706 R/S 4.555535705 R/S Error 0 Feel free to add a check for 0 after line 13 to prevent the error and loop back to label 1 instead of returning in line 19. Cheers Thomas RE: Small Solver Program - Albert Chan - 02-15-2019 09:16 PM (02-15-2019 06:26 PM)Thomas Klemm Wrote:   (02-15-2019 12:07 AM)Albert Chan Wrote:  It would be nice if we can temper the oscillation, or slow convergence. We can also use Aitken's delta-squared process to accelerate the speed of convergence: $$a_{n}=x_{n+2}-\frac {(x_{n+2}-x_{n+1})^{2}}{(x_{n+2}-x_{n+1})-(x_{n+1}-x_{n})}$$ Amazingly, my rate formula is same as Aitken extrapolation formula ! Assuming we have 3 values, x0,x1,x2 and tried to guess where it should end up. My rate analysis: r = (x2-x1)/(x1-x0) = Δx1 / Δx0 We need this for the proof: (Δx0)² = ((Δx0 - Δx1) + Δx1)² = (Δx0 - Δx1)² + 2 * Δx1 (Δx0 - Δx1) + (Δx1)² If same trend continue, where it will ends up = x0 + Δx0 * 1/(1-r) = x0 + (Δx0)² / (Δx0 - Δx1) = x0 + (Δx0 - Δx1) + 2 * Δx1 + (Δx1)² / (Δx0 - Δx1) = x0 + (x1-x0) + (x1-x2) + 2*(x2-x1) − (Δx1)² / (Δx1 - Δx0) = x2 − (Δx1)² / (Δx1 - Δx0) RE: Small Solver Program - Thomas Klemm - 02-16-2019 03:58 AM (02-15-2019 09:16 PM)Albert Chan Wrote:  Amazingly, my rate formula is same as Aitken extrapolation formula ! That's not really a surprise, is it? It's essentially the same idea: geometric series $$1\,+\,r\,+\,r^{2}\,+\,r^{3}\,+\,\cdots \;=\;{\frac {1}{1-r}}$$ RE: Small Solver Program - Csaba Tizedes - 02-16-2019 12:24 PM (02-14-2019 05:25 AM)Gamo Wrote:  a small but versatile solver program specifically for HP Voyager Series calculator. Example for X^X = Y // X to the power of X equal to Y Yes, the classic problem: which number you can power to itself to get the largest number on your calculator: x^x=10^100. Do LOG() of both sides: x×LOG(x)=100, then do the iteration: x_0:=50, then x_i+1:=100/log(x_i). The simple fixed point iteration and converges. By the way, if you interested, here is my SOLVE(i) for HP-12C, which can solve any equation for ANY variable like on 15C, without rewriting the equation. Looks like indirect addressing without indirect addressing. 28 steps only. Yes, my English is terrible, but more convenient than read the full material in Hungarian. At the end of pdf you can find lots of applications which solved with this little solver (like ODE solve with EULER+SOLVE(i) on 12C (!!!), pressure vessel pneumatic conveying pipe sizing program, humidity and dew point solving, outer temperature calculation of a concrete silo wall or comparison of present values of educated and not educated people salaries with the little SOLVE combined with 12C's cash flow calculation capabilities - I guess this is shows clearly the power of this little routine possibilities and why I say everywhere, the SOLVER is MUST to implement ALL current calculators, I hope somebody read and understand what I want to say). The small secant method in this paper for CASIO fx-50F meantime reduced to 13 steps, I guess the smallest solver for these blind machines. You can find attached. Csaba RE: Small Solver Program - Thomas Klemm - 02-16-2019 01:42 PM (02-16-2019 12:24 PM)Csaba Tizedes Wrote:  By the way, if you interested, here is my SOLVE(i) for HP-12C Steps of secant method from 03 to 18: Code: 03-   45 2 :    RCL 2 04-   45 0 :    RCL 0 05-     30 :    - 06-  43 36 :    LSTx 07-   44 2 :    STO 2 08-     35 :    CLx 09-   45 3 :    RCL 3 10-   45 1 :    RCL 1 11-     30 :    - 12-  43 36 :    LSTx 13-   44 3 :    STO 3 14-     33 :    R↓ 15-     10 :    ÷ 16-   45 1 :    RCL 1 17-     20 :    × 18-44 30 0 :    STO+ 0 This can be shortened to: Code: 03-   45 2 :    RCL 2 04-   45 0 :    RCL 0 05-   44 2 :    STO 2 06-     30 :    - 07-   45 3 :    RCL 3 08-   45 1 :    RCL 1 09-   44 3 :    STO 3 10-     30 :    - 11-     10 :    ÷ 12-   45 1 :    RCL 1 13-     20 :    × 14-44 30 0 :    STO+ 0 Don't hesitate to publish your programs in this forum. You can use Google to translate to English. Cheers Thomas RE: Small Solver Program - Csaba Tizedes - 02-16-2019 03:24 PM (02-16-2019 01:42 PM)Thomas Klemm Wrote:  Steps of secant method from 03 to 18 can be shortened to. Thanks, this was eons before, I am sure this is only the first iteration and required to fine tuning. (02-16-2019 01:42 PM)Thomas Klemm Wrote:  Don't hesitate to publish your programs in this forum. You can use Google to translate to English. OK, sometimes I feel myself as Don Quixote... Frankly speaking, who want to discuss about calculators and calculator programming now? 15 years before I wrote something useful, today I am only (an angry) user... Csaba RE: Small Solver Program - Gamo - 02-17-2019 02:57 AM I just fine tune this solver program for HP-11C Added the iterations limit to 28 loops count Added the condition test to end program when root is found. If iteration is over 28 loops count this indicate that user must use a closer guess. ------------------------------------------- Procedure: Use R1 for unknown X in your equation. End your formula with [RTN] X1 ENTER X2 [A] display Answer Answer is in R1 ------------------------------------------------------------------------ Program: Code:  LBL A STO 0 Rv STO 1 28 STO I  // Store Loop Count Limit CLx --------------------- LBL 0 GSB B GSB 1 STO 1 GSB B GSB 1 RCL 1 X=Y  // End if Root is found GTO 2 DSE  // Loop Limit Counter GTO 0 CLx FIX 9  // 0.000000000 indicate that Maximum Loops is use up. PSE PSE FIX 4  GTO 2 --------------------- LBL 1  // Solver Routine RCL 0 ÷ CHS RCL 1 + RTN -------------------- LBL 2 RCL 1  // Answer  RTN -------------------- LBL B .  // f(x) equation start here . .  // X = Register 1 [R1] . . RTN Gamo RE: Small Solver Program - Thomas Klemm - 02-17-2019 09:06 AM (02-17-2019 02:57 AM)Gamo Wrote:  I just fine tune this solver program for HP-11C Your program doesn't work for this simple example: $$x^2+x=x(x+1)=0$$ I've entered the following program: Code: LBL B 1 RCL 1 + LSTx × RTN Example: -2 ENTER -0.5 A 9.999999 99 The reason is that the function value is divided by -0.5, making it larger and larger: Code: LBL 1  // Solver Routine RCL 0 ÷ CHS RCL 1 + RTN Another thing is this in this loop: Code: LBL 0 GSB B GSB 1 STO 1 GSB B GSB 1 RCL 1 X=Y  // End if Root is found GTO 2 DSE  // Loop Limit Counter GTO 0 When you return to label 0, the same value in register 1 will be used as before. This function evaluation can be avoided if you save the result of the previous evaluation. You could move this block: Code: LBL 2 RCL 1  // Answer  RTN Over here and avoid the jump GTO 2: Code: CLx FIX 9  // 0.000000000 indicate that Maximum Loops is use up. PSE PSE FIX 4  LBL 2 RCL 1  // Answer  RTN Take a look at Csaba's solution for the HP-12C using the secant method. Cheers Thomas RE: Small Solver Program - Gamo - 02-17-2019 02:33 PM Thanks Thomas Klemm I try single guess and it work. CHS 2 [A] display -1 If using two negative guesses and result is Overflow then use Single guess instead. Formula Routine is RCL 1 1 + RCL 1 X Gamo RE: Small Solver Program - Thomas Klemm - 02-17-2019 04:57 PM (02-17-2019 02:33 PM)Gamo Wrote:  If using two negative guesses and result is Overflow then use Single guess instead. Then try this function $$f(x)=4x(x+1)$$: Code: LBL B RCL 1 1 + RCL 1 × 4 × RTN Examples: -0.5 ENTER 0.5 A 9.999999 99 -0.5 ENTER A 9.999999 99 0.5 ENTER A 9.999999 99 Even initial guesses very close to the solutions lead to the same result: -1.00001 ENTER -0.99999 A 9.999999 99 -0.00001 ENTER 0.00001 A 9.999999 99 Your algorithm is fundamentally broken. Cheers Thomas RE: Small Solver Program - Gamo - 02-18-2019 03:49 AM Your're right that only work for certain situations. Gamo RE: Small Solver Program - Thomas Klemm - 02-18-2019 05:20 AM If you compare your algorithm with Newton's method: $$x_{n+1}=x_{n}-\frac {f(x_{n})}{f'(x_{n})}$$ you may notice that the number in register 0 should in fact be $$f'(x_{n})$$ or at least a good approximation thereof. And with this in mind indeed the solution can now be found: $$f'(x) = \frac{d}{dx}4 x (x + 1) = 8 x + 4$$ $$f'(-1) = -4$$ $$f'(0) = 4$$ Examples: -1.5 ENTER -4 A -1.00000 0.5 ENTER 4 A 0.00000 Thus your 2nd parameter is an estimate of the slope at the root. And then you don't really have to call the function twice with the same value: Code: LBL A STO 0 Rv STO 1 28 STO I  // Store Loop Count Limit LBL 0 GSB B RCL 0 ÷ STO - 1 X=0  // End if Root is found GTO 1 DSE  // Loop Limit Counter GTO 0 CLx FIX 9  // 0.000000000 indicate that Maximum Loops is use up. PSE PSE FIX 4 LBL 1 RCL 1  // Answer  RTN LBL B .  // f(x) equation start here . .  // X = Register 1 [R1] . . RTN So the question remains how to approximate the derivation at the roots. If you keep record of the previous function evaluation you can get an estimate of the slope using the secant method. Cheers Thomas RE: Small Solver Program - Dieter - 02-18-2019 07:46 PM (02-18-2019 05:20 AM)Thomas Klemm Wrote:  If you compare your algorithm with Newton's method: $$x_{n+1}=x_{n}-\frac {f(x_{n})}{f'(x_{n})}$$ you may notice that the number in register 0 should in fact be $$f'(x_{n})$$ or at least a good approximation thereof. Yes. Gamo already posted this method some time ago in the General Software Library. I have mentioned several times that the first input it NOT a guess for the root but an estimate for the function's DERIVATIVE. For more details take a look at this post: http://www.hpmuseum.org/forum/thread-12344-post-111680.html#pid111680 Gamo, you really should correct and update the instructions. Here it still says: Quote:Enter 2 guesses closest to the root. But again: that's NOT true. The first input is NOT a guess for the root. You now have seen a variety of cases that prove it. Dieter RE: Small Solver Program - Thomas Klemm - 02-18-2019 10:22 PM (02-18-2019 07:46 PM)Dieter Wrote:  Yes. Gamo already posted this method some time ago in the General Software Library. I have mentioned several times that the first input it NOT a guess for the root but an estimate for the function's DERIVATIVE. You could think him stubborn. But it made Csaba post his program which made me translate some of his papers from Hungarian to English and stumble totally unrelated upon one of his previous posts: Quadratic fit without linear algebra (32SII) - NOT only for statistics lovers - LONG Which I still do not understand in detail. So all in all I'm happy with the outcome. Cheers Thomas RE: Small Solver Program - Albert Chan - 02-19-2019 01:10 AM Quadratic fit without linear algebra (32SII) - NOT only for statistics lovers - LONG Linear regression on slope had another benefit, which my Casio cannot do. Doing Quadratic Regression on Casio, you don't really know if it is any good. There is no "correlation coefficient" to lookup. (Why?) Linear regression on slope can show whether a Quadratic curve is warranted. I think the code add slope data, like this: Data points (0,90) and (30, 90.2), add to linear regression: (30+0)/2, (90.2-90)/(30-0) Next point (60, 86.4), add to linear regression: (60+30)/2, (86.4-90.2)/(60-30) ... N points (sorted) thus generate N-1 slopes to regress. RE: Small Solver Program - Csaba Tizedes - 02-19-2019 08:39 AM (02-18-2019 10:22 PM)Thomas Klemm Wrote:   (02-18-2019 07:46 PM)Dieter Wrote:  Yes. Gamo already posted this method some time ago in the General Software Library. I have mentioned several times that the first input it NOT a guess for the root but an estimate for the function's DERIVATIVE. You could think him stubborn. But it made Csaba post his program which made me translate some of his papers from Hungarian to English and stumble totally unrelated upon one of his previous posts: Quadratic fit without linear algebra (32SII) - NOT only for statistics lovers - LONG Which I still do not understand in detail. So all in all I'm happy with the outcome. Cheers Thomas I feel myself like a star (at least for 15 minutes...) Yes, this is my "linear regression on slope", as Chan wrote. When I was in my "Summer Practice" at a company I have modelled lots of water network flows between pump houses (sources), consumers and reservoirs, and sometimes the characteristics curves given only in paper form. Thats why I wrote this little routine. As you can see, if the measured points has significant error, the slope will be very inaccurate then the regression also inaccurate and the coefficients can be meaningless. But as a programming question this was an interesting "game" for me. To the original question from Gamo: I will check carefully the lists here, but I have another idea: what if, if we just simple walking along on the function's curve until the sign is changing, then we turning back, decreasing the stepsize and walking again. If the sign of the function changed, we are turning back again, decreasing the size of the steps and so on... If the stepsize is less than a small value, we found the root. This "marching method" not so elegant, not so fast, but good to play with it, just for a try - maybe it can be shorter than the other methods. (I have an exactly same method for finding a minimum point of a function for CASIO fx-3650P/Sencor SEC103; perfectly works, like a ball at the bottom of a hole, slowly find the deepest point, when swings around the extremum). Csaba