There's a somewhat-builtin facility for detecting submatrices that are all zeros: ArrayTools:-IsZero.

RemoveZeroRowsAndColumns:= (A::Matrix)->
   A[remove(i-> ArrayTools:-IsZero(A[i,..]), [$1..op([1,1], A)]),
     remove(j-> ArrayTools:-IsZero(A[..,j]), [$1..op([1,2], A)])
    ]
:

I haven't done a time test, but I suspect that this'll be significantly faster than the other Answers on your 1000x1000 matrix.

Suppose that an Array A is declared and initialized with a statement such as

A:= Array(1..9, 1..9);

Now A has two dimensions, usually called rows and columns. In the language of your error message, the number of dimensions is called the rank. (This has nothing to do with the concept of rank of a matrix or linear transformation from linear algebra.) The elements of A can be accessed as A[i,j] where i and j are integers in the appropriate ranges, in this example 1 to 9. The i and j are called indices. If you were to attempt an access of A as A[5,3,7], then there'd be three indices, and the number of indices would exceed the rank, like the error message says. You have something akin to this in your code.

Generally you should post your code. If the code is more than, say, 100 lines, and you are able to isolate the erroneous part, just post that. If you can't isolate it, just post the whole thing. Please post code as either plaintext or a Maple worksheet file. Screenshots or cut-and-paste of 2D input is not helpful.

Your first part is incorrect. Since the two distributions being added are identical, your statement is equivalent (as far as Maple is concerned) to

w:= 2*RandomVariable(BetaDistribution(1, m-1));

and thus the two distributions are not independent. (This is due to Maple's automatic simplification.) To get around the automatic simplification and achieve an i.i.d. (independent identically distributed) sum of RVs (Random Variables), you need to do something like

(X,Y):= 'Statistics:-RandomVariable(BetaDistribution(1, m-1))' $ 2; #Note carefully the quotes.
W:= X + Y;

To achieve your second goal, the simplest imaginable thing works:

X2:= W/2;

I can confirm the bug that you're encountering.

You can use unapply instead of makeproc:

unapply(Myfun, [a,b]);

Assuming that you're restricting the analysis to real numbers, ln(dy/dx) = 0 simply implies that dy/dx = 1. Using the ln just introduces unnecessary complexity.

On the other hand, perhaps you are recalling something that you learned about the derivative of the logarithm rather than the logarithm of the derivative? This is a useful concept.

The command op doesn't "get the variable" in an expression. Please read ?op. To get the variable, assuming that there's only one, use

martin:= indets(f, And(name, Not(constant)))[];

plot(Statistics:-PDF(Weibull(6.2, 2), x), x= 0..26);

Using Statistics:-EmpiricalDistribution, a procedure that generates your parameterized random variable can be written in one line:

MyRV:= (i::posint, x::posint)->
   Statistics:-RandomVariable(EmpiricalDistribution(irem~([$i-x..i-1, $i+1..i+x]  -~ 1, 100) +~ 1))
:
 Example:   
X:= MyRV(98,5):
map(trunc, Statistics:-Sample(X,25));

      [93, 2, 97, 93, 95, 2, 96, 3, 96, 96, 97, 94, 1, 3, 100, 2, 99, 95, 100, 1, 94, 2, 100, 1, 97]

Statistics:-Mean(X);

      68

If your goal is simply to draw samples from the distribution rather than to analyze it, then a much more efficient procedure can be used. On the other hand, if your goal is to express the distribution's statistics (such as Mean) in terms of the parameters i and x, it is significantly more complicated, but still I think doable. Either way, let me know.

You are trying to create a Matrix with datatype= float[8]. Such a Matrix is not allowed to contain variables which have not been assigned real numeric values. It looks like you may have two distinct "Delta bar" variables, which look confusingly similar. I need to see a worksheet to say for sure.

The convert command is unnecessary. Once all the variables have been assigned numeric values, you can say

M:= Matrix(K, datatype= float[8]);

You did not directly call evalhf, but it is used in the conversion to float[8], so it was called indirectly.

Maple's rules of operator precedence are a bit weird for mod (see ?operators,precedence): It has different precedence on its left and right sides! You no doubt think that a + b mod c is parsed as a + (b mod c), but it's actually parsed as (a + b) mod c. Furthermore, you can't tell whether it's using mods or modp (see ?mod). Finally, mod is a high-level symbolic operator which probably shouldn't be used in numeric code. If you want the integer remainder, then use the command that's designed for that: irem.

sieve[i+j*j]:= 2 + irem(sieve[i+j*j], 2);

No, but there's enough prime infrastructure built-in such that writing such a procedure is trivial. And here it is:

PrimesBetween:= proc(a::posint, b::posint)
local L:= table(), p:= nextprime(a-1);
   while p <= b do
      L[p]:= [][];
      p:= nextprime(p)
   end do;
   {indices(L, 'nolist')}
end proc:
   
PrimesBetween(100, 200);

     {101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199}

By the way, Maple commands do not generally "print" their output; they return it. Whether the output is printed and where it is printed is determined by a variety of means under the user's control.

Perhaps you don't realize that Psi is an already-defined special function? See ?Psi. If this is the case, then you should use a different name than Psi. But if you really have a good reason for wanting to use the name Psi, then first execute

local Psi;

Since subs works, you should use it.

The set of cases where algsubs works is not a superset of the set of cases where subs works. Thus, "It is a generalization of the subs command" should perhaps be removed from its help. Rather, think of algsubs as something that will sometimes work for cases where subs fails. Please read carefully the paragraph of ?algsubs that deals with how negative powers of the substitution target are handled. That is what is happening in the case that you present.

Algsubs is a high-level command, written in Maple. It consists of 194 lines (*footnote 1) of highly recursive (*footnote 2) code split over 11 procedures. On the other hand, subs is a low-level builtin command, and it is easy to fully understand what it does, even though you can't read the code. The job of replacing one name with another is clearly a job for subs.

(*footnote 1): Lines are as counted by showstat, which doesn't count proc, local, elif, else, end, blank lines, or comments as separate lines.

(*footnote 2): On a quick read through, I counted 13 points where the code calls itself recursively.

If A is the array, and it has n elements, then the animation can be produced by

plots:-display([seq(A[k], k= 1..n)], insequence);

It is not your fault. This is a bug. The issue is that there are two ways to enter an inert double integral (using 1D input), but only one of those ways is understood by dchange. The two ways are

Int(Int(f(x,y), x= a..b), y= c..d);

That's the old way, which is still 100% acceptable. The newer way is

Int(f(x,y), [x= a..b, y= c..d]);

but this isn't understood by dchange and some other commands.

The workaround is to use value(UP) instead of simply UP in the dchange command.