Re: [racket] Plot woes
Neil Toronto
On 05/22/2012 02:41 AM, Jens Axel Søgaard wrote:
> Hi Neil,
>
> I am attempting to draw 1/x in using plot.
> Ideally I'd like to be able to write:
>
> (lambda (x) (if (zero? x) #f (/ x))
>
> in order to plot a function with singularities.
I've thought of such things before. Unfortunately, plot doesn't know
where any singularities are, so it can't specifically sample the
function at those points. I could probably allow the user to tell plot
where they are, but I can't think of a good API for it. I'd want
something general enough to handle 1D, 2D and 3D functions, which allows
the user to specify countably many missing values, and direction of
(dis)continuity.
> Since the contract for plot requires I must use numbers,
> I attempted to split the domain in two from -1 to -epsilon
> and from +epsilon to 1.
>
> That trick didn't work (see attached image).
>
> Am I missing something?
You're not missing anything. There's currently a problem with the
renderers returned by `function', `surface3d' and others: they sample
functions on the entire plotted domain, not just within the functions'
stated domains. To split a function in half and not have it draw a line
between the halves, you currently have to do something like this:
(plot (list (function (λ (x) (if (x . < . 0) (/ x) +nan.0)) -2 0)
(function (λ (x) (if (x . > . 0) (/ x) +nan.0)) 0 2)))
This also causes performance problems when there are many small function
renderers. Fixing it is next on my to-do list.
> Apropos something completely different:
>
> At the page "Introduction" it would be nice if
> the first example included (require plot).
>
> http://docs.racket-lang.org/plot/intro.html?q=plot
Duly noted.
Neil ⊥
____________________
Racket Users list:
http://lists.racket-lang.org/users
Matthias Felleisen
On May 22, 2012, at 8:05 AM, Neil Toronto wrote:
> I've thought of such things before. Unfortunately, plot doesn't know where any singularities are, so it can't specifically sample the function at those points. I could probably allow the user to tell plot where they are, but I can't think of a good API for it. I'd want something general enough to handle 1D, 2D and 3D functions, which allows the user to specify countably many missing values, and direction of (dis)continuity.
-- Matthias
Jens Axel Søgaard
> On May 22, 2012, at 8:05 AM, Neil Toronto wrote:
>
>> I've thought of such things before. Unfortunately, plot doesn't know where any
> singularities are, so it can't specifically sample the function at those points. I could > probably allow the user to tell plot where they are, but I can't think of a good API
> for it. I'd want something general enough to handle 1D, 2D and 3D functions,
> which allows the user to specify countably many missing values, and direction of
> (dis)continuity.
It is in no way a simple matter to plot functions with singularities. Especially
not if points are not known in advance. A simple heuristic is better than none
though.
Consider the one dimensional case. Let's for the sake of argument say
that we sample the function f in two values a and b such that f(a) is
below y-min
and f(b) is over y-max.
If there is a singularity between a and b, then the points P(a,f(a))
and Q(b,f(b)) should not be connected, and if there is no singularity
(but a pretty steep ascent)
then the points should be connected.
By sampling points between a and b it ought to be possible to
give a heuristic of whether a singularity has been found or not.
Thinking aloud: Would a comparison of the rate ascent between a and b to the
mean rate of ascent of all sampled pairs of neighbors give anything usable?
If you know any references to papers and/or algorithm
on the subject, then do give a pointer. I was only able to find
this paper of Fateman:
http://www.cs.berkeley.edu/~fateman/papers/intervalplot.ps
which present a some nice examples of things to look out for.
(Just so I don't forget where I found it)
Matlab code using succesive refinement:
http://www.stat.colostate.edu/~estep/paov/plotfuns.m
--
Jens Axel Søgaard
Eli Barzilay
>
> On May 22, 2012, at 8:05 AM, Neil Toronto wrote:
>
> > I've thought of such things before. Unfortunately, plot doesn't
> > know where any singularities are, so it can't specifically sample
> > the function at those points. I could probably allow the user to
> > tell plot where they are, but I can't think of a good API for
> > it. I'd want something general enough to handle 1D, 2D and 3D
> > functions, which allows the user to specify countably many missing
> > values, and direction of (dis)continuity.
>
> Isn't this where a well-specified exception interface could help? It
> won't be fast but at least you can specify singularities as
> exceptions.
numbers? IOW, this problem can be fixed with
(lambda (x) (/ (exact->inexact x)))
Maybe introduce that as some keyword option, but it looks like for
plots floating point numbers should be the defaults. Without that,
you get nasty surprises like:
(plot (function / -249 249)) --> works
(plot (function / -249 250)) --> error
Also, it looks like the bounds are rounded and exacted, since
(plot (function / -249.0 250.0))
throws an error too. If that get dropped, then it can be an indirect
way to ask for floats.
--
((lambda (x) (x x)) (lambda (x) (x x))) Eli Barzilay:
http://barzilay.org/ Maze is Life!
Neil Toronto
> On 05/22/2012 02:41 AM, Jens Axel Søgaard wrote:
>> Since the contract for plot requires I must use numbers,
>> I attempted to split the domain in two from -1 to -epsilon
>> and from +epsilon to 1.
>>
>> That trick didn't work (see attached image).
>>
>> Am I missing something?
>
> You're not missing anything. There's currently a problem with the
> renderers returned by `function', `surface3d' and others: they sample
> functions on the entire plotted domain, not just within the functions'
> stated domains. To split a function in half and not have it draw a line
> between the halves, you currently have to do something like this:
>
> (plot (list (function (λ (x) (if (x . < . 0) (/ x) +nan.0)) -2 0)
> (function (λ (x) (if (x . > . 0) (/ x) +nan.0)) 0 2)))
>
> This also causes performance problems when there are many small function
> renderers. Fixing it is next on my to-do list.
I just pushed a fix for this. The following now works as you'd expect:
#lang racket
(require plot)
(plot (list (function (compose / exact->inexact) -1 0)
(function (compose / exact->inexact) 0 1)))
(plot (list (function / -1 -0.01)
(function / 0.01 1)))
Eli Barzilay
>
> I just pushed a fix for this. The following now works as you'd expect:
>
> #lang racket
>
> (require plot)
>
> (plot (list (function (compose / exact->inexact) -1 0)
> (function (compose / exact->inexact) 0 1)))
>
> (plot (list (function / -1 -0.01)
> (function / 0.01 1)))
(plot (list (function / -249.0 250.0)))
Neil Toronto
> On Thursday, Neil Toronto wrote:
>>
>> I just pushed a fix for this. The following now works as you'd expect:
>>
>> #lang racket
>>
>> (require plot)
>>
>> (plot (list (function (compose / exact->inexact) -1 0)
>> (function (compose / exact->inexact) 0 1)))
>>
>> (plot (list (function / -1 -0.01)
>> (function / 0.01 1)))
>
> This still doesn't work:
>
> (plot (list (function / -249.0 250.0)))
That's a separate issue that I haven't decided how to resolve yet.
Neil ⊥
Jens Axel Søgaard
Attached is my version of a 2d-plotter that handles singularities.
Give it a spin with your favorite ill-behaved functions and
report back with ones badly handled.
/Jens Axel
Neil Toronto
For approaches like this, I think the definition of "well-behaved" is
"locally Lipschitz continuous with Lipschitz constant < K", where K
depends on the constants used to detect discontinuities (e.g.
machine-epsilon, number-of-regions). For everywhere-differentiable
functions, K is the maximum absolute derivative on the domain plotted,
scaled by various plot- and detection-specific factors.
(Yes, "Lipschitz" is like, the most unfortunate name for a mathematical
property EVER. It takes undergraduate math majors a whole semester to
get over it, if they ever do.)
Knowing that, it doesn't take long to find an ill-behaved function. We
just need one that's differentiable but not locally Lipschitz:
(λ (x) (* (expt (abs x) 3/2) (sin (/ x))))
Things like this will always be ill-behaved when plotted on some domain,
no matter how fine the discontinuity detection is.
Further, there are locally and globally Lipschitz functions whose K is
too large. (I didn't spend time finding one, though.) For each one, some
constant assignment in the discontinuity detection will make them plot
nicely. But for any constant assignment, there are infinitely many of them.
</analysis_lesson>
That doesn't mean I won't try something like this. The definition of
"functions without discontinuities that are well-behaved when plotted
non-adaptively" is also "locally Lipschitz continuous with Lipschitz
constant < K", where K depends on different constants like the
discretization step size. I would be perfectly happy with an adaptive
sampler or just discontinuity detection if I could show that its set of
ill-behaved functions is no larger than those for the current sampler.
That'll take some time.
Also, I'm still considering letting users specify discontinuities. That
would allow them to tell plot what kind they are: removable, removable
singularity, right-step, left-step, etc. Those kinds of properties are
pretty much undetectable.