Solving for an Integral: 28S solution « Next Oldest | Next Newest »

 ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 02-28-2008, 04:42 AM I'm referring here to Chuck's arc length problem which can be found here: Link to archive. The 28S solver and integrator work differently than in the later machines. Notably, integration can't be used as an algebraic function because it returns two values to the stack, the solution and the error. Failed algebraic attempt My first attempt was to create a wrapper to make integration an algebraic function: ```<< << -> fn var xl xh eps << fn var xl xh 3 ->LIST eps S DROP >> >> >> 'INTFN' STO ``` The local variables fn, var, xl and xh aren't neccessary if you want to execute INTFN from the stack, because the values are fetched from the stack to fill the variables and are than put back in the same order. But the overhead is required by the algebraic syntax. INTFN can be called like this: ```'INTFN(v/(1+4*SQ(T)),T,0,2.3,.0001)' EVAL ``` Due to the algebraic overhead, the function is relatively slow. Now enter the solver menu: ```'INTFN(v/(1+4*SQ(T)),T,0,X,EPS)-6' STEQ SOLVR ``` The resulting variables screen shows 'EPS', 'T' and 'X' as variables. In the parameterlist of INTFN, 'T' is a formal parameter, but the solver logic fails to see this and complains about an unset variable 'T' when I try to solve for 'X'. If I assign a value to 'T', integration fails bacause T is no longer accepted as a formal variable. So this is a dead end. :-( RPL solution I decided that I had to go the RPL way. Integration works without a formal variable if the integrand is a program. Here is my solution: ```<< -> x eps << << SQ 4 * 1 + v/ >> { 0 x } eps S DROP >> >> 'FN' STO ``` FN can be called as an algebraic: ```'FN(2.3,.0001)' EVAL ``` This is faster than the original INTFN version. Now enter the solver menu: ```'FN(X,EPS)-6' STEQ SOLVR ``` The resulting variables screen shows 'EPS' and 'X' as a variables. You can now set the accuracy for the integration and solve for 'X'. The result appears within seconds. :) I didn't try the ROOT command but it should work likewise. Marcus Edited: 28 Feb 2008, 4:43 a.m. ▼ Karl Schneider Posting Freak Posts: 1,792 Threads: 62 Joined: Jan 2005 02-28-2008, 03:25 PM Hi, Marcus -- Again, good work! You provided a solution for the HP-28C/S, which has different protocol and lesser capabilities for integration than the successor HP-48 series. I'd mentioned this in a post from 18 Feb 2004 within the thread "RPL challenge -- area between two curves" (Archive #14): Quote:Arnaud's program can't even be run on a 28C because its integration accepts arguments in a different format (dummy variable and limits enclosed in a list), which I don't know how to automate. It also didn't accept value-loaded variable names for the limits. You also demonstrated RPL-based integration of a program instead of an algebraic expression -- something I'd previously wondered how to do. That capability is very important, as some things aren't easy to incorporate into algebraic expressions. On the RPN side, the HP-32SII and its successors offer a choice of defining the integrand function by keystroke program or by algebraic expression; all other RPN models require a keystroke program. I trust that you noticed my solution for the TI-82 (14 Feb 2008), which I based upon your TI-92/Voyage solution. BTW, I find those one-liners clearer than the HP-71 solution, mainly for the absence of the latter's required variable names "IVAR" and "FVAR". However, the TI-82 (at least) cannot perform nested multiple integrations, while the older HP-71 can. -- KS Edited: 28 Feb 2008, 3:48 p.m. ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 03-02-2008, 10:49 AM Karl, thank you for the nice words. I was recently thinking about how the algorithm could be made a little speedier and came to a solution which works in RPL as well as for the 35s in RPN. The most time is spend in the numerical integration. The first improvement in this special case is to make use of the fact that the derivative of an integral is simply the integrand. Valentin has already jumped in and has modified his generic Newton solver for the 35s. Numerical solvers iterate in a way that the X value approximates the solution in ever smaller intervalls. This lead me to another simplification: Integral from 0 to Xn f(t)dt is the same as the sum of Integral from 0 to Xn-1 f(t)dt and Integral from Xn-1 to Xn f(t)dt. So, instead of computing the numerical integral between 0 and Xn in each iteration, the program can simply save a copy of the last argument and integral value and add the difference in each step. The integration is faster for smaller intervalls as long as the required accuracy is given as an absolute (not relative) value. It works out for both the 35s and 28S. My 35S-Implementation: ``` LBL F CHK=188F, LN=29 "Integrand RCL T x² 4 * 1 + v/ RTN LBL I CHK=A2B9, LN=19 "Standard integration FN= F I003: 0 x<>y S FN d T RTN LBL J CHK=CA57, LN=48 "Modified integration x!=0 GTO J008 STO Y "Initialize for X=0 STO Z FN= F RTN J008: FS? 1 "Flag 1 set: use standard intgration GTO I003 RCL Y "Last iteration for X as lower boundary x<>y STO Y "Save new lower boundary S FN d T STO+ Z "Summation RCL Z RTN LBL S CHK=1954, LN=80 "Solver 0 XEQ J001 "Initialize summation INPUT E "Absolute error in X for solver INPUT X "First guess S006: RCL X STO T "Argument to integrand XEQ F001 "Derivative of integral STO A "Intermediate storage RCL X "Current guess XEQ J001 "Integration, upper boundary on stack 6 - "Subtract right hand side of equation ENTER x<> A "Save result and get derivative 1/x * STO- X "Xn+1 := n - f(n+1) / f'(n+1) ABS RCL E "Maximum absolute error in X x x eps << x IF THEN << SQ 4 * 1 + v/ >> { XL x } eps S DROP RS + ELSE 0 END x 'XL' STO DUP 'RS' STO >> >> 'FN2' STO ``` If you solve 'FN2(X,EPS)-6' you'll see X, EPS, RS and XL on the menu. RS is the last ReSult, XL ist the Last X value. You need to initialize these by storing 0 to both or storing 0 to X and evaluating the expression once. The speed improvements are noticable: Depending on the given starting value and accuracy (by display format or variable EPS) the running times can be as little as 25% of the original version. Edited: 2 Mar 2008, 4:54 p.m.

 Possibly Related Threads... Thread Author Replies Views Last Post [HP-PIRME] BUG Pretty Print, Solution: abs => |...|, ABS => ||...|| CompSystems 2 1,323 12-13-2013, 09:36 AM Last Post: CompSystems [HP Prime] Tips for Solving Differential Equations More Efficiently Chris Pem10 8 1,923 11-21-2013, 08:25 PM Last Post: Chris Pem10 HP-PRIME CAS SOLVING fabrice48 8 1,945 10-19-2013, 01:21 PM Last Post: Han HP Prime Solving Nonlinear System of Equations for Complex Results Helge Gabert 11 2,697 09-30-2013, 03:44 AM Last Post: From Hong Kong Integral Richard Berler 32 5,306 09-22-2013, 08:59 PM Last Post: Les Koller [HP-Prime xCAS] Arrays: Programming Commands (solution to the ambiguity) CompSystems 6 1,693 08-30-2013, 06:32 PM Last Post: CompSystems HP-Prime: Solution to ambiguity with [ ... ] & another BUGs CompSystems 12 2,394 08-28-2013, 08:48 AM Last Post: CompSystems Advanced User's Manual and solution to the ambiguity of data types [HP-Prime xCAS] CompSystems 15 3,761 08-20-2013, 03:37 PM Last Post: Thomas Klemm Riemann's Zeta Function (HP-28S) Gerson W. Barbosa 8 1,837 02-03-2013, 03:23 PM Last Post: Gerson W. Barbosa [WP34s et al.] Solving the TVM equation for the interest rate Dieter 24 3,844 12-01-2012, 05:53 AM Last Post: Paul Dale

Forum Jump: