Integrals.jl icon indicating copy to clipboard operation
Integrals.jl copied to clipboard

Specify initial evaluation points [Feature Req.]

Open agerlach opened this issue 3 years ago • 4 comments

For DiffEqUncertainty expectation() applications, distributions that are narrow relative to the support can lead to incorrect results via Quadrature as the integrand is not sampled at points w/ non-zero (numerically) joint pdf values.

B/c we know the pdf, it would advantageous to "seed" or initialize any adaptive quadrature methods w/ the mean and/or random samples from the distribution.

E.g. for a simple linear system, u'=p*u, with uncertain IC,

u0_dist = [truncated(Normal(3.0,2.0),-1000,1000)]

will produce expectations of 0

while

u0_dist = [truncated(Normal(3.0,2.0),3-1000,3+1000)]

produces the correct result. The reason being that the midpoint of the integration domain is used as initial quadrature points in most algorithms supported.

agerlach avatar Jul 20 '20 14:07 agerlach

@agerlach I agree that localized integrands can be very difficult for quadrature schemes to resolve without any prior information, however I don't think this issue can be fixed at the level of Integrals.jl. In particular, not all libraries used by Integrals.jl support breakpoints for the interval of integration. Did DiffEqUncertainty.jl end up making a fix, such as breaking up the integration domain into cubes refined about the peaked distribution?

I would like to attempt fixing related issues like #41 and #160, possibly with a solution like AsDomain from DomainSets.jl. However, I don't think there is ever going to be a universal fix, since not all algorithms and domain discretizations are compatible.

lxvm avatar Oct 30 '23 22:10 lxvm

I agree completely. When Chris and I started this package we had some discussions about develop new scheme to support some other work. Hence this issue. PR #188 extends the interface to add support for this w/ CUBA.

agerlach avatar Oct 31 '23 12:10 agerlach

Interesting! I didn't know that Cuba supported features like that. I wonder if the information about initial points could be specified with more expressive domain types that have "special" points. For example, QuadGK.jl has something similar where the user provides an interval with breakpoints (a, b, c...) that is used to partition the interval. In that case, QuadGK clusters points about the breakpoints.

lxvm avatar Nov 01 '23 02:11 lxvm

One issue with these Cuba features is it is not readily obvious of how to use them. Cuba.jl just mimics the Cuba C library docs and those are sparse on details. E.g., when trying to use the peakfinder option with Divonne, the best resource I could find to figure out how to actually use it was the pythons bindings for Cuba.

agerlach avatar Nov 01 '23 13:11 agerlach