The reason is that sqrt(200) is the standard deviation of the discrete uniform distribution on {$1..49}.

Slightly more generally, suppose you're taking random combinations of k numbers out of a finite set S (with at least k elements).  Let X[i,j] be the i'th number in the j'th random combination, for i = 1 .. k and j = 1 .. n.  Then the sample variance is

1/(k*n) * sum(sum(X[i,j]^2, j=1..n), i=1..k) - (1/(k*n) * sum(sum(X[i,j], j=1..n),i=1..k))^2

For convenience, I'm not using the version with Bessel's correction, but that makes no difference in the limit n -> infinity.  Also for convenience, I'll suppose the members of each combination are presented in random order, so that each X[i,j] has a uniform distribution on S.  Using the Strong Law of Large Numbers, as n -> infinity, 1/(k*n) * sum(sum(X[i,j], j=1..n), i=1..k) and 1/(k*n) * sum(sum(X[i,j]^2, j=1..n),i=1..k) converge with probability 1 to E[X] and E[X^2] respectively, where X has a uniform distribution on S, and so your sample variance converges to the variance of X.

The limit as n -> infinity should be sqrt(200) = 10*sqrt(2).  But it's late, so I'll show how to get this result tomorrow...

@wschrabi: It is a discontinuity in cm(x,t) (at least for x near 0), but I don't see how it could be a pole.

@wschrabi: It is a discontinuity in cm(x,t) (at least for x near 0), but I don't see how it could be a pole.

The only thing really wrong in your code is the line

c[i] := [seq(op(a[i]), i = 1 .. sets)]: 

First of all, you shouldn't really use the same variable i for two different purposes in the same command.  In this case it's not so bad because of the special evaluation rules for seq, but it is confusing, and it means that you can't access the local variable i of the procedure inside the seq: it would be better to use, say, j for the index variable of the seq.  But note that at this point you have only defined a[1] to a[i], not all the way to a[sets]: when a[j] is unassigned, op(a[j]) is j.  The result is extra elements i+1, i+2, .. sets in c[i].  You can notice this if you look closely at the histograms on the right in the first few frames.  What you want is

c[i]:= [seq(op(a[j], j = 1 .. i)]:

I don't know what you mean by "intergrate from the one side only".  I only used the PDE's and initial and boundary conditions that were given.  I did need to specify the timestep, because with default values the interval 0 < t < t0 would be missed completely.

As for Mathematica, this is the wrong site for questions about it.

I don't know what you mean by "intergrate from the one side only".  I only used the PDE's and initial and boundary conditions that were given.  I did need to specify the timestep, because with default values the interval 0 < t < t0 would be missed completely.

As for Mathematica, this is the wrong site for questions about it.

I wouldn't say that it is not invertible, rather that (assuming the parameters are all positive) the inverse is not defined on the whole real line, only on the range of the function, which is the interval (0, g+h). The function is increasing, and therefore one-to-one, under these assumptions.

Given real values of the parameters, the best way to actually calculate the inverse function is using fsolve.

I wouldn't say that it is not invertible, rather that (assuming the parameters are all positive) the inverse is not defined on the whole real line, only on the range of the function, which is the interval (0, g+h). The function is increasing, and therefore one-to-one, under these assumptions.

Given real values of the parameters, the best way to actually calculate the inverse function is using fsolve.

What Matrix are you taking eigenvalues of?  Do you have any reason to believe these eigenvalues should be real?
If your Matrix is symmetric, you should construct it as a symmetric Matrix, e.g.

Matrix(..., shape=symmetric);

otherwise Maple's numerical methods may produce very small imaginary parts of some eigenvalues. 

What Matrix are you taking eigenvalues of?  Do you have any reason to believe these eigenvalues should be real?
If your Matrix is symmetric, you should construct it as a symmetric Matrix, e.g.

Matrix(..., shape=symmetric);

otherwise Maple's numerical methods may produce very small imaginary parts of some eigenvalues. 

There are several issues here.  The various facilities for handling names (print, printf, lprint, convert(..., string), ...) handle trailing tildes in different ways.

showassumed affects print, but not the others.  It doesn't matter whether or not there actually are assumptions on the variable, showassumed controls what print does with it.  The help file ?interface could be clearer on these points.

printf (and its relatives) and lprint always remove trailing tildes, as far as I can see.  convert(..., string) always preserves them.

There are several issues here.  The various facilities for handling names (print, printf, lprint, convert(..., string), ...) handle trailing tildes in different ways.

showassumed affects print, but not the others.  It doesn't matter whether or not there actually are assumptions on the variable, showassumed controls what print does with it.  The help file ?interface could be clearer on these points.

printf (and its relatives) and lprint always remove trailing tildes, as far as I can see.  convert(..., string) always preserves them.

The "variable defined in Maple" was A[sl].  The "Equation pasted from MathType" had A[s*l].  It's rather strange, because there doesn't appear to be any space between the s and the l, but you can see it by copying and pasting into a 1D Maple input region.

The "variable defined in Maple" was A[sl].  The "Equation pasted from MathType" had A[s*l].  It's rather strange, because there doesn't appear to be any space between the s and the l, but you can see it by copying and pasting into a 1D Maple input region.