Another tough integral « Next Oldest | Next Newest »

 ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-13-2012, 03:09 PM ``` 1 cos-1(cos2(pi*x)) Let f(x) = --- + -------------------, 2 3*sin(pi*x) ``` a fair approximation to the square wave function with an offset of 1/2 below: ``` y | | ______ 1 |______ ______ ______ | | | | | | | | | | | | | | | | | | | | | | | | ----------------+------------------------------------- x -2 -1 0 | 1 2 3 4 5 | ``` Let us find the integral of f(x) from x=0 to x=1. WolframAlpha returns only a numerical result, 0.995555, which suggests there is no known closed-form solution. Let's see how some calculators handle the task: ```WP 34S: 001 LBL 00 002 # pi 003 * 004 COS 005 x^2 006 ACOS 007 x<>y 008 # pi 009 * 010 SIN 011 3 012 * 013 / 014 . 015 5 016 + 017 RTN 0.9955554977636525 wp34s (v2) 0.9955554973850958 wp34s (v3.0 2458) ``` HP-15C LE: ```0.9955555102 (FIX 7) ``` HP-50g: ```'Integral(0,1,1/2+ACOS(SQ(COS(pi*X)))/(3*SIN(pi*X)),X)' 0.995555497542 (Fix 9) ``` Emu71: ```>RADIANS >STD >INTEGRAL(0,1,1E-10,1/2+ACOS(COS(PI*IVAR)^2)/(3*SIN(PI*IVAR))) 0.995555497542 >INTEGRAL(0,1,1E-11,1/2+ACOS(COS(PI*IVAR)^2)/(3*SIN(PI*IVAR))) 0.995555198305 ``` It looks like our calculators arsenal cannot give us more than 8 correct decimal places, 0.995555497. (I haven't tried the latest WP 34S release yet). Edited for grammar. Edited: 13 May 2012, 3:22 p.m. ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 05-13-2012, 03:59 PM The latest WP 34S release should return the same result but you can try in double precision and compare the figures. Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-13-2012, 06:57 PM Quote: ``` 1 cos-1(cos2(pi*x)) Let f(x) = --- + -------------------, 2 3*sin(pi*x) ``` a fair approximation to the square wave function with an offset of 1/2 below: ``` y | | ______ 1 |______ ______ ______ | | | | | | | | | | | | | | | | | | | | | | | | ----------------+------------------------------------- x -2 -1 0 | 1 2 3 4 5 | ``` Let us find the integral of f(x) from x=0 to x=1. ``` ``` Well, it doesn't seem that tough at all, a simple integration routine I wrote for high-precision returns: ``` 0.99555549776236782298677656111175223160096264596209... ``` essentially instantly (0:00:00, to the second), which is correct to all places shown. A second or so will get you 200 decimal places beginning with: ``` 0.995555497762367822986776561111752231600962645962092224364259066155312626099042699264864294947084302054786866956732685735494339225156047080022695529860042854... ``` This might be useful to check for accuracy the results you get from calculators. Best regards from V. ``` ``` ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-13-2012, 09:30 PM Hello Valentin, Quote: A second or so will get you 200 decimal places beginning with: ``` 0.995555497762367822986776561111752231600962645962092224364259066155312626099042699264864294947084302054786866956732685735494339225156047080022695529860042854... ``` I've checked the digit frequency of your result and compared it against that of Pi for the first 156 significant digits: ``` your Pi result 0 15 14 1 8 14 2 24 20 3 9 17 4 16 16 5 18 14 6 26 13 7 11 10 8 10 21 9 19 17 ``` I would expect the result of the integral to have a more even digit distribution. Of course there is no solid reason behind this. Often integrals of periodic transcendental functions, depending on the limits of integration, are not transcendental numbers themselves, much less normal. The repetition of the digit '5' in the beginning might be an indication the result is not normal (as I may have wrongly assumed). Not that I don't trust it, but is there another integration method that may confirm your result? Best regards, Gerson. ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-14-2012, 02:02 AM Quote: Not that I don't trust it, but is there another integration method that may confirm your result? It might be the case that the meager 156 digits I gave aren't enough for an statistically significant digit frequency calculation so I let my routine run for a minute and there you are, about 1,000 decimal places for you to fiddle with to your heart's content ``` 0.99555549776236782298677656111175223160096264596209222436425906615531262609 90426992648642949470843020547868669567326857354943392251560470800226955298600428 54061223230260667542990343571102929104397030422418588620801216306986461891226713 12458538366567533281556903416297515134295492756873107350286487265699127726709068 02866195453154345434350558071871988003701209930033784114664742830891982441996720 14204397003759667135107370215878694197174112228164086756751024279765262148865802 14426869108178025469920516054579901976743735865432183534421777147289493851489995 34900344460945637494513818059645286910899589401623806133301582292846669833025305 75559408429314709927714148536378809730554519179542358963419925913618096446986160 89707018534727889549943105208309888722681503586468722935816711860533157811264171 46873750916264977276849274618602746297307576605027226287396682942279959930995857 63544972773656734777298255881221987614643146784249414791351811025549311074079353 91817654555749920736361474456647892716652144730718950692026594512267763789970029... ``` For reliable digit frequency results it's frequently the case that many thousands or even millions of decimal digits are necessary, and even more. Consider for instance this simple square root computation: ``` sqrt(308642) = 555.55557777777733333335111111022222227199999... ``` Seeing that 44-digit result for its 6-digit integer argument you might be induced to believe that it is anything but normal, with obviously biased digit counts. However, the mirage quickly fades away the moment you look further: ``` 555.555577777777333333351111110222222271999997013333521066654464000813511055792 35937757839912506782260202536969737338174042982693335359822074938022460559677134 80453936671895444061294508193099456043131900500139679957197210080367577711228344 88033492927605660571893841517753635572851565452015965971830849488612380140811581 65086640947765975127949499079439218226790840268539816815091386804202414469605745 72371836881741155667721985539461145895339203043493051583779961546365025268947444 55853751944622774941855495715419031133120938831434794450531321168458752080467791 63014094929372714096921707984170197616439799857571395157854920325300060444559270 50627657849491643773882305857452761242714492527784638035394215448725969098842689 41210190584508455675940399110409463416619952140752221724065975226358154445769047 20955493675153265579204818884769350972126051641884022716404454461648737349442853 71015481254270598108264574726782545171053067106881548003259525083330892109951136 99664941708194946290983595047955472874158265042111361146447247573233473458118379 10816149125024610778070772206907264334574601215626802573125232732663840608892935 67078985853834062943931460643680720460050842175328769283886123614901626961378958... ``` which surely looks much more "normal" and about as normal as an irrational square root is going to be. Making judgments based on scarce data is what Wikipedia calls "Hasty generalization, a logical fallacy also known as 'the law of small numbers'" Best regards from V. ``` ``` ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-14-2012, 06:03 AM Thank you, Valentin, for taking the time to compute a staggering one thousand digits of the integral and providing yet another interesting example and explanations. Best regards, Gerson. ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-14-2012, 06:55 AM Quote: Thank you, Valentin, for taking the time to compute a staggering one thousand digits of the integral and providing yet another interesting example and explanations. ``` ``` You're welcome, Gerson, but this does not deserve any thanks, it was just a 1-min computation plus the time to copy & paste the result to post it. What's weird about it all is that I can perform such computations, i.e., multi-digit integrals, matrix operations, general programming, etc., not just on a desktop PC, not just on a laptop, not just on a tablet or even a mobile phone, but also on a device as simple as my eReader ! So I think handheld calculators are really, really a thing of the past and can't compete at all with modern technology. Just imagine, doing pretty convoluted math on a humble eReader which weights even less than a typical advanced graphing calculator, doesn't tax your eyes, consumes extremely little power, and which you can carry with you at all times with no burden whatsoever. We old-timers are now living in the future ... :D Best regards from V. ``` ``` ▼ Oliver Unter Ecker Member Posts: 239 Threads: 9 Joined: Dec 2010 05-14-2012, 12:24 PM Hi Valentin, Quote: Just imagine, doing pretty convoluted math on a humble eReader which weights even less than a typical advanced graphing calculator, doesn't tax your eyes, consumes extremely little power, and which you can carry with you at all times with no burden whatsoever. Yes, I second that. Were you extrapolating, or do you have an actual implementation on your eReader? If so, native code or JavaScript? If latter, I'd be very interested to add your integration routine to my recent BigFloat implementation, if you'd be willing to share it. Cheers. ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-14-2012, 06:42 PM Quote: Were you extrapolating, or do you have an actual implementation on your eReader? If so, native code or JavaScript? No, no, I wasn't extrapolating, my high-precision numerical integration routine actually does run flawlessly on my Sony eReader (10.4 ounces ~ 295 grams) under DosBox for Android, as I wrote it many years ago for an MS-DOS environment. Matter of fact lots of pretty useful MS-DOS software runs as well (Turbo Pascal, for instance, DOS-based emulators, etc.), as well as lots of native Android applications. As long as the particular app isn't heavy with animations (Chess, Sudoku, crosswords, word-search puzzles, comics readers, audio players, etc), it'll run fine on its beautiful E-Ink Pearl 800x600 screen. Thanks a lot for your interest and Regards from V. ``` ``` ▼ Oliver Unter Ecker Member Posts: 239 Threads: 9 Joined: Dec 2010 05-15-2012, 01:58 AM Ah, I see. Thanks for the clarification. I thought of JavaScript because it's the one language that comes built-in with these devices, as part of their browsers. I share your enthusiasm about what could be done with eReaders and calcs. Edited: 15 May 2012, 2:00 a.m. Crawl Senior Member Posts: 306 Threads: 3 Joined: Sep 2009 05-14-2012, 01:02 PM Quote: So I think handheld calculators are really, really a thing of the past and can't compete at all with modern technology. Just imagine, doing pretty convoluted math on a humble eReader which weights even less than a typical advanced graphing calculator, doesn't tax your eyes, consumes extremely little power, and which you can carry with you at all times with no burden whatsoever. We old-timers are now living in the future ... :D Best regards from V. I'm not sure why the moral here is "handheld calculators are a thing of the past" rather than "handheld calculators of the future have something to learn from existing ereaders". ▼ bill platt Posting Freak Posts: 2,448 Threads: 90 Joined: Jul 2005 05-15-2012, 07:56 AM It's the "convegence theorem" ▼ Crawl Senior Member Posts: 306 Threads: 3 Joined: Sep 2009 05-15-2012, 09:17 PM I love the idea of a capacitive touch screen for a calculator. For one, I'd like to use it literally like scratch paper for "manual" calculations. I'd also like to see a way of inputting formulas through handwriting recognition. And I'd like to see these two capacities work together. DeboT Member Posts: 59 Threads: 3 Joined: Oct 2010 05-14-2012, 06:57 AM Thank you for your insightfulness. However, are the 1000 digits generated from a different method than to 200 digits? If not, then you have not answered Gerson's question (you have merely refuted the basis of his query). Quote: Making judgments based on scarce data is what Wikipedia calls "Hasty generalization, a logical fallacy also known as 'the law of small numbers'" The same goes for the number of answers from fundamentally different sources to affirm the accuracy. So far, yours is the only one more than 16 digits. Quote: This might be useful to check for accuracy the results you get from calculators. It would be very stupid indeed for anyone to take your single source as a definitive answer to check all others against. Perhaps as a comparison, but certainly not a check. I do not even use Wolfram Alpha answers as definite checks, but merely as comparisons.So far from this thread, the best we can come up with is the first 6 digits (0.995555) and then only with a very low confidence due to the small sample. I.e. it could all be horribly wrong... (unlikely, but possible). Steve Perkins Member Posts: 104 Threads: 0 Joined: Dec 2007 05-14-2012, 07:48 PM I can confirm the result given by Valentin to the first 30 digits of precision, using a high-precision method I programmed a few years back. Adaptive integration using double exponential quadrature. My method was designed for tough integrals, but doesn't achieve results very quickly in high precision applications. It took about a minute to get (and confirm) 30 digits. As always I admire Valentin's abilities and the tools he creates. Edited: 14 May 2012, 7:51 p.m. ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-14-2012, 11:36 PM Quote: I can confirm the result given by Valentin to the first 30 digits of precision, using a high-precision method I programmed a few years back. Adaptive integration using double exponential quadrature. My method was designed for tough integrals, but doesn't achieve results very quickly in high precision applications. It took about a minute to get (and confirm) 30 digits. As always I admire Valentin's abilities and the tools he creates. First of all, thank you very much for your appreciation and kind words. It would be simple to check my 1000+ digit result using Mathematica, Maple, or PARI/GP for instance but I have none of them at hand right now, just my eReader with 894 books, 1200+ (mostly math) papers, 550+ audio files, and about two dozen assorted applications installed on it which, regrettably, don't include any of those. However, I have no doubts whatsoever that my results are correct, as I've thoroughly tested the routine innumerable times using much more troublesome cases than Gerson's fairly mild integral, with reasonable success so far. Just for instance, this simple, elementary integral is probably a harder nut than Gerson's if dealing with it straightforwardly: ``` / 10 | x2*sin(x3).dx /0 ``` To some 750+ decimal places, my quadrature routine produces: ``` 0.14587364123643233630725025779820134374806272608726769409905827138475545678 91873252258156604930609502083829281765371738264413094048636339249999429802062379 57961603563346425312582138935758511981714815713913326990354637466454659373462989 23856880365641819094067543287366102194594287500142135303138814301671482938831794 27974049757460366605942774949622910049520201352146447612765559727085467322830992 54379599727222565251775009728250768063382128194525269251261178725736672752044626 44859088056137118266467542972367110856139562847936891992188412495430195619961309 13352677329362058906113949895974412408123236712255902834497437414519914257696481 63439748933708957007997706981128617570867623590746630295371101538752591300416917 574665836311221415826565867926678856523022042718456... ``` which is demonstrably correct to all digits shown. On the other hand, the naive approach using the HP-71B and asking for full precision gives: ``` >INTEGRAL(0,10,0,IVAR^2*SIN(IVAR^3)) .145873640643 ``` An interval approach fares somewhat better: ``` >S=0 @ FOR I=0 TO 9 @ S=S+INTEGRAL(I,I+1,0,IVAR^2*SIN(IVAR^3)) @ NEXT I @ S .145873641554 ``` You might want to try assorted calculators and your own high-precision quadrature to see how they cope with this simple integral. Best regards from V. ``` ``` Edited: 14 May 2012, 11:44 p.m. ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-15-2012, 12:19 AM Quote: However, I have no doubts whatsoever that my results are correct, as I've thoroughly tested the routine innumerable times using much more troublesome cases than Gerson's fairly mild integral, with reasonable success so far. I should have known better. Sorry! I have yet to see a wrong result coming from you. Well, at least this has elicited some interesting discussions and insights. Best regards, Gerson. Mike Morrow Posting Freak Posts: 758 Threads: 9 Joined: Jul 2007 05-15-2012, 10:10 AM The Casio fx-115ES Plus ("The Quadra-nator") confirms the accuracy of the first 11 significant digits shown above. :-) But seriously...when programmed in LineIO mode (using the 115's default integration tolerance value) as S(x2sin(x3),0,10) it returns 0.145873641232327 (accurate to 11 digits). The 10 digits shown in the display are correct. Computation time is 1003 seconds...just making it before the 115's 1050-second time-out halt. Better have some bright light on the solar cell to do many of these! The equivalent integral (Ssin(u)du)/3 evaluated from 0 to 1000 produced 0.145873641237629, a result accurate only to 11 digits. It took 737 seconds to compute. When solved and computed analytically as -(cos(1000)-cos(0))/3 on the 115, the result is 0.145873641236197. This is accurate only to 12 digits. I was expecting a better result. Edited to add results from 115's analytic solution. Edited: 16 May 2012, 1:01 a.m. Steve Perkins Member Posts: 104 Threads: 0 Joined: Dec 2007 05-15-2012, 02:47 PM That is certainly a numerically challenging function to integrate over the suggested range. It increases while oscillating more and more as x increases. And it has the benefit of being analytically easy to solve with a simple substitution. My high-precision routine had the grace to give up, since it quickly hit built in limits. Other double precision integrators in my toolkit had mixed results. Some gave close results for lower precision demands, but hung (taking too much time) when trying for precision better than 1E-14 or so. Other integrators gave close answers, but with error estimates that were off by a factor of 10 or more. A nice test case indeed. And motivation to attempt to make my tools more robust. ▼ Mike Morrow Posting Freak Posts: 758 Threads: 9 Joined: Jul 2007 05-15-2012, 10:00 PM Quote: A nice test case indeed. And motivation to attempt to make my tools more robust. I wish there were at least one HP calculator with a numerical integration function that could even remotely match the performance found in my cheap non-programmable calculator that is sold in mass-market department stores like Walmart and Target in the USA. It's embarrassing that all of the HP line, even the 50g, performs so poorly when executing its native quadrature firmware. Like Casio, the TI line has utilized Gauss-Kronrod quadrature for many years. It would be interesting to try this test case on my TI-89 Titanium HW version 4, but I don't have it nearby at present. It is HP that needs its tools to be more robust. Oliver Unter Ecker Member Posts: 239 Threads: 9 Joined: Dec 2010 05-16-2012, 03:33 AM I confirmed your 750 digit result with the following W|A command-line: NIntegrate[x^2*sin(x^3), {x, 0, 10}] 1000 digits Strangely enough, this command-line does not work with Gerson's function ("1/2+(acos(cos(pi*x)^2))/(3*sin(pi*x))" or "1/2+(ArcCos[Cos[Pi*x]^2])/(3*Sin[Pi*x])" for Mathematica syntax). Instead, this NIntegrate[1/2+(acos(cos(pi*x)^2))/(3*sin(pi*x)), {x, 0, 1},WorkingPrecision->100] will work, but be limited to 200 digits. While NIntegrate[x^2*sin(x^3), {x, 0, 10},WorkingPrecision->100] also works, it doesn't work with WorkingPrecision->1000 (though the first command proves that W|A can compute 1000 digits without timing out). If someone can make rhyme or reason out of this, I'd love to hear it. I'm presently working on a CAS for ND1 that uses W|A to do the hard work. This thread inspired me to sort out high-precision numerical integration, and the result of a little bit of integration work is that RPL code like this ```\<< 200 setBigFPrecision 0 toBig 10 'x^2*sin(x^3)' 'x' integrate \>> ``` will numerically integrate a function to a desired number of digits. The "integrate" function, when fed BigFloat values, will produce a BigFloat result. A type in ND1, it can be used for subsequent computation. The internal implementation of "integrate" looks like this ``` "integrate": function(a, b, expr, name) { return BigFloat.fromString(calculator.callWA("", "NIntegrate[" + calculator.expressionToWAexpression(expr) + ", {" + calculator.unquote(name) + ", " + a + ", " + b + "},WorkingPrecision->" + BigFloat.precision + "]")); }, ``` The idea is that users can write their own adapter functions and thereby compute anything they like in Mathematica/WolframAlpha, while using a stack-based calc that operates very much like a 50g. ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-16-2012, 05:32 AM Quote: I confirmed your 750 digit result with the following W|A command-line: NIntegrate[x^2*sin(x^3), {x, 0, 10}] 1000 digits ``` ``` Thanks, you can get as many digits as you want for checking purposes or otherwise by simply evaluating: ``` (1-Cos(1000))/3 ``` which is the exact result. Regards from V. ``` ``` Steve Perkins Member Posts: 104 Threads: 0 Joined: Dec 2007 05-16-2012, 07:58 PM Valentin, I'd be quite interested in learning any details you could share about the workings of your quadrature routine, as well as your extended precision setup. I've programmed adaptive quadrature using Simpson's and Gauss+Lobatto (which both require endpoint evaluation). Also programmed Gauss+Kronrod, and Gauss+Bond routines which don't require endpoint evaluation. For these I have fixed weights and abscissa in the code making variable precision hard. I also programmed double exponential, which works much differently than the above methods, but lends itself to variable precision calculation admirably. I use MAPM by Michael C. Ring for extended precision since it's freeware and fairly easy to use and quite reliable. It can also be slow, especially if you're not careful with rounding results. Anyway, if you can share anything more about your methods, I'd appreciate it. If not, I certainly understand. ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-17-2012, 06:15 AM Hi, Steve: Quote: Valentin, I'd be quite interested in learning any details you could share about the workings of your quadrature routine, as well as your extended precision setup. Thanks for your interest, Steve, and congratulations on your many and varied implementations of assorted quadrature algorithms, I see you're really keen on numerical integration. As for details on my own quadrature routine which I mentioned a few posts ago, I'm sorry but I'll keep them private for the time being, mainly for the following reasons: 1) my routine uses a simple algorithm, nothing fancy at all, and thus I don't think you'd learn anything worthwhile by having a look at it, in fact I think you'd be somewhat disappointed, a likely case of unfulfilled high expectations. 2) I wrote this routine some years ago for an article I was going to submit to Datafile for publication as part of my "Boldly going ..." regular series of sizable articles featuring novel techniques, as it was purpose-made to achieve a "bold" goal. After a month or so I eventually finished the article, which includes the routine proper with full comments and, as usual, fully-worked-out examples galore, with lots of relevant comparisons and tests, caveats, etc. Alas, things went south with Datafile, regrettably, and the 16-page article remains unpublished to this day so I'd rather not disclose its fine details till I eventually try and publish it somewhere else, together with about six or seven additional unpublished ones. A pity indeed because the finished article looks awesome, even if I say so myself. As a policy, I always write the kind of articles I'd love to read myself. Thanks for your understanding and Best regards from V. ``` ``` ▼ W. Bruce Maguire II Member Posts: 180 Threads: 23 Joined: Apr 2012 05-17-2012, 01:32 PM Valentin wrote: Quote: As for details on my own quadrature routine which I mentioned a few posts ago, I'm sorry but I'll keep them private for the time being ... I'd rather not disclose its fine details till I eventually try and publish it somewhere else... You have got to be kidding! What sort of reward do you envision for "publishing" your programmable-calculator code (outside of this forum)? And where might you "publish"? Do you really believe that there is a big, receptive audience outside of this forum? Bruce. Steve Perkins Member Posts: 104 Threads: 0 Joined: Dec 2007 05-17-2012, 04:27 PM Since you've gone to a lot of work, and hope to publish, I certainly respect your choice to keep it private for now. I was kind of hoping it might be appropriate for use in calculator circles, and it appears that this is the case! In the meantime, my interest in numerical integration has been re-kindled. I've already found ways to speed up my variable precision routine by a factor of 3, so that has been worthwhile. It's not so much finding test cases that a given routine can't handle that bothers me, but finding cases which the routine claims to handle with a reported error estimate, but is way off. My code handles cases with singularities (integrated from 0 to 1) like: 1/sqrt(x) log(x)*sqrt(x) without much problem. Also (integrating say from 0.01 to 1): sin(1/x)/(x*x) also fairly readily. I'm always looking for more interesting test cases as have been presented here. Hopefully we can soon have calculator routines at least as good as Casio's available to us on our HP machines. As always, I greatly appreciate your input to the forum. ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-17-2012, 05:07 PM Quote: I'm always looking for more interesting test cases as have been presented here. Have you tried Kahan's integral already? ``` /1 | I = | (sqrt(x)/(x-1) - 1/ln(x))dx /0 ``` It has the advantage of having a closed-form solution, as shown here: ```I = 2 - ln(4) - EulerGammaI = 0.036489973978576520559023667001244432806840395339565892952872746128345029282945897851326282715415875401 ``` An alternate form is ```I = Psi[3/2] ``` which can be checked at W|A. ▼ Mike Morrow Posting Freak Posts: 758 Threads: 9 Joined: Jul 2007 05-18-2012, 02:14 PM At the risk of boring folks with more fx-115ES Plus results, I report that Dr. Kahan's classic HP-34C test integral produces for once slightly mixed results. Adding 1E-12 offsets to both limits of integration to eliminate singularity at those points during the Gauss-Kronrod process, the 15-digit result after 121 seconds is: 0.0364899739857302 The actual displayed value is 0.03648997399, so only the first 9 significant digits of the result are correct...the last should have been 8. Still, this is a much better result in far less time than any HP calculator I've used. This includes the HP 50G, which has much faster hardware producing results for other problems about 20 times faster than does the 115ES Plus. This is a trivial-cost high school student's calculator. I can't help being impressed with the algorithm engineering at Casio, and a little disappointed by HP's in this area. It's way past time for Gauss-Kronrod quadrature to be an option in HP scientific calculator firmware. Edited: 18 May 2012, 2:27 p.m. ▼ Steve Perkins Member Posts: 104 Threads: 0 Joined: Dec 2007 05-18-2012, 04:30 PM I'm surprised you had to modify the end-points. My understanding is that Gauss-Kronrod does not need to evaluate the function at the extremes. That's part of it's usefulness. My implementation of Gauss-Kronrod also comes up with a reasonable approximation of the integral, but the error estimate isn't as good as I would hope. Overly optimistic I guess I might say. A good test case. Steve Perkins Member Posts: 104 Threads: 0 Joined: Dec 2007 05-18-2012, 04:19 PM My double exponential quadrature routine with variable precision handles that one quite nicely too. I incorporated it into my test suite when I read about it in this forum. Thanks for the suggestion. That function does have interesting behavior at certain spots. Kiyoshi Akima Senior Member Posts: 325 Threads: 18 Joined: Jul 2006 05-17-2012, 06:36 PM Another simple integral with a closed-form solution is the integral of FRAC(x) from 0 to 6.4 . Every HP handheld calculator to date with a built-in integrator produces a >50% relative error (or >100% depending on how you define relative error). The Casio fx-115ES Plus gets it exact. Gene Wright Posting Freak Posts: 1,545 Threads: 168 Joined: Jul 2005 05-17-2012, 07:47 PM Submit your articles to HP Solve. I know it is PDF, but it would reach a great deal more people than Datafile ever has. Be good if more people would learn about Valentin! ▼ Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-18-2012, 07:22 PM Quote: Submit your articles to HP Solve. I know it is PDF, but it would reach a great deal more people than Datafile ever has. Be good if more people would learn about Valentin! Thanks for your always kind words and your continued appreciation, Gene. The problem is, I already have the article written as a PDF document formatted per the mandatory Datafile template in use at the time and reformatting it per HP Solve guidelines would surely be a lot of work and thus a lot of additional time I simply don't have. Also, I have two versions of the routine, one of them is written for and runs in a non-HP multiprecision environment and thus is likely to be unsuitable for HP Solve. The other version is a standard-precision HP-71B program but I think that HP Solve will probably be much more interested in publishing materials for newer, current, easily available models than for hard-to-get deprecated old timers like the HP-71B. Anyway, thanks again for your interesting, well-meant advice and Best regards from V. ``` ``` Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-13-2012, 11:41 PM Quote: Well, it doesn't seem that tough at all Well, it is tough enough for all HP calculators I tried. For instance, on the ultra fast HP-15C LE it apparently will run forever when FIX mode is set to 8 or 9. Mike Morrow's Casio fx-115ES Plus below has 15 internal digits, but will get only the first 11 digits right. It is very fast and all 10 digits in the display will be correct, however. Best regards, Gerson. ▼ Mike Morrow Posting Freak Posts: 758 Threads: 9 Joined: Jul 2007 05-14-2012, 01:32 AM The interesting thing about the speed of the results from the fx-115ES Plus (aka 991ES Plus C) is that the equivalent of a 2500-iteration Savage benchmark takes 23 times longer on the Casio than on the HP 15C-LE. The magic in the Casio's quadrature result comes from the algorithm, not the hardware. Edited: 14 May 2012, 1:51 a.m. Mike Morrow Posting Freak Posts: 758 Threads: 9 Joined: Jul 2007 05-13-2012, 09:49 PM The Casio fx-115ES Plus uses a Gauss-Kronrod quadrature method and requires the function to be evaluated at the endpoints. Adding 1E-12 to each endpoint and setting a tolerance value of 1E-14 yields the following result after 63.0 seconds: 0.995555497761582 What shows up in the 10-digit display is 0.9955554978, which is accurate for every digit. That's not too shabby for an \$18 calculator. Edited: 13 May 2012, 10:17 p.m. ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-14-2012, 06:17 AM Quote: That's not too shabby for an \$18 calculator. Agreed! One of these days I had the band of my wristwatch replaced at a Casio representative and saw one of these on the shelf. Coincidence or not the price tag was exactly R\$ 115 (~US\$ 57.50), what made me lose the interest. Bart (UK) Posting Freak Posts: 850 Threads: 10 Joined: Mar 2009 05-14-2012, 09:59 AM www.numberempire.com gives: .9955554977623585 ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-14-2012, 05:41 PM Thanks for the link! Previously I had tried this link, which offers various integration methods, but I wasn't able to make it work with this example. I tried the WP 34S (v2) again and divided the integration interval into three subintervals, [0,0.1], [0.1,0.9] and [0,9,1], and obtained 0.9955554977624464. Only slightly closer though. P.S.: The following is better and will give 14 correct digits: ```WP 34S v2: 001 LBL 00 002 # pi 003 * 004 COS 005 x^2 006 ACOS 007 x<>y 008 # pi 009 * 010 SIN 011 / 012 RTN 0 ENTER .5 g Integral 00 2 * 3 / .5 + --> 0.9955554977623607 ``` Edited: 14 May 2012, 9:30 p.m. ▼ Bart (UK) Posting Freak Posts: 850 Threads: 10 Joined: Mar 2009 05-16-2012, 03:24 AM That is an interesting link, I had not come accross it. I like the choice of different methods too. Les Wright Posting Freak Posts: 1,368 Threads: 212 Joined: Dec 2006 05-14-2012, 10:46 PM Gerson, you seem to have a glitch in the code for your functions program. Once you consume your x with the cos(Pi x) calculation, etc, where are you getting a second x for the denominator? This is what I used: ```LBL'SS' # Pi * FILL COS STO* X ACOS SWAP SIN 3 * / . 5 + RTN ``` I use my own integrator on the 34S, which is a simple port of the PPC IG routine. My routine requires that only the last two extrapolated estimates agree, whereas the internal routine requires three and this means things take much longer. With my routine, on a WP34S with the latest firmware, I get 9.95555497791e-1 at SCI 11. This is off just one ULP in the 12th digit. The integral isn't that pathological and deserves another look with a routine that computes it correctly. Les ▼ Les Wright Posting Freak Posts: 1,368 Threads: 212 Joined: Dec 2006 05-14-2012, 10:59 PM PS: My batteries are getting sluggish so the WP34S built-in integration routine in the first surge and shuts the calc off and I am too lazy to dig out fresh ones at the moment. However, in the most recent emulator set to ALL 00 and DP, I not only get the first 12 digits exactly, but internally the routine gives us 6 more, to 18 in totally. And the emulator does so instantaneously, which means this will not take forever on the real calc. And considering that the calc's routine picks the accuracy according to display mode (ALL 00 is treated as SCI 11), getting 6 more digits is darned impressive. Les Edited: 14 May 2012, 11:00 p.m. Valentin Albillo Posting Freak Posts: 1,755 Threads: 112 Joined: Jan 2005 05-14-2012, 11:08 PM Quote: With my routine, on a WP34S with the latest firmware, I get 9.95555497791e-1 at SCI 11. This is off just one ULP in the 12th digit. Nope, it's 3 units off in the 11th digit, see my 1,000+ decimal places result above which begins with: ``` 0.995555497762367822986776561111752231600962645962... ``` Regards. V. ``` ``` ▼ Les Wright Posting Freak Posts: 1,368 Threads: 212 Joined: Dec 2006 05-15-2012, 12:12 AM Sorry, typo. 11th digit in mine should be a 6, not a 9. My claim is right. I transcribed my result incorrectly. Les Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-14-2012, 11:16 PM Hello Les, Quote: Once you consume your x with the cos(Pi x) calculation, etc, where are you getting a second x for the denominator? From the stack, which is automatic filled with the content of the X register every time the function to be integrated is called. From the manual: `"Lower and upper integration limits mustbe supplied in Y and X, respectively. Otherwise,the user interface is as in HP-15C. Please refer tothe HP-15C Owner’s Handbook, Section 14 andAppendix E, for more information about automaticintegration and some caveats."` Quote: With my routine, on a WP34S with the latest firmware, I get 9.95555497791e-1 at SCI 11. This is off just one ULP in the 12th digit. By simplifying the integrand prior to integration I was able to obtain 14 correct digits. Please see my reply to Bart above. Quote: The integral isn't that pathological and deserves another look with a routine that computes it correctly. I think you're right. The transitions at integer values of x appear to be fast enough to be a problem. Regards, Gerson. Alexander Oestert Senior Member Posts: 429 Threads: 31 Joined: Jul 2011 05-15-2012, 03:55 PM Please, someone give me a hint on how long this takes on the 34s and the 15c LE. It seems to take forever and I don't want to wear out my batteries. ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-15-2012, 04:34 PM Hello Alexander, On the WP 34S (v2) it takes only a couple of seconds. On the HP-15C you should set FIX mode to 7 at most. Alternatively use the following program: ```001 LBL A 002 pi 003 * 004 COS 005 x^2 006 ACOS 007 x<>y 008 pi 009 * 010 SIN 011 / 012 RTN 0 ENTER .5 g Integral A 2 * 3 / .5 + f --> 0.995555504 ``` It will take less than 3 seconds on the HP-15C LE, even when FIX mode is set to 9. Cheers, Gerson. ▼ Alexander Oestert Senior Member Posts: 429 Threads: 31 Joined: Jul 2011 05-15-2012, 04:52 PM Duh! Switching to RAD did it on the WP-34s V2 and the 15c LE, but on the WP-34s V3 2931 I can't get it to work. EDIT: Solved on the WP-34s V3 by explicitly setting an END at the end instead of terminating only by RTN. I didn't know that. Edited: 16 May 2012, 2:10 a.m. ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 05-17-2012, 01:20 AM This sounds strange. Can you send me an email with details so that I can have a closer look when I find the time? ▼ Alexander Oestert Senior Member Posts: 429 Threads: 31 Joined: Jul 2011 05-18-2012, 02:02 AM I can't remember what state the calculator was in when I encountered this, so I took a recently converted machine (VERS 3049, only calc.bin) and did a RESET. Set it to [g][RAD] Enter the WP-34s program from posting #1, exactly the 17 steps, no less, no more 0 [ENTER] 1 [g][2] 0 0 For a short time I see 'Running Program' then the integral solving symbol appears and the 9.95555497761E-1 converges. So far so good. Because the RCL annunciator is still flashing, I [EXIT]. I key in again 0 [ENTER] 1 [g][2] 0 0 Now I see 'Running Program' and the RCL flashing, no integral symbol. After about 2 minutes RCL flashing is off, the matrix display shows 'Reset' and 9.96943604201E-1 is in the lower display. Another 0 [ENTER] 1 [g][2] 0 0 has the machine running for 5 minutes before I forcibly [EXIT] with 'Stopped' and 4.72450120691E-1 in the display. Now I do [GTO][.][.] and enter the following program 001 LBL A 002 [g][x^2] 003 5 004 - 005 [g][RTN] Then [f][SLV] brings up -2.2360679775 in the display. Again 0 [ENTER] 1 [g][2] 0 0 does a fast conversion, with RCL still flashing when I [EXIT]. The next 0 [ENTER] 1 [g][2] 0 0 has to be aborted after 5 minutes with [EXIT]. [SLV]ing program A from above results in 'Running Program'. I end the experiment here. Marcus, I send you this as an email but post this here so everybody can have fun... ;-) ▼ Paul Dale Posting Freak Posts: 3,229 Threads: 42 Joined: Jul 2006 05-18-2012, 02:19 AM I'd be guessing that the return stack isn't being cleared each time and that is confusing things. Try a g-shirt RTN in run mode to clear things. Still, it seems like a bug/annoyance. - Pauli Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 05-18-2012, 02:49 AM Pauli is right, you need to do do g-RTN before you can try again. What happens is that you are effectively nesting calls to he integrator. EXIT stops the process while in your user code routine to integrate. You can R/S to continue. If, instead, you start a new integration operation this will nest. If you stop it again, it will eat more and more of the available memory. There is still something wrong. Showing "Running Program" is not intended. "Reset" is most probably the watch dog kicking in. ▼ Alexander Oestert Senior Member Posts: 429 Threads: 31 Joined: Jul 2011 05-18-2012, 03:34 AM [g][RTN] does it. I just wanted to mention, if I enter program mode after an integral run and without [g][RTN], memory is only about 415 steps after the first try and 320 steps after the second try, 225 after 3rd... ▼ Paul Dale Posting Freak Posts: 3,229 Threads: 42 Joined: Jul 2006 05-18-2012, 03:46 AM About right. The integration code allocates a lot of local registers. We probably should clear the locals/return stack when solve or integrate is run from the keyboard.... - Pauli ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 05-18-2012, 01:58 PM Quote: We probably should clear the locals/return stack when solve or integrate is run from the keyboard.... That's what I came up with on the last kilometers of my latest bike trip. Let me see when I find the time or a clean implementation. The "Running Program" is indeed intentional because the software detects that it is not run from the top level and is therefore quiet (no display of intermediate results). This leaves the original "Running" message in the display. I'm still thinking about the "Reset" case. I have an idea how to fix this (if it's the watchdog) but I must check if I'm right. ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 05-19-2012, 03:21 AM I did something about the issue. Please report any findings! fhub Posting Freak Posts: 1,216 Threads: 75 Joined: Jun 2011 05-18-2012, 04:06 AM Quote: I just wanted to mention, if I enter program mode after an integral run and without [g][RTN], memory is only about 415 steps after the first try and 320 steps after the second try, 225 after 3rd... I guess this happens only if you abort the integral calculation with [EXIT] each time!? I've now tried a few integrals with the emulator (last build) but I don't get this problem. Of course I couldn't interrupt the integration because the emulator is just too fast. Pauli has already posted the solution: clearing the return stack (and local regs) for a manual call should solve the problem. This should also be done for SUM and PROD (not only INTG and SLV), I don't know if there are other such routines which call a user program. Franz Edited: 18 May 2012, 4:07 a.m. Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-19-2012, 09:49 PM Quote: It seems to take forever and I don't want to wear out my batteries. Same here with versions 3.0 2609 and 3.1 2988. Has it been fixed? Thanks! ▼ Marcus von Cube, Germany Posting Freak Posts: 3,283 Threads: 104 Joined: Jul 2005 05-20-2012, 02:25 AM The adaptive integrator in WP 34S V3 takes a few seconds even on the emulator so be prepared for some delay. I have no idea what the END/RTN issue is. ▼ Gerson W. Barbosa Posting Freak Posts: 2,761 Threads: 100 Joined: Jul 2005 05-20-2012, 05:59 PM After 25 minutes it is showing 9.95555498002-1 for five or six minutes now, diverging from a previous closer result. A useful feature here would be the possibility of interrupting the calculation whilst keeping the latest approximation on the display instead of the not meaningful result I get when I eventually abort it, 1.54227330201 (v3.1 2988).

 Possibly Related Threads... Thread Author Replies Views Last Post Integral Richard Berler 32 3,552 09-22-2013, 08:59 PM Last Post: Les Koller That's a tough nut to crack... Matt Agajanian 7 1,121 05-25-2012, 01:53 PM Last Post: Mike Morrow can't calculate integral on hp 41 Hans Holzach 12 1,653 05-16-2012, 02:59 PM Last Post: Hans Holzach Just how tough is a 32S-II ? Matt Agajanian 12 1,517 03-11-2012, 05:36 AM Last Post: Thomas Radtke Polynomial Derivative and Integral for the HP 35S Eddie W. Shore 0 391 05-24-2011, 10:29 PM Last Post: Eddie W. Shore Indefinite integral calculation with 49g+/ 50g ferroburak 7 1,034 11-10-2009, 02:27 PM Last Post: David Hayden Solving for an Integral: 28S solution Marcus von Cube, Germany 2 567 03-02-2008, 10:49 AM Last Post: Marcus von Cube, Germany Handheld Calculator Helps Solving Kahan's Integral Symbolically Gerson W. Barbosa 3 641 02-26-2007, 11:28 AM Last Post: Gerson W. Barbosa HP integral Sara Linn 4 692 01-23-2007, 05:26 PM Last Post: dbatiz Kahan's "weird" integral revisited Les Wright 1 467 01-02-2007, 01:29 PM Last Post: Namir

Forum Jump: