A colleague came to my office and asked if IDL's implementation of QROMO produces accurate results since the example given by IDL is off from the analytic result by 1%. When you look at the example on IDL's QROMO manual page, it is to integrate the function 1/x^4 from 2 to infinity. Now, it is true that when you do exactly what is stated on the manual page, you get the stated results: PRINT, QROMO('hyper', 2.0, /MIDEXP) does yield 0.0412050 Anyone can integrate this by hand and one obtains 1/24 = 0.041666666.... Which differs in the fourth decimal place from the example in the IDL manual (~4 parts in 400, or 1%). IDL is certainly doing what is asked. But reading Numerical Recipes about using midexp says to only use midexp when the function decreases exponentially rapidly. Which the function in the example in the IDL manual pages isn't doing. So, this naive application in this IDL example is a poor example if one wants to show a valid result, and not only example usage. It is better to integrate this function using /midinf - and one can examine what the upper limit should be. Each of the three below is a much more accurate integration than the one listed in the example: IDL> print,qromo('hyper',2.0,5000.,/midinf) 0.0416667 IDL> print,qromo('hyper',2.0,500.,/midinf) 0.0416667 IDL> print,qromo('hyper',2.0,50.,/midinf) 0.0416640 If possible, I would suggest modifying the example on the QROMO manual page. -Marcos Montes
|