On Tue, 23 Mar 2010, Hans W Borchers wrote:
> Barry Rowlingson <b.rowlingson <at> lancaster.ac.uk> writes:
> > On Mon, Mar 22, 2010 at 4:28 PM, Hans W Borchers
> >
> > <hwborchers <at> googlemail.com> wrote:
> > > Still I believe that a clever approach might be possible avoiding the
> > > need to call a commercial solver. I am getting this hope from one of
> > > Jon Bentley's articles in the series Programming Pearls.
> >
> > Is this the 'Largest Empty Rectangle' problem?
> >
> > http://en.wikipedia.org/wiki/Largest_empty_rectangle
>
> Dear Barry,
>
> thanks for this pointer. I never suspected this problem could have a name
> of its own. Rethinking the many possible applications makes it clear: I
> should have searched for it before.
>
> I looked in some of the references of the late 80s and found two algorithms
> that appear to be appropriate for implementation in R. The goal is to solve
> the problem for n=200 points in less than 10-15 secs.
>
How about less than 2 seconds? [And 500 points in less than 15 seconds - on a
2-year-old
DELL Optiplex GX755.]
The implementation below (at end) loops over all 'feasible' pairs of x values,
then
selects the largest rectangle for each pair, subject to your specified
constraints. I
have no idea if it implements a previously published algorithm.
Other constraints are reasonably easily accommodated.
HTH,
Ray Brownrigg
> Thanks again, Hans Werner
>
> > I had a look at some of the references from Wikipedia, but they all
> > follow a similar pattern, one I have noticed in many computer science
> > journal articles:
> >
> > 1. State a problem that looks tricky.
> > 2. Say "We have an efficient algorithm for the problem stated in #1"
> > 3. Proceed to derive, using much algebra and theory, the efficient
> > algorithm. 4. Stop.
> >
> > The idea of actually producing some dirty, filthy, actual code to
> > implement their shiny algorithms never seems to cross their minds.
> >
> > I also found a similar question from 2008 asked on the R-sig-geo
> > mailing list. That didn't get much help either!
> >
> > Sorry.
> >
> > Barry
N = 200
x <- runif(N)
y <- runif(N)
ox <- order(x)
x <- x[ox]
y <- y[ox]
x <- c(0, x, 1)
y <- c(0, y, 1)
plot(x, y, xlim=c(0, 1), ylim=c(0, 1), pch="*", col=2)
omaxa <- 0
for(i in 1:(length(x) - 1))
for(j in (i+1):length(x)) {
x1 <- x[i]
x2 <- x[j]
XX <- x2 - x1
if (XX > 0.5) break
yy <- c(0, y[i:j], 1)
oyy <- order(yy)
yy <- yy[oyy]
dyy <- diff(yy)
whichdyy <- (dyy <= 0.5) & (dyy >= 0.5*XX) & (dyy <= 2*XX)
wy1 <- yy[whichdyy]
if (length(wy1) > 0) {
wy2 <- yy[(1:length(dyy))[whichdyy] + 1]
k <- which.max(dyy[whichdyy])
maxa <- (x2 - x1)*(wy2[k] - wy1[k])
if (maxa > omaxa) {
omaxa <- maxa
mx1 <- x1
mx2 <- x2
my1 <- wy1[k]
my2 <- wy2[k]
}
}
}
lines(c(mx1, mx2, mx2, mx1, mx1), c(my1, my1, my2, my2, my1), col=2)
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
