Maple 18 was wrong because the "closed form" that you give is not valid for all x. In Maple 2018, this error has been corrected. Now you either need to put an assuming clause (such as assuming x > 0) or throw caution to the wind and include the formal option:

sum(1/(1+x)^t, t= 1..infinity, 'formal');
sum(1/(1+x)^t, t= 1..infinity) assuming x > 0;

Here's an example of what you're trying to do. The procedure immediately below takes an expression f and returns the Newton's method iteration that can be used to solve f=0 for "some variable". The user can specify the variable to use as that "with respect to" variable (to bind that variable is the logician's parlance) in the second argument or leave the second argument blank to let the procedure try to figure that out on its own.

Newton:= proc(f::algebraic, var::name:= ())
local x:= var;
   if x=() then #User didn't pass a 2nd argument.
      #The Not(constant) part is to prevent things like 
      #Pi from being considered variables. 
      x:= indets(f, And(name, Not(constant)));
      #indets always returns a set, which we process based on its
      #number of elements, nops.
      if nops(x)=1 then x:= x[] #Extract the one and only member.
      elif x={} then error "No variables found"
      else error "Multiple variables found"
      fi
   fi;
   x - f/diff(f,x) #Return Newton iteration.
end proc: 

Then try

Newton(y^2 - 2);
Newton(x^2 - a);
Newton(x^2 - a, x);
Newton(2);
Newton(2, x); #Intentional divide-by-zero error

Note that the parameter specification var::name:= () has both a type, name, and a default value, (), and that that default value doesn't match the type. Thus the procedure knows that if var=(), then the user definitely didn't pass a second argument. (Actually, it's impossible to pass the particular value () as an argument, but that's beside the point. This type/default mismatch trick works for any default.)

The example that you showed of plot(x^2) is unusual compared to most Maple commands in that the input is an expression, yet the command doesn't require the user to specify the variable to bind. (And that's what I did in Newton.) But the far more common situation is that the procedure assumes that if no variable is specifed then the function/expression argument is a procedure whose (first) argument is the variable to bind (as in the commands plot(x-> x^2, -2..2) and plot(sin, -Pi..Pi)). This is simpler to code. Here's an example:

Newton2:= proc(f, var::name:= ())
local x;
   if var=() then #User didn't pass a 2nd argument.
      #Treat f as a procedure, temporarily bind a local, and return a procedure:
      unapply(x - f(x)/D(f)(x), x)
   else
      #Bind var and return an expression
      var - f/diff(f,var)
   fi
end proc:

Note that binding a variable is mostly a logical/mental process rather than something that happens inside the computer. It is something that happens in the minds of the writer and user of the procedure that gives a particular variable a different status.

When Maple needs to bind a variable that will appear in the output of a user-level procedure but which was not specifed in the input, it almost always chooses one whose name is an underscore, followed by one or two letters, followed by a nonnegative integer. Mostly, this is a convention that is respected by programmers rather than something that is enforced by a compiler.

The relation is not generally true. It would be true if you changed 2*Pi to Pi. The argument of a complex number is defined (in Maple, and perhaps generally) to be in the interval (-Pi, Pi]. I think you're assuming that that interval is [0, 2*Pi).

Your expression doesn't necessarily simplify to 0, even with the real assumption. Let phi be 3*Pi/4, for example. Actually, it's only true for angles in the 1st quadrant.

I think that a dedicated procedure is a bit of overkill for this. Also, calling RandomVariable is redundant. If Lambda is a "container" (such as a Matrix or Vector from Excel, or any other source) of parameter values, and Alpha is either a single percentile level (e.g, Alpha:= 0.98) or a container of percentile levels matched to Lambda, then all that you need to do is

Statistics:-Quantile~(Poisson~(Lambda), Alpha);

The following procedure will extract all products of two or more functions from an expression, with the only exponent allowed being -1 (which in Maple is equivalent to being a factor of the denominator).

FunctionProducts:= e-> 
   select(
      type, 
      map2(select, type, indets(e, `*`), {function, function^identical(-1)}), 
      `*`
   )
:

This procedure extracts a nonnegative integer from the end of any symbol or string, using all available trailing digits from the end.

GetTail:= proc(X::{symbol, string})
local E:= ["A", StringTools:-Explode(X)[]], k;
     for k to nops(E) while StringTools:-IsDigit(E[-k]) do od;
     `if`(k=1, FAIL, parse(substring(X, -(k-1)..-1)))
end proc:

Here's one way to do it. I'm going to name the procedure Apply because it applies its first argument to its second argument.

Apply:= proc(f, x)
   return f(x)
end proc;

g:= x-> x^2;
Apply(g, 2);

The above should return 4 of course. Make sure that it does that for you.

So, it's actually much simpler than what you were trying. Actually, the above is more elaborate than what's needed. It can be shortened to

Apply:= (f,x)-> f(x);

Let me take apart your code line-by-line to explain some of its problems. Your first line includes the parameter declaration var:= x. This makes var the parameter, and its default value will be x if SomeProc is called with less than two arguments. This is not a problem per se; I just want you to be aware of what your syntax means to Maple.

Your next line is f(var):= g. This is a very dangerous line, and I'm surprised and disappointed that the syntax checker allows it. It's dangerous because it allows the modification of something, f, which is external to the procedure and hasn't even been declared.

The command LargeExpressions:-Veil is pretty much tailor made for your problem. See ?LargeExpressions.

You mentioned using Newton's method, so I assume that you can do some iteration in Excel. I am about to present a very simple single-line procedure which when iterated exactly three times in hardware-float arithmetic will compute W(x) with a maximum absolute error of 6.7e-13  over your interval of interest, 0.3..30. The procedure is so simple that the only constants that it contains are 1, 2, and 3 (nary a decimal point!), and the only operations are 6 additions, 5 multiplications, 4 assignments, 3 divisions, and 1 natural logarithm. So that's only 57 operations over 3 iterations.

Here it is:

W:= (x,y)-> ((X,g)-> (S-> x*(1-(S+g/2)/(X*(S/g+1)+g/3)))(X^2))(x+1, x+ln(x/y));

"Whoa!" you say, "What do you mean by 'iterate' when the procedure takes more than one argument?"

It's technically called "folding"; I mean computing W(W(W(x,x),x),x). (I don't recommend that you evaluate that symbolically because the result is massive, as is typical with Newton's methods.) You can verify what I said about the absolute error with this plot:

plot(evalhf@((x-> W(W(W(x,x),x),x)) - LambertW), .3..30);

And this plot shows that the relative error is less than 2.7e-13 in hardware-float arithmetic:

plot(evalhf@((x-> W(W(W(x,x),x),x))/LambertW - 1), .3..30);

Now you say, "How did you derive that?" To which I say "with meticulous hand-craftsmanship, decades of practice, and a little help from Wikipedia and Maple." I'll give more details later. In the meantime, browse the Wikipedia article "Householder's method". My procedure W is a third-order Householder iteration. The fact that it's iterated three times is a coincidence.

I doubt that Excel will accept my terse functional-programming style, which isn't required for the efficiency. It's trivial to deconstruct, but I'll do it for you anyway:

W:= proc(x,y)
local g, X, S;
   g:= x+ln(x/y);
   X:= x+1;
   S:= X^2;
   return x*(1-(S+g/2)/(X*(S/g+1)+g/3))
end proc:

I agree that the infomation for each constant should be stored better. I think Record would be the best format. But even as it is, you do not need to know the position of derive in the sequence (which could vary anyway) in order to extract it. You can do

eval(derive, select(type, [ScientificConstants:-GetConstant(e)], `=`))

This will also act gracefully in cases where there is no derive.

AFAIK, circular derivations are not used in the package, and I'm glad for that.

Builtin procedures sometimes act weirdly when used in "naked" form, i.e., without (immediate) arguments. The most common case for me is eval. You can work around this by putting the builtin in a wrapper:

curry((e,i)-> `?[]`(e,i), f)([1]); 

Try this

funcCoeffList:= subsindets[flat](
   (Funcs:= [select(
      hastype, 
      indets[flat](
         subsindets[flat](expression, function, f-> One*f),
         `*`
      ), 
      function
   )[]]),
   {function, identical(One)}, 1
);
Funcs:= eval(Funcs/~funcCoeffList, One= 1);

In case it's not clear, I make One a dummy coefficient of all functions so that it'll be extracted properly, then I set its value to numeric 1.

Yes, Maple has a command nearly tailor-made for this job: Iterator:-TopologicalSorts (I don't know where that name comes from (and I've taken four graduate courses in topology)).  Here is the procedure code, and I'll provide some commentary after it. Warning, this is quite dense code, but it's super fast.

L:= proc(c::posint, L::And(list, satisfies(L-> nops(L)>=2*c)))
local p:= L[], A:= p[1+c..-c-1], P:= (p[1..c], p[-c..-1]); 
     [seq(
          [P[[seq(p[1..c])]], A, P[[seq(p[-c..-1])]]],
          p= Iterator:-TopologicalSorts(
               2*c, rhs~(op(3, rtable((1..c)$2, (i,j)-> i<c+j, 'storage'= 'triangular'['upper'])))
          )
     )]
end proc:

Examples of usage:

L(3, [a1,a2,a3, a,a, b1,b2,b3]);
       (output omitted)

nops(%); #Count entries of previous output.
       57

Many of your input parameters are redundant (which is dangerous in computer code). All that's needed is the value of c and the list of objects to permute that has at least 2*c entries. The list length m isn't needed. The as in the center are not needed, nor their count, and they don't even need to be the same. Also b is not needed. All that information comes from the input list itself.

Here's a speed test:

Perms:= CodeTools:-Usage( #command that measures time and memory usage
     L(6, [a||(1..6), b||(1..6)])
): #Very important! End this with a colon because the output is huge!
memory used=2.10GiB, alloc change=495.45MiB, cpu time=22.97s, real time=12.52s, gc time=13.78s
nops(Perms);
                            1566813

So, that's more than 1.5 million length-12 permutations in less than13 seconds.
 

Since it's such a common Question, it bears repeating that any number of plots with the same number of dimensions (2 or 3) can be combined with plots:-display as Kitonum has shown. However, for this common situation of showing a curve along with some related points, I prefer to do it with a single plot command, like this

X:= <0, 3.38, 5.21, 6.97, 8.4108, 10.099, 10.9232, 11.8091>:
Y:= <0, 0.760e-1, .1275, .1994, .2286, .3222, .3637, .999>:
h1:= solve(Vdc = 0.15e-2*sqrt(2.53669508e8*u^3-6.06101011e8*u^2+3.46343435e8*u), u):
nC:= nops([h1]): #number of branches
plot(
   [h1, <X|Y>], Vdc= 0..11.5, thickness= 0, color= [magenta$nC, blue], 
   style= [line$nC, point], symbol= asterisk, symbolsize= 24
);