search

maroba / findiff

Python package for numerical derivatives and partial differential equations in any number of dimensions.
MIT License
413 stars 59 forks source link

Analytic inversion of Vandermonde matrix #66

Closed dwierichs closed 1 year ago

dwierichs commented 1 year ago

Hi, new fan of this package! :)

I was wondering: At least within the applications I'm looking at, the equation solving performed in L99 of coefs.py (or L89 for symbolic=True) can also be done analytically, because matrix is a Vandermonde matrix, and computing specific elements of its inverse is doable without any matrix inversion (see for example ProofWiki). Would this be interesting to implement? I noticed that findiff.coefficients becomes unstable at some point (admittedly for quite large values of acc~25).

Thanks for the good work, kind regards!

maroba commented 1 year ago

Hi @dwierichs ,

thanks! Glad that you like it. :-)

In principle, you are right that one could solve the equation analytically. However, I guess the price would be less readability. Since this part of the code is not critical in terms of performance (the equation system is always small for a computer), I would give priority to readability and maintainability.

Regarding the other point you mentioned about stability, you are absolutely right, and this should be fixed. I think the easiest way to do it, would be to use symbolic computation when high accuracy is requested and convert it back to floating point numbers at the end.

Kind regards, Matthias

maroba commented 1 year ago

@dwierichs Now that I'm thinking about it, using the direct solution without solving the system, as you proposed initially, may also solve the stability issue. That would make the readability argument less important. So, I think it's worth a try!

dwierichs commented 1 year ago

Hey @maroba, thanks for the quick reply! Yeah I definitely agree that the equation solver is more readable and maintainable :) Indeed my thought was that the direct computation could be more stable, but if the trick of doing stuff symbolically and then converting back is stable as well, that might actually be the nicer path (even if not optimal computationally, but as you say, that doesn't really matter).

maroba commented 1 year ago

@dwierichs I'm working on some other topic right now. So if you would like to implement it yourself, feel free to do so. If you need help finding orientation in the code, just tell me.

maroba commented 1 year ago

Closed with PR #67