Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add Thiele's interpolation formula #24109

Open
wants to merge 12 commits into
base: master
Choose a base branch
from
Open

Conversation

baruchel
Copy link
Contributor

@baruchel baruchel commented Oct 4, 2022

Brief description of what is fixed or changed

Added a new interpolation formula

Release Notes

  • functions
    • Added a new function for Thiele's interpolation formula

@sympy-bot
Copy link

sympy-bot commented Oct 4, 2022

Hi, I am the SymPy bot (v169). I'm here to help you write a release notes entry. Please read the guide on how to write release notes.

Your release notes are in good order.

Here is what the release notes will look like:

This will be added to https://github.com/sympy/sympy/wiki/Release-Notes-for-1.12.

Click here to see the pull request description that was parsed.
#### Brief description of what is fixed or changed
Added a new interpolation formula

#### Release Notes

<!-- BEGIN RELEASE NOTES -->

* functions
  * Added a new function for Thiele's interpolation formula

<!-- END RELEASE NOTES -->

@baruchel baruchel closed this Oct 4, 2022
@baruchel baruchel reopened this Oct 4, 2022
>>> f.subs(x, 0.5) * 3
3.14159265358979

>>> thiele([1, 2, 3, 4], [10, 12, 14, 16])
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Would it be better to raise an exception in this case?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, but how can I achieve that? Not sure about the cleanest wai...

Copy link
Contributor Author

@baruchel baruchel Oct 5, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@asmeurer I tried to implement it. Did I do it the right way with expr.has(nan, zoo)?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@asmeurer Could you review my code and tell me if it is OK?

@asmeurer
Copy link
Member

asmeurer commented Oct 4, 2022

Please add tests for this.

@oscargus oscargus added the polys label Oct 5, 2022
@baruchel
Copy link
Contributor Author

baruchel commented Oct 6, 2022

Please add tests for this.

What do you mean? Tests for the error cases? How could I do?

@bjodah
Copy link
Member

bjodah commented Oct 10, 2022

@baruchel take a look at a few of the other (merged) pull requests, we require that all new functions are covered by unit tests.

@oscarbenjamin oscarbenjamin marked this pull request as draft November 14, 2022 20:42
@oscarbenjamin
Copy link
Contributor

I'm going to mark this as draft. Feel free to set it as ready for review when some tests have been added.

@baruchel baruchel marked this pull request as ready for review January 25, 2023 17:53
@baruchel
Copy link
Contributor Author

@oscarbenjamin Sorry for the long delay, I added some more tests. Is it OK. Some checks were not successful above, but I can't figure out if it comes from me or from something else in the checking process.

@oscarbenjamin
Copy link
Contributor

I'm not sure what caused the doctests to fail. I've restarted them.

"""
n = len(u)
assert len(v) == n
seterr(divide=True)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh actually this is the problem.

@github-actions
Copy link

github-actions bot commented Jan 25, 2023

Benchmark results from GitHub Actions

Lower numbers are good, higher numbers are bad. A ratio less than 1
means a speed up and greater than 1 means a slowdown. Green lines
beginning with + are slowdowns (the PR is slower then master or
master is slower than the previous release). Red lines beginning
with - are speedups.

Significantly changed benchmark results (PR vs master)

Significantly changed benchmark results (master vs previous release)

       before           after         ratio
     [41d90958]       [857a97e6]
     <sympy-1.11.1^0>                 
-     1.11±0.02ms          695±3μs     0.62  solve.TimeSparseSystem.time_linear_eq_to_matrix(10)
-      3.24±0.1ms      1.30±0.02ms     0.40  solve.TimeSparseSystem.time_linear_eq_to_matrix(20)
-     6.55±0.06ms      1.87±0.04ms     0.29  solve.TimeSparseSystem.time_linear_eq_to_matrix(30)

Full benchmark results can be found as artifacts in GitHub Actions
(click on checks at the top of the PR).

@baruchel
Copy link
Contributor Author

Hi, @oscarbenjamin Is it OK now? Regards.

@bjodah
Copy link
Member

bjodah commented Jan 27, 2023

@baruchel I would expect to find a method named "test_thiele_interpolate" here:
https://github.com/sympy/sympy/blob/master/sympy/polys/tests/test_rationaltools.py

@baruchel
Copy link
Contributor Author

Hi @oscarbenjamin I added some tests.

sympy/polys/rationaltools.py Outdated Show resolved Hide resolved
sympy/polys/rationaltools.py Outdated Show resolved Hide resolved
sympy/polys/rationaltools.py Outdated Show resolved Hide resolved
sympy/polys/rationaltools.py Outdated Show resolved Hide resolved
sympy/polys/rationaltools.py Outdated Show resolved Hide resolved
sympy/polys/rationaltools.py Outdated Show resolved Hide resolved
and their function values in vector v (both vectors must be of equal
lengths).

This algorithm might encounter division by zero
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is an example of this (or reference)?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A ZeroDivisionError is yielded whenever values are equally spaced. For instance,
thiele([0,1,2,3], [0,1,2,3])

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would like to see a better criteria stated. If you interpolate on ([1+z,2,3,4],(1,2,3,5)) and then replace z with 0 you get x which is clearly not right. If you use only 3 points as in (1,2,3),(1,2,3) it doesn't fail. Maybe this should be limited to numeric input whose output is unambiguously ok or is an error indicating problems.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is the criteria that there cannot be 3 colinear points? There is a routine for checking colinearity in geometry.


This algorithm might encounter division by zero
when values are equally spaced. For such cases, a better algorithm is
sympy.polys.polyfuncs.rational_interpolate.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is it always true that the rational_interpolate(dat, n) == thiele_interpolate(*dat) if n is the degree of the Thiele numerator?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

7 participants