Re: Complex Exponential Integral

On Tuesday, 8 July 2025 14:21:01 PDT Erik Luijten wrote:
> Hi Ethan,
> 
> In amos_airy.c, you write:
> 
>  * The cexint_() routine is also from the AMOS collection, but is not part
>  * of the Bessel function subset and not included in libopenspecfun.
>  * I could not find source for a version of the source that splits the
>  * real and imaginary parts of the arguments as with the zxxxx.f Bessel
>  * routines; this one accepts and returns a CMPLX argument rather than
> separate
>  * real and imaginary parts.
>  *             D. E. Amos, ALGORITHM 683  A Portable FORTRAN Subroutine
>  *             for Exponential Integrals of a Complex Argument
>  *             ACM Trans. Math. Software 16:178-182 (1990)
> 
> This confuses me for two reasons:
> - You invoke CEXINT as a double precision function, even though this
> Fortran function is single precision.

> - The algorithm 683 you quote DOES provide a separate double precision
> function, ZEXINT. Moreover, that function does exactly what you were
> looking for, namely taking and returning separate real and imaginary parts.

There is some history behind this due to IMHO over-zealous, or at least
overly pessimistic, concern about license contamination with non-free code.
I cited a publication in ACM/TMS that describes the implementation.
However, people have assserted that code obtained through the ACM
journal is non-free because the journal claims copyright.
I believe this to be factually incorrect; the published article may be subject
to copyright but not the code it describes. In any case, rather than pursue
the argument I instead used copies of the code obtained from US national
lab repositories.  I found source there for cexint.f but not for zexint.f,
so that is what I use myself to build libamos and what I invoked from the
gnuplot code.

Even so that was apparently not sufficient by itself to convince hard-core
"free licenses only" gatekeepers.
There is an explicit statement of the [non]copyright status of the code
in the cited 1985 Sandia National Laboratory report, but that report
predated development of the expint routines (1990?) and I did not find a
separate Sandia report that mentions them specifically. 
That is probably why neither routine is included in libopenspec,
and why inclusion in gnuplot is a configuration option.

> I changed f_amos_cexint to use ZEXINT (just for my own edification) and it
> works properly. However, before submitting it, I wonder what am I missing?

Where did you obtain a copy of the source?
Can you document the relevant claim or non-claim of copyright or
licensing?

And, out of curiosity, is there anything in the zexint.f code that would
differentiate it from a single-precision version other than declaration
of the variables as Fortran DOUBLE PRECISION?

	Ethan


> 
> Thank you,
> 
> Erik
>