@janhardo What long calculation? I get the solution _C1*x*exp(-y) + _C2.

@acer Acer, would you please explain the purpose of

a:= 'a':  b:= 'b':  c:= 'c':

that appears before the gc() that appears before timed code? My guess is that its purpose is to force garbage collection of %, %%, and %%%.

Also, I'm guessing from the way that you recoded ZigZag that compiled procedures can't use op(procname). I used that to avoid the extra procedure layer (i,j)-> ZigZag(i,j,n).

@Rouben Rostamian  Thanks. When you edit a Reply, there's a pull-down at the bottom that says Parent. That's how you move Replies to different Answers. 

@Rouben Rostamian  Here is the inverse of the above procedure, i.e., it returns the indices given the entry:

ZigZagInv:= (k::posint)-> 
local n:= op(procname);
    if k > n*(n+1)/2 then n+1 -~ thisproc(n^2+1-k)
    else
        local d:= ceil((sqrt(1+8*k)-1)/2), r:= k - (d-1)*d/2;
        [r, d+1-r][(-1)^d*[1,2]]
    fi
:
#Example usage
M:= Matrix(8$2):
for k to 64 do M[ZigZagInv[8](k)[]]:= k od:
M;

@acer The design that units are names rather than strings (or records or objects) seems to me to be inherently flawed.

@Sayed From my code above, the set of coefficients alone is simply {C}.

@Carl Love In other words, diff acts as if the part of it that handles unevaluated functions as the first argument were coded like this:

restart:
old_diff:= eval(diff):
unprotect(diff);
diff:= overload([
    proc(f::function, x::name, $)
    option overload;
    local F:= op(0,f), d:= cat(`diff/`, `if`(F::indexed, op(0,F), F)), r;
        r:= `if`(F::indexed, d[op(F)], d)(op(f), x);
        `if`(r::specfunc(d), 'procname'(args),  r)        
    end proc,

    old_diff
]):
protect(diff, old_diff);    
        

As far as I can telll, the part about passing the index isn't documented, and VV must've made an educated guess that it would work.

@janhardo 

Do

showstat(Drawlines);

and post the results here. I'll probably be able to tell you how it works and to fix it.

The Google Drive links that you posted to this code require a password, so I wasn't able to download it.

@mmcdara It must be the case that diff passes f's index as the index to `diff/f`.

@janhardo Can you show examples of Drawlines or Drawplanes not working?

@Umang Varshney 

Axel suggested that you do substitutions and then use fsolve, not solve.

@666 jvbasha @tomleslie

The problem with changing ax is that the endpoints of the solution interval depend on it. Thus, the solution for one value of ax will not be accepted as even an approximate solution for another value of ax.

@gawati2611 Your header says "Maple 2020", but I guess that's not true. You need to change this:

PolygonAsLines:= (P::list([numeric,numeric]))->
local t, k;
    plot(
        [seq([((P[k+1]-~P[k])*~t +~ P[k])[], t= 0..1], k= 1..nops(P)-1)],
        _rest
    )
:

to this:

PolygonAsLines:= proc(P::list([numeric,numeric]))
local t, k;
    plot(
        [seq([((P[k+1]-~P[k])*~t +~ P[k])[], t= 0..1], k= 1..nops(P)-1)],
        _rest
    )
end proc:

@janhardo Assuming that the package is freely distributable, without copyright issues, could you post a link to it?

The fact that the read command didn't give any errors is a very good sign. It's possible that you won't be able to use the with command for this, but that's not a major obstacle. Did help pages come with the package?

Note that Maple has a substantial LinearAlgebra package already, and even a Student:-LinearAlgebra.

@Anthrazit Unfortunately MaplePrimes doesn't correctly attach files with .maple extension. Can you make it a .zip file and upload it again?