[SciPy-Dev] Find points in delaunay triangulation : scipy.spatial vs. scipy.interpolation

Sun May 13 13:32:32 EDT 2012

Le 13/05/2012 13:03, josef.pktd at gmail.com a écrit :
> On Sun, May 13, 2012 at 12:12 PM, Pablo Winant<pablo.winant at gmail.com>  wrote:
>> Le 13/05/2012 11:40, josef.pktd at gmail.com a écrit :
>>
>> On Sun, May 13, 2012 at 10:57 AM, Pablo Winant<pablo.winant at gmail.com>
>> wrote:
>>
>> Le 13/05/2012 08:44, josef.pktd at gmail.com a écrit :
>>
>> On Sat, May 12, 2012 at 8:13 PM, Pablo Winant<pablo.winant at gmail.com>
>>   wrote:
>>
>> Hi,
>>
>> I tried to use interpolation routines in scipy recently and I have found
>> two slight performance issues
>>
>>    - The LinearNDInterpolation object implemented in cython requires a
>> list of points and a list of values to be created. But is is not
>> documented how to change the values of the interpolator without doing
>> the mesh again. This is useful when one is solving the values of a
>> function at the vertices of the mesh : one doesn't want to do the
>> triangulation again and again. Maybe there could be a simple specific
>> method to set the values in this case. In that case it would consist in
>> changing the value of a property but it would be consistent with more
>> general interpolation schemes.
>>
>> - I tried to use the delaunay object from scipy and noticed a strange
>> thing: for a given set of coordinates it takes longer to get the indices
>> of the triangles containing the points than it takes to perform the
>> interpolation using LinearND object. This is puzzling since apparently
>> the implementation of LinearND performs many calls to the qhull library
>> to get this indices. Attached is a simple exampe demonstrating this anomaly.
>>
>> One last thing: I have written an interpolation object on sparse grids,
>> using smolyak product of chebychev polynomials. It is written in pure
>> python (vectorized) and licensed under the bsd license. Currently it
>> lives in another library but I guess it would make more sense to have
>> something like that in a more general scientific lib. Let me know if you
>> are interested. (it is available there anyway:
>> https://github.com/albop/dynare-python/tree/master/dolo/src/dolo/numeric:
>> chebychev.py
>> and smolyak.py)
>>
>> Pablo,
>>
>> what is actually your license for dolo?
>>
>> It used to be GPL, but I changed it to BSD a while ago, precisely in
>> order to integrate better with the python community.
>> As I was using google-code at that time I switched to the new-bsd
>> license which was the only bsd option.
>> If I understand well, it is refered to as BSD-3.
>>
>> Now I must say I am a bit lost in this jungle, so I guess I would follow
>> your suggestion if you say BSD-n is better.
>>
>> A while ago I got confused about all these qualifiers for BSD (for
>> statsmodels), BSD with number is less ambiguous or easier to remember.
>> I don't think which BSD or which MIT doesn't matter as long as it is
>> clearly stated.
>>
>> There several other parts of the library which would fit better outside
>> (such as a nonlinear solver with complementarity constraints) so it is
>> important to me that the license makes a sensible relocation of the code
>> possible.
>>
>> your license file is GPL
>>
>> Thank you for spotting that. I will need to properly do all these legal
>> stuff once I am sure about the good license.
>>
>> https://github.com/albop/dynare-python/blob/master/dolo/LICENSE
>> but setup.py says BSD and and I found another package using some of
>> your code as BSD-2
>>
>> Can you tell me which one it is ?
>>
>> https://github.com/christophe-gouel/RECS/blob/master/LICENSE.txt
>> I only looked briefly, he uses your code to parse the model definition
>> files, otherwise matlab.
>> see bottom of this page https://github.com/christophe-gouel/RECS
>>
>> Ah, I know this one: he is a friend. He has written a software to solve
>> rational expectation models in matlab. In addition to dolo, it uses the
>> compecon toolbox which has some interesting interpolation routines
>> (everything is documented in a book "Applied Computational Economics and
>> Finance, Mario J. Miranda&  Paul L. Fackler, MIT Press" ). They have three
>> kind of one-dimensional interpolation routines (linear, splines, chebychev)
>> and they provide a flexible way to produce multidimensional interpolation as
>> a product of these one-dimensional routines.
> I looked at Miranda Fackler's code a few years ago, and thought they
> have a restrictive license. Now, I cannot find anything except the
> standard disclaimer on liability.
I've been wondering many time too and Christophe, the author of 
compecon, has asked one of the authors. Apparently they were not really 
aware of the various opensource license and simply meant the code to be 
"usable by anyone" but wanted to remain in control of their library. 
However, I guess they would have no problem with letting python code of 
part of the code be released under the BSD license, provided they are 
cited and the derived work si clearly distinct from their.

>
> I haven't quite figured out (no code) yet how to go in general from 1d
> to nd. (I was and am more interested in smoothing noisy versions, then
> pure interpolation.)
You can basically do a tensor product of the functions in each base. ( 
f(x,y) = (a11 f1(x) + a12 f2(x)) * (a21 f1(y) + a22 f2(y)) )
 From a numerical point of view, I think it is dominated by direct nd 
representation but it is relatively easy to implement and fast. For 
instance, multilinear interpolation where (f1 and f2 are linear) is very 
fast.


>>
>>
>> I was just browsing after your link. There is some interesting code
>> and it would be good if we can share some of it.
>>
>> I would be very glad to do so. Until now, I have been very busy doing
>> research (thesis, post-doc,...) but I am now trying to turn the code from
>> dolo into something useful for others.
>>
>> I also didn't know about a python - octave bridge.
>> and was trying to see how easily we could create qnwnorm with scipy.
>>
>> qnwnorm comes from the compecon toolbox. It is based on a very well known
>> algorithm for the one-dimensional quadrature and it should be fairly easy to
>> rewrite/port it to python. I can do that if you are interested (it is also
>> something I tried a while ago).
> scipy provides the tools for 1d quadrature (points, weights), but I
> didn't know that it's so easy to extend to multivariate integration.
> Do you know if this extends to quadrature with respect to other
> distributions than normal.
Frankly, no.

>
> Josef
>
>>
>> Cheers,
>>
>> Josef
>>
>> Best,
>>
>> Pablo
>>
>> Thanks,
>>
>> Josef
>>
>> Best regards,
>>
>> Pablo
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-dev
>>
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-dev
>>
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-dev
>>
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-dev
>>
>>
>>
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-dev
>>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-dev