Re: Implicit plots/level curves possible?

[Hans Hagen <j.hagen@xs4all.nl>]

On 10/8/2018 11:25 PM, Aditya Mahajan wrote:
> On Mon, 8 Oct 2018, Alan Braslau wrote:
> 
>> On Mon, 8 Oct 2018 16:00:10 -0400 (EDT)
>> Aditya Mahajan <adityam@umich.edu> wrote:
>>
>>> On Sun, 7 Oct 2018, Hans Hagen wrote:
>>>
>>>> On 10/7/2018 7:14 PM, Alan Braslau wrote:
>>>>> On Sun, 7 Oct 2018 17:25:35 +0200
>>>>> "Mikael P. Sundqvist" <mickep@gmail.com> wrote:
>>>>>
>>>>>> ContourPlot[2 x^5 + x y + y^5 == 0, {x, 0, 2}, {y, -2, 1/2}]
>>>>>
>>>>> Brut force:
>>>>> [...]
>>>>>
>>>> as this takes some time here's a cheat:
>>>>
>>>> \starttext
>>>>
>>>> \startbuffer[demo]
>>>> [...]
>>>> \stopbuffer
>>>>
>>>> \startTEXpage
>>>>     \typesetbuffer[demo]
>>>> \stopTEXpage
>>>> \stoptext
>>>>
>>>> a next run the already prepared buffer will be taken unles it has been
>>>> changed.
>>>
>>> I thought that this will also be a good usecase of showing Lua+MP
>>> interaction. I wrote the code below following the metafun manual, but 
>>> I cannot
>>> get it to compile. What am I missing?
>>>
>>>
>> You need to put it into the mp namespace:
> 
> Thanks. The metafun manual is confusing in this regard and I got the 
> impression that any lua namespace could be used.
> 
>> (then eps should be made a linear function of xi)
> 
> Oh, I missed that. Thanks.
> 
>> (and, indeed this is much faster than calculating in MP)
> 
> Yes. Iterating over 10^6 values on a 1GHz computer should roughly take 

hm, my 3.4 gig laptop cpu needs 0.18 sec for the slightly more that i 
million steps ... it depends on the function too

> 1ms in any reasonable programming language. Metapost for loops work with 
> macro expansion, which can be very expensive for large loops.

generalized:

\startluacode
----- abs = math.abs

local contour = { }
local data    = { }
local origin  = { 0, 0 }
local length  = 0

local function generate(f, x_min, x_max, y_min, y_max)
    local points = { }
    local length = 1000
    local eps    = 1e-3
    local spe    = -eps
    local n      = 0
    local code   = "return function(x,y) return " .. f .. " end"
    local action = load(code)
    if action then
         action = action()
    end
    for xi = x_min, x_max, (x_max - x_min)/length do
        for yi = y_min, y_max, (y_max - y_min)/length do
          -- if abs(action(xi,yi)) < eps then -- 10% gain with:
             local v = action(xi,yi)
             if v < eps and v > spe then
                 n = n + 1
                 points[n] = { xi, yi }
             end
        end
    end
    return points, n
end

function mp.Countour(...)
     data, length = generate(...)
end

function mp.ContourN()
     mp.print(length)
end

function mp.ContourPoint(i)
     mp.pair(data[i] or origin)
end

function mp.ContourPath(f,...)
     if f then
         data, length = generate(f,...)
     end
     mp.path(length > 0 and data or origin)
end

\stopluacode

\starttext
     \startMPpage[instance=doublefun]
       % lua.mp.Countour("2*x^5  + x*y + y^5", 0, 2, -1, 0.5) ;
       % draw lua.mp.ContourPath() withpen pencircle scaled .01 ;
         draw lua.mp.ContourPath("2*x^5  + x*y + y^5", 0, 2, -1, 0.5) 
withpen pencircle scaled .01 ;
         setbounds currentpicture to (0,-2) -- (2,-2) -- (2,.5) -- 
(0,.5) -- cycle ;
         currentpicture := currentpicture xsized 5cm ;
         picture pic ; pic := currentpicture ;
         drawarrow llcorner pic -- lrcorner pic ;
         drawarrow llcorner pic -- ulcorner pic ;
         label.rt ("$x$", lrcorner pic) ;
         label.top("$y$", ulcorner pic) ;
         for x=0 step .5 until 2 :
             label.bot(decimal x,(x/2)[llcorner pic,lrcorner pic]) ;
         endfor ;
         for y=0 step .5 until 2.5 :
             label.lft(decimal (y-2),(y/2.5)[llcorner pic,ulcorner pic]) ;
         endfor ;
     \stopMPpage
\stoptext


-----------------------------------------------------------------
                                           Hans Hagen | PRAGMA ADE
               Ridderstraat 27 | 8061 GH Hasselt | The Netherlands
        tel: 038 477 53 69 | www.pragma-ade.nl | www.pragma-pod.nl
-----------------------------------------------------------------
___________________________________________________________________________________
If your question is of interest to others as well, please add an entry to the Wiki!

maillist : ntg-context@ntg.nl / http://www.ntg.nl/mailman/listinfo/ntg-context
webpage  : http://www.pragma-ade.nl / http://context.aanhet.net
archive  : https://bitbucket.org/phg/context-mirror/commits/
wiki     : http://contextgarden.net
___________________________________________________________________________________