-
-
Notifications
You must be signed in to change notification settings - Fork 30
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
Added Gauss-Turán quadrature implementation #254
Added Gauss-Turán quadrature implementation #254
Conversation
""" | ||
Gauss Turan Quadrature using xgt51 and agt51 | ||
""" | ||
struct GaussTuran{B} <: SciMLBase.AbstractIntegralAlgorithm |
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.
export GaussTuran
is missing from src/Integrals.jl
I'm wondering how this ties in with the package that I'm working on that does this more generally: https://github.com/SouthEndMusic/GaussTuranQuadrature.jl. We've had some discussions about this in the issue (#245) and in my PR (#252) and on slack: https://julialang.slack.com/archives/C66NPKCQZ/p1721823029663799. As suggested by @ChrisRackauckas, higher order derivatives can efficiently be computed with TaylorDiff.jl, but Gauss-Turán quadratures are probably most useful when you can compute the higher order derivatives cheaply (analytically). |
Although I would say this PR is largely superseded by #252, whose goal of computing the Gauss-Turán rules is broader, there are still a few aspects of this PR that would be good to include in GaussTuranQuadrature.jl, including:
I may make a PR to GaussTuranQuadrature.jl with these changes |
That would be much appreciated (: |
struct GaussTuran{B} <: SciMLBase.AbstractIntegralAlgorithm | ||
n::Int # number of points | ||
s::Int # order of derivative | ||
ad_backend::B # for now ForwardDiff | ||
end |
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.
This struct can be deleted since it is already defined in Integrals.jl
maxiters=nothing) | ||
integrand = prob.f | ||
a, b = domain | ||
return gt51((x) -> integrand(x, p), a, b) |
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.
This function should return a SciMLBase.jl solution, such as SciMLBase.build_solution(prob, alg, val, nothing, retcode = ReturnCode.Success)
|
||
a, b = 0.0, π | ||
|
||
result = solve(IntegralProblem((x, p) -> sin(x), (0, π)), GaussTuran(5, 2, ForwardDiff)) |
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.
result = solve(IntegralProblem((x, p) -> sin(x), (0, π)), GaussTuran(5, 2, ForwardDiff)) | |
result = solve(IntegralProblem((x, p) -> sin(x), (0, π)), GaussTuran(5, 2, ForwardDiff)).u |
@@ -0,0 +1,61 @@ | |||
module IntegralsGaussTuran |
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.
module IntegralsGaussTuran | |
module IntegralsGaussTuranExt |
This file should be renamed to IntegralsGaussTuranExt.jl and the Project.toml
file is missing a line IntegralsGaussTuranExt = "ForwardDiff"
in the [extensions]
section
@AbhayBhandarkar I've left some additional comments to get this pr to work. I hope this helps you understand how this works. I will also open a pr to your branch. I will close this pr because I have tested that the nodes and weights here do not give the full theoretical accuracy of a Gauss-Turán rule whereas those in #252 appear to. |
This draft pull request implements the Gauss-Turán quadrature method to the Integrals.jl package. The implementation includes:
gt51
):xgt51
) and weights (agt51
).While the code works within the Julia environment in Visual Studio Code IDE, testing within Integrals.jl is still an issue that needs tending to. Feedback and suggestions are welcome to improve this implementation. Thank you so much for this opportunity!