Test suite failures with cvs version?



Hello Ray,

It seems to me that I am a bit tired. Again I was too fast. I did the check for
conjugate values with different algorithm and not with conjugate values. 

Now I arrived at this routine:

(defun expintegral-ei (z)
  (+
    (- (expintegral-e 1 (- z)))
;    (- (* 0.5 (- (log z) (log (/ 1 z)))) (log (- z)))
    (cond 
      ((> (imagpart z) 0)
       ;; Positive imaginary part. Add phase %i*%pi.
       (complex 0 (float pi)))
      ((< (imagpart z) 0)
       ;; Negative imaginary part. Add phase -%i*%pi.
       (complex 0 (- (float pi))))
      ((> (realpart z) 0)
       ;; Positive real value. Add phase -%i*pi.
       (complex 0 (- (float pi))))
       ;; Negative real value. No phase factor.
      (t 0))))

Again the results I obtained with GCL:

(%i9) expintegral_ei(1.0);
(%o9) 1.895117816355937
(%i10) expintegral_ei(-1.0);
(%o10) -0.21938393439552
(%i11) expintegral_ei(1.0 + %i);
(%o11) 2.387769851510522*%i+1.764625985563854
(%i12) expintegral_ei(1.0 - %i);
(%o12) 1.764625985563854-2.387769851510522*%i
(%i13) expintegral_ei(-1.0 + %i);
(%o13) 2.962268118550434*%i-2.8162445198177954E-4
(%i14) expintegral_ei(-1.0 - %i);
(%o14) -2.962268118550434*%i-2.8162445198177954E-4

This results including the conjugate values are correct for the expintegral_ei.
The testsuite includes some further tests for conjugate values which know work
as expected.

Dieter Kaiser

-----Urspr?ngliche Nachricht-----
Von: raymond.toy at ericsson.com [mailto:raymond.toy at ericsson.com] 
Gesendet: Montag, 27. Oktober 2008 21:06
An: Dieter Kaiser
Cc: maxima at math.utexas.edu
Betreff: Re: [Maxima] Test suite failures with cvs version?

>>>>> "Raymond" == Raymond Toy <raymond.toy at ericsson.com> writes:

    Raymond> Dieter Kaiser wrote:
    >> 
    >> (defun expintegral-ei (z)
    >> (+
    >> (- (expintegral-e 1 (- z)))
    >> ;    (- (* 0.5 (- (log z) (log (/ 1 z)))) (log (- z)))
    >> (if (and (= (imagpart z) 0) (> (realpart z) 0))
    >> ;; Positive real value. Add a phase factor -%pi*%i
    >> (complex 0 (- (float pi)))
    >> ;; For all other values. No phase factor.
    >> 0)))

    Raymond> I was getting ready to check in all the changes we had done last
week,
    Raymond> but I was looking at this expression.  I must be stupid.  I don't
see
    Raymond> how  f(z) = 1/2*(log(z)-log(1/z)) - log(-z) simplifies that way.
For
    Raymond> real z, yes, I see f(z) is -%i*%pi for positive z and 0 for
negative z.
    Raymond>  But for complex z, not on the real axis, I always get %i*%pi (or
    Raymond> -%i*%pi, depending on how you want to write -1).

As best as I can tell f(z) = -%i*%pi for all z except the negative
real axis where it is 0.

However, if we use this version, we get weird results, assuming
expintegral_ei is continuous in the right half plane:

expintegral_ei(1.0)         -> 1.895 - 6.283185307179586*%i
expintegral_ei(1.0+1d-8*%i) -> 1.895 - 3.14159*%i
expintegral_ei(1.0-1d-8*%i) -> 1.895 + 3.14159*%i

If we use the actual expression f(z), we get something that seems
reasonable:

expintegral_ei(1.0)         -> 1.895
expintegral_ei(1.0+1d-8*%i) -> 1.895 - 2.71828d-8*%i
expintegral_ei(1.0-1d-8*%i) -> 1.895 + 2.71828d-8*%i

But this depends on expintegral-e returning the correct complex
values, and I'm not quite sure about that yet.

Ray