|Extended Complex Number Calculations in SysRPL|
Message #1 Posted by Les Wright on 25 June 2007, 7:50 p.m.
This follows up a back and forth I was having with James Prange in an earlier thread.
I have completed my SysRPL code to compute the sine and cosine integrals by complex continued fraction. The original Modified Lentz algorithm and formulae can be found in Numerical Recipes 2nd ed or later. Internally, the complex numbers are treated as two part entities, where each part is an extended real. The parts can be either the real or imaginary parts of rectangular form, or the modulus and argument of polar form. I have used NULLLAMs and stack manipulation. The complex arithmetic really is normal extended real number arithmetic on the constituent parts--I work in rectangular form for additions and subtractions, and polar form for multiplication, division, roots, and powers. The SysRPL functions %%R>P and %%P>R allow ready switch back and forth without getting bogged down in the trigonometry.
I have managed to keep track of everything and it works fine, but this is one case where the SysRPL code, albeit not yet optimized is perceptibly slower than my UserRPL version. However, the trade of can be construed as worth it--since the internal calculations are done to 15 digits, the 12 digit output is completely accurate provided that the input argument is of suitable size.
This leads me to a practical question, that I believe someone asked here earlier. When working in extended reals in SysRPL, what would you chose as a "machine epsilon" for the purpose of testing convergence? In double precision on my PC, I think machine epsilon is about 2.2e-16. In my program, the convergence test is when a particular multiplicative factor becomes unity and doesn't change in a subsequent iteration. In UserRPL, converging to actual equality works fine, but in extended reals in SysRPL, this expectation drags the computation out unnecessarily for some input. I have by trial and error settled on 2e-15 as my EPS. Any smaller and you might as well just converge to equality unnecessarily. So if variable del is the factor I want to be unity, the loop exits when LAM del %%1 %%- %%ABS %%2E-15 %%< returns True. Requiring LAM del %%1 EQUAL to be true consumes cpu time and offers no added benefit.
Any thoughts on this?