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
manualintegrate nested pow #23716
manualintegrate nested pow #23716
Conversation
✅ Hi, I am the SymPy bot (v167). 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.
Update The release notes on the wiki have been updated. |
Benchmark results from GitHub Actions Lower numbers are good, higher numbers are bad. A ratio less than 1 Significantly changed benchmark results (PR vs master) Significantly changed benchmark results (master vs previous release) before after ratio
[41d90958] [6afa478c]
<sympy-1.11.1^0>
- 966±3μs 620±2μs 0.64 solve.TimeSparseSystem.time_linear_eq_to_matrix(10)
- 2.80±0.01ms 1.16±0ms 0.41 solve.TimeSparseSystem.time_linear_eq_to_matrix(20)
- 5.67±0.01ms 1.70±0ms 0.30 solve.TimeSparseSystem.time_linear_eq_to_matrix(30)
Full benchmark results can be found as artifacts in GitHub Actions |
@oscarbenjamin is this good to go? |
8fcb040
to
7065541
Compare
7065541
to
942629f
Compare
if e.has_free(x): | ||
raise ValueError() | ||
base_, sub_exp = _get_base_exp(b) | ||
return base_, sub_exp * e |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Does this assume that (a**b**c) == a**(b*c)
? If so then that identity does not always hold.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Basic power rule: to calculate integral of x**n
, n != -1
, we multiply it by x
and divide it by n+1
, so we get x**(n+1)/(n+1)
.
Basic log rule: to calculate integral of x**-1
, we multiply it by x*log(x)
, so we get log(x)
.
Now let's extend these basic rules to more complicated forms like (x**2)**(-1/2)
and x*cbrt(1/x)
.
I know (x**2)**(-1/2) != x**-1
, but to calculate the integral, we should decide whether to use power rule or log rule. The decision is made by treating (x**2)**(-1/2)
as x**(2*(-1/2))
, so the "nominal" exponent is -1, thus we apply log rule: multiply it by x*log(x)
to get x*log(x)/sqrt(x**2)
. The correctness is easily verified: diff(x*log(x)/sqrt(x**2)) == 1/sqrt(x**2)
.
Similar thing happens to x*cbrt(1/x)
as well. The "nominal" exponent is 1+(1/3)*(-1) == 2/3
which is not -1, so we apply power rule: multiply it by x
and divide it by 5/3
. The result is 3*x**2*cbrt(1/x)/5
, and diff(3*x**2*cbrt(1/x)/5) == x*(1/x)**(1/3)
.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Oh, I see. So we just use sub_exp*e
in the divisor.
try: | ||
base, exp_ = _get_base_exp(integrand) | ||
except ValueError: | ||
return |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It is generally better not to raise and catch ValueError
as part of normal flow control because various bugs can result in raising ValueError
and then this would obscure those bugs by silencing the exception. Here a special purpose exception should be used like class NoMatch(Exception)
.
I think this looks good. |
I'm not sure that this would work if the base is not positive. |
This aims to solve negative base issue. See #23712 (comment) which is wrong for negative base. |
References to other Issues or PRs
Fixes #23712
Fixes #23505
Brief description of what is fixed or changed
This is an extension to power rule.
Nested pow like
sqrt(x**(S(5)/3))
can be treated asx**(S(5)/6)
, but the structure should be kept in the end.However, if the exponent is equivalent to -1, the integrand should be multiplied by
x*log(x)
instead ofx/(exponent+1)
Other comments
Release Notes