int(-2*log((1+sqrt(s))/(1-sqrt(s)))/((-s^2+1)*(s-1)*sqrt(s)), s = 0..z) assuming 0 < z,z < 1;

First, get rid of all evalm. It is never needed in modern code, and it does nothing at all in the way that you are using it.

Assign the result of dsolve(..., numeric) to a variable:

Sol:= dsolve(..., numeric);

Do not use the output option for dsolve.

Then the following procedure evaluates the Matrix AA:

AA:= proc(S,t)
local Cu_inv:= eval(Cu, S(t))^(-1), tCu_inv:= eval(tCu, S(t))^(1);
     eval(Avv - Avu.Cu_inv.Cv + tCv.tCu_inv.Auu.Cu_inv.Cv - tCv.tCu_inv.Auv, S(t))
end proc:

To use it, invoke AA(Sol, t) where t is an actual numeric value.

I notice that your Auv and Avu are identically 0. I assume that you were going to change those. If not, two of the terms in the procedure are trivial.

The difference is the order that the operations are performed. If you use the gcd, then the gcd is done first, ignoring the mod. The result of that is then passed to mod. If you use Gcd, then it recognizes that it needs to take the modulus into account while finding the gcd. The same is true for all of the capitalized operators listed at ?mod .

For the vast majority of Maple procedures, arguments are evaluated before being passed: In f(g(x)), g(x) is evaluated before being passed to f. For gcd(x^2+1, x+1), the result is 1 and that is all that `mod` receives. The Gcd avoids this because it does not nothing by itself; it is an "inert" operator.

You can't use square brackets as parentheses---as an outer level of grouping, for example. You can only use round parentheses for algebraic grouping. So, your expression should be like this:

solve({((alpha[1]-alpha[2]*lambda)*sqrt(x)+p[2]*lambda)*(k[1]*(1-lambda^2)+2) = p[1]*(2*k[1]*(1-lambda^2)+2), ((alpha[2]-alpha[1]*lambda)*sqrt(x)+p[1]*lambda)*(k[2]*(1-lambda^2)+2) = p[2]*(2*k[2]*(1- lambda^2)+2)}, [p[1], p[2]]);

P:= Matrix([''[1,4],[1,5],[2,6],[3,7]'' $ 3]):
x:= 3*Statistics:-MovingAverage(P[..,1], 3):
X:= 3*Statistics:-MovingAverage(P[..,1]^~2, 3):
y:= 3*Statistics:-MovingAverage(P[..,2], 3):
Y:= 3*Statistics:-MovingAverage(P[..,2]^~2, 3):

select(has, f, u);

select(has, f, x);

Use this procedure:

MatrixToStrings:= proc(M::Matrix(integer))
uses LA= LinearAlgebra;
local k,x;
     [seq(cat(seq(sprintf("%d",x), x= M[k,..])), k= 1..LA:-RowDimension(M))]
end proc:

Use add, not sum.

add(add(add(a_(i,j,k)*d(i, j, k), k = 0 .. 1), j = 0 .. 1), i = 0 .. 1);

The reason that it doesn't work with sum is that i, j, and k have the value 1 from some previous computation. Using add, it will ignore previous values of the index variables. But sum---like the vast majority of Maple procedures---evaluates its arguments before they are passed.

f:= n-> a^n + b^n + c^n:
Sol:= {solve({f(1)=1, f(2)=2, f(3)=3}, {a,b,c})}:
nops(Sol);
                               1
evala(eval(f(5), Sol[1]));
                               6

Here's a solution technique that tends to give a large number of zeros in the solutions vector, and a large number of integer entries.

I will assume that your 500x1300 matrix contains floating-point data. If it contains integer data, then you can skip the Integerize step below.

The first is convert the problem to an integer problem by multiplying each row by its least common denominator. The procedure Integerize does that.

Download Integerize.mw

This is probably not more efficient or direct, but it does look neater to me and easier to understand:

tmp:= seq(phi[j-1]*~[1-p[j], p[j]], j= 5..2, -1);
Vector(op~([tmp]));

A *~ B

See ?elementwise

It takes several short steps to remove all of Maple's knowledge about gamma. In Maple 16, do these:

restart:
constants:= ({constants} minus {gamma})[]:
`evalf/gamma`:= proc() end proc:
`evalf/constant/gamma`:= proc() end proc:
unprotect(gamma);

Note that this will affect certain special functions, such as Psi, that depend on gamma. The effect will ensue even when those functions are being used (unbeknowst to the user) in the background, rather than being explicitly invoked by the user.

Consider the following plot:

plot([x, 1/3 - (1/3)*x], x= 0..1);

Since the two lines intersect between x=0 and x=1, there is no meaningful way to order the expressions; the is will return false no matter which direction you put the inequality.

sys:= {diff(y(t),t) = 2 - t^2/2 - 1/4 * Y(t), diff(Y(t),t)= y(t), Y(0)=0}:
eval(y(t), dsolve(sys, {Y(t), y(t)}));