The command is map(limit, M, [x1= 0, x2= 0]), not the map2 command that you gave. You must respect the order of the arguments when using map and its relatives.

I disagree with Rouben's method. I think that the coeffs command is needed for this:

V:= [u[2], u[2,2], u[2,2,2]]:
C:= coeffs(collect(f, V, distributed), V, 'terms'):
Coeffs:= table(sparse, [terms]=~ [C]);

The following should be faster if you're just as likely to be checking residues as nonresidues:

if NumberTheory:-QuadraticResidue(x,n) then 
    a:= NumberTheory:-ModularSquareRoot(x,n);
    b:= b+1 #or b++ in Maple 2019+
else
    c:= c+1 #or c++ in Maple 2019+
fi;

And if you want the last returned value rather than the last computed value and also use Joe's "last" trick, then do this:

Foo:= ()-> (foo("last"):= foo(args)): #Wrapper procedure
foo:= proc(x) option remember; x^2 end proc:
foo("last"):= "foo not yet used":
Foo(2);
                               4
foo("last");
                               4

Foo(3);
                               9

foo("last");
                               9

Foo(2);
                               4

foo("last");
                               4

Of course you can't apply 1/z or Joukowsky to regions containing z=0.

Here's an example of Joukowsky applied to a triangle:

restart:
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
    )
:
J:= z-> z + 1/z; #Joukowsky
Jw:= evalc([Re,Im](J(x + y*I)));
commonopts:= view= [1..2.5, 0..1], size= [500$2], scaling= constrained:
R:= plots:-display(
    PolygonAsLines([[1,0],[2,0],[1,1],[1,0]], color= blue, thickness= 3),
    seq(plot(i, x= 1..2-i, color= black, thickness= 0), i= 0..1, 0.1),
    plot([seq([i, t, t= 0..2-i], i= 1..2, 0.1)], color= black, thickness= 0),
    commonopts
);
F:= plottools:-transform(unapply(Jw, (x,y))):
plots:-display(F(R), commonopts);

By far the easiest way to use parallel programming is through the Seq and Map commands of packages Grid or Threads. You simply take a seq or map command from your sequential code and replace it with the package command:

foo:= i-> i^2:
Threads:-Seq(foo(i), i= 1..9);
Grid:-Seq(foo(i), i= 1..9);
X:= [$1..9]:
Threads:-Map(foo, X);
Grid:-Map(foo, X);

So the complexity of Threads:-Task is not needed here, whereas it may be needed to parallelize code that can't be expressed with seq or map such as recursive code.

Code that writes to global variables (without locking them) is not "thread safe". The vast majority of Maple's high-level symbolic commands do this, thus separate invocations are not truly "independent" (as you said), thus Threads can't safely be used. If you try to use Threads with code that isn't threadsafe, the results are very unpredictable. You may get incorrect results, an error, or a crash.

On the other hand, usage of global variables doesn't affect Grid, but you pay for that with a much higher startup time and memory usage.

@nm z^(-1) obviously has a z in it!  

It helps to think of breaking expressions down into trees:

z^(-1):   Root node:  ^
2 branch nodes:  z  and  -1

There is a z node.

1/(x*y):   Root node:       *
2 branch nodes: x^(-1) and y^(-1)

There's no x*y node.

You could substitute x = t/y.

The command ToInert(1/(x*y)) ​​​​​​or ToInert(1/z) will show you the trees (in 1D function form).

The purpose of this Answer is not to address the relationship that you report between kernel crashes and input mode. While that's a very serious issue (which you've reported clearly), I don't feel like going down that rabbit hole at the moment. Rather, the purpose of this Answer is to advocate that sum (with a minor adjustment) is better than add for this computation. So, I disagree with the other respondents on this point.

Your sum is a geometric sum that can easily be put into a very simple closed form. This was known even to medieval bankers, centuries before the invention of calculus. So, it's a shame that Maple's sum automatically converts this to direct addition of terms (sum does that when the number of terms is small). We can avoid that by making the upper limit of summation temporarily symbolic:

restart:
RR:= 4: #return rate PV/C 
n:= 12:
simplify(sum(1/(1+r)^(155/10/365 + k/%n), k= 1..%n)) assuming %n::posint;
                                      r                 
                     -----------------------------------
                                      /            /1 \\
                                      |            |--||
                            (761/730) |            \%n/|
                     (1 + r)          \-1 + (1 + r)    /

fsolve(value(%) = RR, r);
                                 8.161501788

The "Length of output exceeds limit" message is not an error or even a warning. Maple has computed your result. It just chooses to not display it in prettyprinted form. If you want to see the result in plaintext form, use lprint(%).

The presence of this message and whether you see the output has no effect whatsoever on any subsequent command that uses that output.

You can increase the output limit to any value that you want. I recommend that you don't do this. 

The best way to convert an expression sequence into a vector or an expression sequence of vectors into a matrix is to use angle-bracket operators:

f:= x-> <x | x+1>: # *row* vector
test:= <seq(f(i), i= 0..1, 0.1)>; #column of row vectors, i.e., a matrix

#You can plot a 2-column matrix without first extracting its x- and y-vectors:
plot(test, style= pointline);

#extraction of sequence of column vectors (don't know why this would be useful, but you asked):
ColVecs:= seq(test[..,k], k= 1..upperbound(test)[2]);

In your worksheet, you've repeatedly made the mistake of attempting to index a structure S with parentheses as in S(1). Square brackets should usually be used for indexing, as in S[1] and test[.., k]. With Vectors, Matrices, and Arrays (known collectively as rtables), in some cases you may use parentheses; this is known as programmer indexing. You used it correctly, albeit haphazardly, in your procedure f in your worksheet. I'd recommend that you stay away from this for the time being.

The Maple structures for which square-bracket indexing can be used for member extraction are expression sequences, lists, sets, tables, rtables, modules (which includes Records), and objects that overload `?[]`. Modules may also be indexed as A:-B rather than A[B] if B is explicitly an export of the module (in which case I recommend using the :- operator).

Implies(And(d1::even, d2::odd), d1+d2 <= 9);

The prettyprinted output of this will be in standard logical notation.

Try this:

TypeTools:-AddType(
    special_polynomial,
    (p, x::name)-> 
    local 
        X, k:= 0,
        T:= subsindets(
            p, 
            identical(x)^{constant, And(name, Not({identical(x), X[posint]}))}, 
            ()-> X[++k]
        )
    ;
        type(T, polynom(anything, {x} union indets(T, X[posint])))
);

Usage:

type(expr, special_polynomial(x));

I don't think that your way is good, for two reasons:

• I don't think that you understand what Maple's type algebraic is. Look up ?type,algebraic. Surely you don't want x^sin(x) + 1 to be considered a special polynomial.
• The scanning by remove is only through the top level of an expression--the terms in this case. Looking within the terms requires the deeper scan provided by subsindets (or evalindets). For example, your procedure would not call 3*x^(1/2) + 1 a special polynomial because it doesn't see the x^(1/2) because of the 3.

Also note that Pi, sqrt(3), etc., are type constant but not type numeric.

The operators numer and denom are mathematical, not purely syntactic like lhs and rhs. It's easy enough to write a syntactic operator to do what you want:

Apparent_Numer_Denom:= (e::`*`)-> 
    ((D,N)-> (N,1/D))(selectremove(type, e, anything^negative))
:
Apparent_Numer_Denom(expr);
 

ProductOfDiagonal:= (M::Matrix)-> mul(ArrayTools:-Diagonal(M));

Trying to guess your next Question, the sum of the k-wise products of the diagonal entries is

SumKwiseProducts:= (M::Matrix, k::nonnegint)-> 
    add(mul~(combinat:-choose([seq](ArrayTools:-Diagonal(M)), k)))
:

plots:-odeplot(sol, [t, diff(x(t), t)], t= 0..20);