What do you mean by "multiple roots"? A sequence can only converge to one value. Yes, an iterative root finder can be used to find multiple roots, but only if given multiple initial values, and there'll be separate convergent sequences for each root.

We'll define the absolute error of the kth term of a sequence as the absolute value of the difference between this term and the limit. If a sequence converges, we'll say that the order or convergence is r if the absolute error of the kth term of the sequence can be approximated by A^(r^k) for some 0 < A < 1 and r > 0. The values of A and r can be estimated by linear regression. Here are some Maple procedures for doing that:

restart
:
(*>>>-------------------------------------------------------------------------
This procedure approximates the order of convergence of a convergent sequence
S. It outputs 3 values: the order of convergence, the error-approximating 
function, and the coefficient of determination (aka r-squared) of the 
regression.

This syntax requires Maple 2018 or later.
-------------------------------------------------------------------------<<<*)
ConvergenceOrder:= proc(S::~Vector);
local Err:= (-ln@abs)~(S[..-2] -~ S[-1]), n:= numelems(Err), k, j, R, r;
   Digits:= 15;

   #Check integrity of sequence:
   for k to n do until Err[k] > 0;
   if k > n-2 then
      if k > n then 
         error "sequence doesn't converge sufficiently rapidly"
      fi;
      error "need more terms in sequence"
   fi;
   for j from k to n-1 do 
      if Err[j] >= Err[j+1] then 
         error "absolute errors must be strictly decreasing" 
      fi
   od;

   R:= Statistics:-LinearFit(
      [1, _k], <($k..n)>, ln~(Err), _k, 
      ':-output'= [':-rsquared', ':-parametervector']
   );
   Digits:= 5;
   evalf([(r:= exp(R[2][2])), _e = exp(-exp(R[2][1]))^(r^_k), R[1]])[]
end proc
:
(*>>>----------------------------------------------------------------------
This procedure takes a function f whose root we want and a template T for a
root-finding method that uses only one initial value and returns the 
iteration function.
-----------------------------------------------------------------------<<<*)
MakeIterator:= proc(f::procedure, T::procedure)
local x;
   unapply(simplify(T(x,f)), x)
end proc
:
(*>>>----------------------------------------------------------------------
This procedure applies the iteration J, starting with initial value v0, and
returns the sequence of iterates as an Array.

This uses Maple 2019 syntax.
-----------------------------------------------------------------------<<<*)
Sequence:= proc(
   J::procedure, #unary iteration function 
   x0::complexcons, #initial value
   {  #keyword parameters:
      abserr::And(realcons, positive):= 10.^(1-Digits), #convergence criterion
      max_iters::posint:= 99 #max numbers of iterations
   }
)
local x:= x0, x1, X:= Array([x]), n;
   for n to max_iters while abs((x1:= evalf(J(x))) - x) > abserr do
      X,= (x:= x1)
   od;
   if n > max_iters then error "did not converge" fi;
   X
end proc
: 

The iterative method that you gave is close to---but not quite the same as---Newton's method. If you provide a function and initial value for which your method converges, I'll apply these procedures to it. In lieu of that, I'll just show an example using regular Newton's method:

#Example usage:
Digits:= 999: #Large Digits is needed to get sufficient data.

M:= (x,f)-> (x - f(x)/D(f)(x)): #method--Newton's method in this case.
f:= x-> x^3 - x - 2: #function to find root of
x0:= 1: #initial guess of root

ConvergenceOrder(Sequence(MakeIterator(f,M), x0));
                                /      _k\         
                                \2.0391  /         
            2.0391, _e = 0.77243          , 0.99721

So, the order of convergence is about 2.04, the error-approximating function is as shown, and the coefficient of determination is > 99.7%, indicating that the regression was highly accurate.

I'll give you two methods---the first is easy and lazy and the second guarantees uniform selection from the available primes (meaning each prime is equally likely to be chosen):

EasyRandPrime:= (R::range(realcons))->
   nextprime(rand(ceil(lhs(R))-1 .. prevprime(floor(rhs(R)))-1)())
:
#This syntax requires Maple 2018 or later:
UniformRandPrime:= proc(R::range(realcons))
local rnd:= rand(ceil(lhs(R))..floor(rhs(R))), r;
   do until isprime((r:= rnd()));
   r
end proc:

The usage is the same for both:

EasyRandPrime(3..30);
UniformRandPrime(3..30);

If you require a large number of random selections from the same range, the above procedures can be made more efficient fairly easily.

The Answer by dharr exploits the low-level structure of a procedure in a way that won't work if f's parameters have type decalrations and is potentially dangerous if f's parameters have default values.

The Answer by Christian Wolinski does something that could be done more simply with the commands curry and/or rcurry.

What you want can be achieved much more robustly by simply declaring the procedure f such that its parameters have default values that are global names:

f:= proc(a:= ':-a', b:= ':-b', c:= ':-c') a*b*c end proc:
f(3);
                             3*b*c

Another option, which is even more robust are eliminates the correspondence between parameters and the order that they're declared--but which requires a little more typing--is to use keyword parameters: 

f:= proc({a::algebraic:= ':-a', b::algebraic:= ':-b', c::algebraic:= ':-c'})
   a*b*c 
end proc:
f(b=7);
                     7*a*c

Now f's parameters can be given in any order.
 

The values of keyword parameters of type name and the keywords themselves (the word output in this case) are always global. When there is a conflict between a local (your local Q) and a global (your desired Q setting of the keyword) with obstensibly the same name, the local name takes precedent. In those cases, it's necessary to distinguish the global with the :- prefix, thus, ':-Q'.

I give this Answer just to explain the reasoning behind the other Answer.

Some Maple programmers always use the ':-...' syntax for these unevaluated global names.

It's as simple as 

diff(BesselI(alpha, x), x);

Note that the last letter of the function name is uppercase I ("eye"), not lowercase l ("ell").

I can't get Maple to give the Fourier transform of your function, even by direct integration in the a=1 case.

However, regarding your other question, about algsubs: As discussed on its help page ?algsubs, this command treats denominators differently than non-denominators. To substitute x=ex (where x is a name and ex is an expression) into another expression LE, there are only two commands that should be used: either eval or subs---reserve algsubs for more-complicated substitutions, such as when x is an expression.

• If you want to make the substitution regardless of whether it's mathematically valid, use subs(x= ex, LE).
• If you want Maple to consider the mathematical validity, use eval(LE, x= ex).

Either of these commands will be more robust and far more efficient than algsubs.

For this reason among others, SearchText should be deprecated. StringTools:-Search and StringTools:-SearchAll are externally compiled, so they're adequately efficient replacements; probably more efficient when multiple searches are specified in a single command. Making a wrapper for a case-insensitive version of these to replace searchtext is trivial (if we say that only ASCII characters can have a "case").

While I agree with your premise in theory, in practice you may need to make a distinction between different kinds of infinity. The code that you show returns +infinity, somewhat by happenstance. It could just as well have been -infinity or some complex flavor of infinity. Also note that once you allow infinities, the usual field rules of real and complex arithmetic no longer apply. In general, Maple's ability to navigate these issues on its own is pretty good.

As Acer's Answer shows, mod is complicated, and so I usually reserve its use for "higher" algebraic and number-theoretic computations. For numeric arithmetic, I use irem, which automatically returns unevaluated for non-numeric arguments--suitable for your purpose.

iremp:= (x,m::posint)-> irem(m+irem(x,m),m):
plot(iremp(ceil(x),2), x= -9..9);

The first of these commands is to make the correct adjustment for negative arguments, similar to the distinction between modp and mods.  

A not-necessarily-integer parameter---nominally called "degrees of freedom"---is commonly used for the t and chi-square distributions. See "Welch's t test" and "Welch-Satterthwaite equation". This is taught even in second-semester statistics. Perhaps "degrees of freedom" wasn't the best choice for the name of this parameter, but this usage is unfortunately well established, and the mathematical concept is solid.

Unfortunately (AFAIK) there's no way (either with or without assumptions) that you can use the commands solve or isolve to achieve your goal, which I think is to extract 100, 20, 3 from 123. But there is an easy way:

b:= 10:  convert(123, base, b) *~ [seq(b^i, i= 0..floor(log[b](123)))]  

The 10 could be replaced by any integer base greater than 1. The results are returned with the least-significant digits on the left.

In order to implement the Path psuedocode that you show, it is necessary that the AllPairsDistance procedure record the path information in the matrix Next (note that lowercase next is a reserved word in Maple). So, I wrote an object-oriented module that implements both procedures.

This uses Maple 2019 syntax. If you're using an earlier Maple, let me know, and I can change it.
 

Download FloydWarshall.mw

It just assumes that all unweighted edges have weight 1. That is explained in the first paragraph of the help page ?GraphTheory,AllPairsDistance.

I'm not sure if that answers your Question, since I can't find a direct reference to "Betweeness centrality matrix".

I think that you've conflated the mostly typographical notion of the first or left-most coefficient with the algebraic concept of the leading coefficient. And it seems that you simply want the sign of the left-most. I think that this procedure will do it:

Sign:= (e::algebraic)-> 
   if e::{`+`,`*`} then thisproc(op(1,e))
   elif e::{`^`,function} then 1
   else try sign(e) catch "invalid argument": thisproc(op(1,e)) end try
   fi
:

Edit: Added type `^` to code above.

Your

= "inverse of tan (x)"

contains an invisible bad character. Just backspace over this part of the command and retype it manually; don't copy and paste.