adaptive sampling, call for collaborators

559 views
Skip to first unread message

Andrei Berceanu

unread,
Jun 19, 2014, 7:34:23 AM6/19/14
to juli...@googlegroups.com
I think a really nice adition to Julia (either in Base or as a separate package) would be the possibility of adaptive sampling. I am thinking in terms of a stand-alone routine for n-dimensions, which could then be called for plotting, numerical integration, etc.
I am currently doing research in theoretical physics and in my day to day work I often encounter functions which are very expensive to compute and plotting them would go so much faster with the help of a good sampling scheme.
As a starting point, see this question on the Mathematica(R) stack exchange community:

http://mathematica.stackexchange.com/questions/216/adaptive-sampling-for-slow-to-compute-functions-in-2d

I think this should not be hard to implement in Julia, but we could easily do better, for instance with support for parallel evaluations.

So, first I would like to know who would be interested in having such functionality in Julia, and then we decide how to proceed.

Stefan Karpinski

unread,
Jun 19, 2014, 10:35:29 AM6/19/14
to juli...@googlegroups.com
Yes, absolutely this would be good to have. Lots of packages would probably need adaptive sampling so having this as it's own low-level package might be good.

Andrei Berceanu

unread,
Jun 20, 2014, 5:22:01 PM6/20/14
to juli...@googlegroups.com
Ok, great. So how do I start? Should I take it that you are willing to collaborate? I am quite new to Julia btw.

Stefan Karpinski

unread,
Jun 20, 2014, 9:41:45 PM6/20/14
to Julia Dev
No, I'm afraid I don't have time to work on this, but I wouldn't be surprised if someone might be interested. The best way is probably to start trying to write something and ask questions as you encounter problems. You could try porting code from another language like Python, R, Matlab or Mathematica. Keep in mind if you do this, however, that the resulting code will probably have to retain the same license as the original code.

Aerlinger

unread,
Jun 21, 2014, 5:18:42 PM6/21/14
to juli...@googlegroups.com
Not sure how much time I have, but would love to see something like this. 

Alireza Nejati

unread,
Jun 22, 2014, 6:40:12 PM6/22/14
to juli...@googlegroups.com
I'm also a bit strapped for time but I'd like to see something like this and might be able to help out a little. Have you made a choice on the algorithm to use?

Andrei Berceanu

unread,
Jun 23, 2014, 8:20:15 AM6/23/14
to juli...@googlegroups.com
The 1D case is more or less straight-forward and implementations already exist. Therefore I propose we start by focusing on the 2D case. The plan would be

* Make a triangulation on a set of points
* Implement optimization moves as in this paper[1]
* Implement a greedy algorithm for minimizing the curvature of the resulting triangulation

[1] http://graphics.berkeley.edu/papers/Narain-AAR-2012-11/Narain-AAR-2012-11.pdf

Is there any package for tessellation/triangulation in Julia?

Alireza Nejati

unread,
Jun 23, 2014, 6:03:49 PM6/23/14
to juli...@googlegroups.com
I had a look at that paper and it doesn't seem to be very well suited for this problem. What Mathematica does seems to be far simpler. It just starts with a rectangular mesh then subdivides this rectangular mesh by simply adding a vertex at the centroid.


On Thursday, June 19, 2014 11:34:23 PM UTC+12, Andrei Berceanu wrote:

Stefan Karpinski

unread,
Jun 24, 2014, 8:12:19 AM6/24/14
to Julia Dev
If it's good enough for Mathematica, then I suspect it's good enough for us.

Andrei Berceanu

unread,
Jun 27, 2014, 1:54:04 PM6/27/14
to juli...@googlegroups.com
Some ammends are in order, I do not claim to have come up with the plan from the previous post, it was actually Anton Akhmerov's idea (http://antonakhmerov.org).

About the triangulation, as far as I understood, the main advantage of implementing a triangulation as the one used in that paper is that it behaves much better around sharp edges.
However, it may be too ambitious as a starting point, so I suggest to start from Delaunay triangulation. That should already give an improvement over Mathematica's performance.

Andrei Berceanu

unread,
Jun 27, 2014, 2:04:32 PM6/27/14
to juli...@googlegroups.com
The Bowyer–Watson algorithm looks simple enough to implement in Julia.
http://en.wikipedia.org/wiki/Bowyer%E2%80%93Watson_algorithm
http://paulbourke.net/papers/triangulate/

Alireza Nejati

unread,
Jun 29, 2014, 9:02:20 PM6/29/14
to juli...@googlegroups.com
I agree with Stefan here. Using complex triangulations would make it much harder to calculate curvature for subdivision. I don't know how Mathematica subdivides regions, but it makes sense to first try a simple subdivision of a square grid. Not a delaunay triangulation. At any rate, I wouldn't be able to help with a delaunay-based algorithm.

If you feel delaunay is the way to go, why don't you write some initial code and put it up on git.


On Thursday, June 19, 2014 11:34:23 PM UTC+12, Andrei Berceanu wrote:

Jason Knight

unread,
Jul 2, 2014, 11:24:28 AM7/2/14
to juli...@googlegroups.com
I would be very interested in a package like this. Although I have a somewhat specialized use case and have already implemented the following:

I'm interested in moderate dimensional summations over large discrete spaces (Poisson based Bayesian statistical models), so I implemented first an adaptive cube division technique (somewhat hard to read code here), but then moved to a KD-tree division approach to handle large aspect ratio distributions (slightly cleaner code here). Because of the specific nature of my problem, my stopping criteria is somewhat specific (I know elements of the summation vector need to be unity), but the division criteria are somewhat more generic (subdividing regions with high 'discrepancy' using a naive finite difference approach). 

I wasn't aware this fell under the umbrella of adaptive sampling (my original inspiration came from Steven Johnson's Cubature.jl package, which (unsurprisingly) is inefficient for my discrete domains).

The code could be adapted to continuous regions quite easily I believe. 

Jason
Reply all
Reply to author
Forward
0 new messages