You have four specific x values in you IBCs, but Maple can handle at most two. You have -5000, -1000, 1000, 5000.

I don't know why the surd doesn't work, but limit(a^(1/k), k= infinity) does work.

There seems to be a bug when there's only one element in the list being mapped over. The bug is not about a difference between `^` and op; but rather about a difference between [1,2] and [[a,b]]: The former has two elements while the latter has only one. The bug is also manifested by Map(`^`, [2], 2).

Here's a very basic (although structured, functional, and extensible) implementation. I didn't include any error checking of user input. So the user needs to make sure that there's the right number (25 by default) of characters in the key, split doubletons in the plain text, make sure only alphabet characters are in the plaintext or ciphertext, and make sure that texts have an even number of characters. However, I did generalize it so that you are not restricted to (the traditional) 5x5 alphabet: You can use any rectangle (just change R and C in the module locals).

(*
A basic implementation of the Playfair cipher (see "Playfair cipher" in Wikipedia).
There is no error checking of user input, and the user has to split doubletons on
their own.
                                                                                   *)
Playfair:= module()
option `Written 19-Apr-2013 by Carl Love.`;
local
    R:= 5, C:= 5,  #For a 5x5=25 alphabet
    
    #Convert number in 1..25 into row-column coordinates
    to_row_col:= proc(n)  local r,c;  r:= iquo(n-1, C, 'c');  (r+1, c+1)  end proc,

    #Convert row-column coordinates to number in 1..25
    from_row_col:= (r,c)-> C*(r-1)+c,

    rotf:= k-> k mod C + 1,  #mod C, shifted by 1 (rotate forward)
    rotd:= k-> k mod R + 1,  #mod R, shifted by 1 (rotate down)
    
    #For a pair in 1..R*C, return the corresponding Playfair pair. Note that if
    #we just consider the indices in the array, then the mapping from index-pair
    #to index-pair is a constant, independent of the cipher key.
    Playfair_pair:= proc(a,b)
    option remember;
    local r1, r2, c1, c2;
        if a=b then  error "expected unequal arguments"  end if;
        (r1,c1):= to_row_col(a);
        (r2,c2):= to_row_col(b);
        if r1=r2 then  (from_row_col(r1,rotf(c1)), from_row_col(r2,rotf(c2)))
        elif c1=c2 then  (from_row_col(rotd(r1),c1), from_row_col(rotd(r2),c2))
        else (from_row_col(r1,c2), from_row_col(r2,c1))
        end if
    end proc,

    EnDig:= table(),  #Encrypt: Digraph to Digraph
    DeDig:= table()   #Decrypt: Digraph to Digraph
;
export
    SetKey:= proc(Key::string)
    local i, j, i2, j2, DigPlain, DigCrypt;
         #Encrypt and decrypt every possible digraph.
        for i to R*C do  for j to R*C do
            if i=j then  next  end if;
            (i2,j2):= Playfair_pair(i,j);
            DigPlain:= cat(Key[i],Key[j]);
            DigCrypt:= cat(Key[i2],Key[j2]);
            EnDig[DigPlain]:= DigCrypt;
            DeDig[DigCrypt]:= DigPlain
        end do  end do;
        [][]
    end proc,

    Encrypt:= proc(plaintext::string)
    local k;
        cat(seq(EnDig[plaintext[2*k-1..2*k]], k= 1..iquo(length(plaintext),2)))
    end proc,
    
    Decrypt:= proc(ciphertext::string)
    local k;
       cat(seq(DeDig[ciphertext[2*k-1..2*k]], k= 1..iquo(length(ciphertext),2)))
    end proc
;
end module;

Download Playfair.mw

The package linalg, the command evalm, and the operators &. and &* (when used in the linalg/matrix context) should not be used... ever. These commands are ancient and deprecated. Even experts are forgetting the nuances of their usage. If you have a professor giving you code to work with, please tell that person that I said to not use these commands... ever.

Maple's solve does not return anything when there are no solutions. I agree that that is confusing. In this case, some trivial algebra "by hand" shows that there is indeed no solution. The -2s and rs cancel and you're left with

abs(1+sqrt(1+4*r)) < 1

On the left, you have 1 plus something nonnegative, so the inequality is never true.

The assuming has nothing to do with this: There are no solutions no matter what you assume.

A Matrix can be an element of an Array, just like anything else. A Matrix can even be an entry in another Matrix. There's nothing special that you need to do. Just try it like you were writing above,

MatrixArray:= Array(1..9);

for i to 9 do MatrixArray[i]:= V(i) end do;

I'm curious as to why you begin variable names with &. Usually these are only used for binary and unary operators. But if that's what you want to use, I guess that there's no harm in it. If you are trying to emulate C syntax, then I think that it is a bad idea: the & symbol has meaning in both languages, but that meaning is completely different.

Let's say that your three procedures are named s1, s2, and s3. And let's suppose that you wish for your users to refer to these methods as "QR", "LU", and "RREF" respectively. In the module below, I made dummies for the three procedures, with just a print statement in each. You'll put your actual procedures there, and you'll want to have actual parameter declarations such as proc(A::Matrix, B::Vector, ...) where I have proc().

Solve:= module()
local
    s1:= proc() print("s1", args) end proc,    
    s2:= proc() print("s2", args) end proc,
    s3:= proc() print("s3", args) end proc,
    
    Methods:= table(["QR"= s1, "LU"= s2, "RREF"= s3]),

    ModuleApply:= proc({method::identical("QR", "RREF", "LU"):= "LU"})
        Methods[method](_rest)
    end proc
;
end module;

Now the users will call this via

Solve(A,B, ...other arguments..., method= "QR");

If no method argument is given, the method will default to "LU" (that's what the :="LU" does). The position of the method argument in the calling sequence does not matter. All arguments other than the method= ... are passed on to the solving procedure (that's what the _rest does). If the user gives a method other than the three listed, an appropriate error message is given. Note that this module does not have, and does not need, any exports; and it does not have, nor need, any main body (other than the NULL body provided by that isolated semicolon).

If you have Maple 17, check out ?SignalProcessing,Engine, ?SignalProcessing,Engine,FFT, etc.

It's not clear what you mean by "radius from 0..1". The radius of a cylinder is constant. I'll assume that you mean that the radius is 1. If the radius is 1, and the volume is Pi/2, then the height is 1/2.

The basic plot of the cylinder is

plot3d([1,theta,z], theta= 0..2*Pi, z= 0..1/2, coords= cylindrical);

Note the part in square brackets: [1,theta,z]. That gives the [r, theta, z] cylindrical coordinates, but, as I said above, r is constant 1 in this case.

If you want to get fancy and put endcaps on your cylinder,

plots:-display([
    plot3d([1,theta,z], theta= 0..2*Pi, z= 0..1/2, coords= cylindrical), #Same as above
    plot3d(
        [[r,theta,0], [r,theta,1/2]], r= 0..1, theta= 0..2*Pi, coords= cylindrical,
        thickness= 2, grid= [10,25], transparency= .35
   )
]);

From the File menu, select Print Preview.

Why not just write your own version?

Mod:= (a,m)-> a mod m + `if`(m < 0, m, 0):

Your range of e is -500..500. Over this range, 1/cosh(e) is effectively zero for the vast majority of the range. (By "effectively", I mean within the constraints of hardware floating-point arithmetic.) I think that you may have mixed up some units or missed some powers of 10 when you formulated the problem and got that range for e.

Or, if this is a homework problem, could it be that the professor has given you a trick question? You say that you want a plot of U(e,0). This can be deduced immediately from the initial conditions with no PDE solving: U(e,0) = 1/cosh(e).

I wrote the code to do what you want. Please use it and test it. There are two slight differences:

1. Maple's range operator `..` is non-associative, and hence there's no way to write (a..b)..c without the parentheses. I don't think that it would be possible for Maplesoft to make it associative without affecting existing code because it works with lhs and rhs, just like `=`, and hence it inherently has two operands (either or both of which may be NULL).
2. My implementation requires inserting &.. between the rtable and the index specification if you want to use the new "skip" (i.e. "every nth") index specification that you requested. Unlike (1) above, this is not a "hard" requirement. It is possible to implement this new indexing without using any new operator between the rtable and the square brackets. But it is a little more complicated. One would need to ?overload the `?[]` operator.

Any standard-Maple index specification, no matter how complicated, should "pass through" my &.. operator and work the way it used to, so please test this aspect.

(*
 An implementation of "skip" or "every nth" rtable indexing: An attempt
 to implement some of the more-sophisticated Matlab-style indexing.
                                                                       *)
`&..`:= module()
option `Written 15-Apr-2013 by Carl Love.`;
local
     #Dimension bounds, set as each dimension is processed
    lo :: integer, hi :: integer,

     #Process the new type of index spec that this module introduces:
    #(a..b)..skip (also handles the case (..)..skip). Returns a list of the
    #actual indices to use.
    Skip:= proc(IndSpec :: range..posint)
        local a, b, k;
        a:= lhs(IndSpec);
        if a = (..) then (a,b):= (lo,hi) else (a,b):= op(a) end if;
        seq(k, k= `if`(a < 0 and lo = 1, hi+a+1, a)..`if`(b < 0 and lo = 1, hi+b+1, b), rhs(IndSpec))
    end proc,

    #Handle the dimension/index-spec pair for 1 dimension.
    DimProcess:= proc(Ind :: specfunc(anything, ``))
        local
            ind:= ``(op(2,Ind)),  #Index spec
            ret  #return value
        ;
        (lo,hi):= op(op(1,Ind));  #Bounds for this dimension
        ret:= expand(subsindets(ind, range..posint, Skip));
        if ind ::``(range..posint) then  [ret]  else  ret  end if     
    end proc,

     #Main code for &.. operator. The index spec for each dimension is paired
     #with that dimension's bounds.
    ModuleApply:= (A :: rtable, Inds :: list)->
        A[(DimProcess ~ (zip(``, [rtable_dims](A), Inds)))[]]
;
end module:     

IndexSkip.mw

The symbol e is not pre-defined in Maple for the purpose of input, although it does appear in Maple's output. You need to replace e^(-x) with exp(x) and replace e^(-10*x) with exp(-10*x). If you also solved the ODE with Maple, you'll need to resolve it using exp(-x) instead of e^(-x). Even though it may appear that Maple has given you a solution, that solution is just treating your input e as if it were any other variable.

Regarding your second plot, it is not clear what is meant by "plotting an equation", especially when that equation has only one variable. You could move the 300 to the left side of the equation, and then just plot the left side. The basic plot command does not handle equations.