scott-maddox / fdint

Precise and fast Fermi-Dirac integrals of integer and half integer order.
BSD 3-Clause "New" or "Revised" License
17 stars 7 forks source link

for full numpy compatibility: vectorized functions #16

Open marklundeberg opened 8 years ago

marklundeberg commented 8 years ago

Thanks for a great plugin! This is very useful.

For ultimate smoothness with the numpy philosphy, it would be great if the fd* functions would operate elementwise over the input arrays including broadcasting, much in the way of functions in the scipy.special module.

For example, this passes a column 2-vector and row 2-vector to the order-n Bessel function, resulting in a matrix of J_5(1), J_5(2), J_6(1), J_6(2)::

>>> from scipy.special import jn
>>> jn(array([[5],[6]]),array([1,2]))
array([[  2.49757730e-04,   7.03962976e-03],
       [  2.09383380e-05,   1.20242897e-03]])

Currently, this is not directly implemented for the fdint functions, however it can be done by wrapping them with numpy's "vectorize" class.

>>> fdk_vec = numpy.vectorize(fdint.fdk)
>>> fdk_vec(array([[5],[6]]),array([1,2]))
array([[  314.66805418,   817.16312139],
       [ 1920.6214911 ,  5088.14241086]])

This is unfortunately a bit inefficient in that it essentially performs a for loop in python, rather than a for loop in C.

scott-maddox commented 8 years ago

Thanks for the issue submission, I'm glad you find the package useful!

It would definitely be good more closely match scipy's standards, especially if it can be done without reducing performance, or if performance can be maintained by calling a separate function.

I'm no longer actively using this package, so I can't guarantee a speedy followthrough, but I would be happy to accept a pull request if you (or anyone) would like to take it on themselves

marklundeberg commented 8 years ago

Totally understandable, for now I'll use the vectorize wrapper but if I see performance issues you may see a pull request from my end.

By the way, another thought: in case you don't plan on maintaining the package long-term, it might make sense to integrate it into scipy.special, which lacks F-D integrals and the related polylogarithm. You certainly seem to have implemented a very fast code and I think many python users may be unaware that your package exists. If you're supportive of this idea I can get a feature request going and maybe even implement it myself (would be a good excuse to learn the innards of C+Python+cython). Cheers-

scott-maddox commented 8 years ago

Sounds great on both counts. I'm happy to support your efforts in any way I can.