Like this:
 

Download Markov_chain.mw

P:= <0.6199, 0.1450, 0.0636, 0.1549, 0.0166;
     0.6550, 0.0002, 0.0870, 0.1827, 0.0751;
     0.2990, 0.1294, 0.0010, 0.4925, 0.0781;  
     0.2773, 0.7184, 0.0043, 0,      0;
     0.0009, 0.6229, 0.3762, 0,      0
>:
V:= [S, I__1, I__2, I__3, R]:
n:= nops(V):
GT:= GraphTheory:  LA:= LinearAlgebra:
MC:= GT:-Graph(
    V,
    {seq(seq([[V[i],V[j]], P[i,j]], i= 1..n), j= 1..n)}
);
GT:-DrawGraph(MC);
#Find long-term probability of being in each state:
V =~ LA:-LinearSolve(
    <P-Matrix(n, shape= identity) | <seq(1, 1..n)>>^+,
    <seq(0, 1..n), 1>
);

[S = 0.535156546280454, I__1 = 0.215609042590436, I__2 = 0.0648612328168834, I__3 = 0.154231678262430, R = 0.0301415000497959]

Like this:

plot3d(
    [[r, theta, r], [r, theta, 0]], r= 0..1, theta= 0..Pi, 
    coords= cylindrical, scaling= constrained
);

The 3D version of polar coordinates is coords= cylindrical. Making the 3rd coordinate 0 in the 2nd plot restricts it to the xy-plane, so you get the region of integration.

Explicit returns of true or false are almost always redundant. The following is equivalent to VV's procedure, but in boolean-operator form:

restart:
TypeTools:-AddType(
    ElementaryMatrix,
    A-> A::'Matrix'('square') and (
            nops((local Z:= op(2, A-1))) > 1 implies
                nops((local j:= [lhs~(Z)[]])) = 2 and
                Z = {j[]=1, j[[2,1]][]=1, (j[1]$2)=-1, (j[2]$2)=-1}
        )
):
M:= <0,1,0,0; 1,0,0,0; 0,0,1,0; 0,0,0,1>:
type(M, ElementaryMatrix);
                              true

I included the identity matrix as an elementary matrix. If this is not desired, change  > 1 to <> 1.

Here's a much shorter code to do the check that you want, where M is the matrix:

remove(
    type, 
    op(2, M), 
    anything= And(complexcons, satisfies(x-> is(x, ComplexRange(-90,90))))
)[];

To get just the numbers, you can do

rhs~(cozum)[]

but they will not be in the same order.

To retain the order, do

rhs~([cozum[]])[]

To get the number associated with a particular label, do

eval(g[-1], cozum)

Like Mariusz, I recommend that you build upon Maple's own laplace command. But unlike Mariusz, I recommend that unevaluated laplace transforms in the result (which will come from transforming unknown functions such as f(x) and its derivatives) be converted to unevaluated Elziki transforms. Like this:

Elziki:= (f::algebraic, t::name, v::name)->
    (expand@evalindets)(        
        v*inttrans:-laplace(f, t, 1/v),
        specfunc(:-laplace),
        L-> 'Elziki'(op(1..2, L), 1/op(3,L))*op(3,L)
    )
:
Elziki(diff(f(x),x$2) - diff(f(x),x) + 2*f(x), x, v);

Thus the use of laplace is transparent to the end user.

The problem is that the summand is evaluated before k has a value. This doesn't happen if you use command add:

add(diff(M[1,1], x[k]) + diff(M[1,k], x[1]) - diff(M[k,1], x[1]), k= 1..2)

When you've seen LinearAlgebra algorithms labeled as FractionFree, that refers to the internal workings of the algorithm, not to its final output. If you have a matrix of rational numbers, the probability that its eigenvectors can be expressed as rational numbers is infinitesimal, so the example that you show is very special. And, as you've noted, it's easy to put a vector into the fraction-free form. Indeed, it's trival:

FF:= v-> (rationalize@expand)~(mul(denom~(v))*v):
M:= LinearAlgebra:-RandomMatrix(3);
(FF@op@op)~(3, LinearAlgebra:-Eigenvectors(M, output= list))[];

The nesting depth to which statement output is printed in the absence of an explicit print statement is controlled by the environment variable printlevel (see ?printlevel). Its default value is 1.

The command parse will interpret a string as if it was Maple code. So, if S is the string you gave, then U:= parse(S) gives the equivalent expression, with units.

Like this:

``(4/5) %* (2*x-5)

Now that I understand what you want, you can do it with inert operators combined with the thisproc and 'procname' from my previous Answer. Inert operators are infix arithmetic operators preceded by %. The code below uses %^ and %*. You may also change + to %+ if you choose. Depending on your Maple version and input mode (1D or 2D), you may need to specify the inert operators in prefix functional form, e.g., `%^`(x, y). The 2D output display should remain as usual, except that multiplication will appear as a gray * (this can be changed if you want).

hy:= proc(n, x, y)
    if n=0 then y + 1
    elif n=1 then x + y
    elif y=0 then 1
    elif y=1 then x
    elif n=2 then x %* y
    elif n=3 then x %^ y
    elif n=4 then 
        if y::posint then x %^ thisproc(4, x, y-1)
        else x %^ 'procname'(4, x, y-1)
        fi
    elif [n,y]::list(posint) then thisproc(n-1, x, thisproc(n, x, y-1))
    else 'procname'(n-1, x, 'procname'(n, x, y-1))
    fi
end proc
:
hy(4, x, 4);
                     %^(x, %^(x, %^(x, x)))

hy(4, 4, y); 
                     %^(4, hy(4, 4, y - 1))

The above should display in Maple as exponentiation usually does in 2D output.

Note that x is a "dummy" variable in your function: The ability to fully recurse does not depend on whether x is symbolic or numeric, so I now just check n and y.

Like this:

H:= proc(n, x, y)
    if n=0 then y+1
    elif n=1 and y=0 then x
    elif n=2 and y=0 then 0
    elif y=0 then 1
    elif not [n,x,y]::[nonnegint$3] then 
        'procname'(n-1, x, 'procname'(n, x, y-1))
    else
        thisproc(n-1, x, thisproc(n, x, y-1))
    fi
end proc
:
H(n, x, y);
                  H(n - 1, x, H(n, x, y - 1))

H(2,x,y);
                    H(1, x, H(2, x, y - 1))

H(2, x, 2);
                      H(1, x, H(2, x, 1))

If you wish, 'procname' could be replaced by 'H' and thisproc could be replaced by H.

The above may be what you said you didn't want. If so, sorry, and please show some explicit examples of what you do want.

A space curve has one-dimensional input and three-dimensional output. So there can be only one free variable in all 3 components. In nm's Answer, the first space curve has free variable x; the second has y.

On the other hand, a surface (such as created by plot3d) has two free variables.

In formal usage in Maple (such as when used in a type declaration), a function is different than a procedure. (Maple's usage is somewhat idiosyncratic on this point.) What you want is a procedure. Like this:

Newton:= proc(f::procedure, x0::complexcons, n::posint)
local Df:= f/D(f), x:= x0; 
    to n do x:= evalf(x - Df(x)) od;
    x
end proc:
Newton(x-> x^2-2, 1, 9);
                          1.414213562

As an example, sin is a procedure and sin(1) and sin(x) are functions. A function is a name followed by an expression sequence of arguments in parentheses.

Both procedures made with proc() ... end proc and the arrow -> have type procedure.

In informal usage, such as on help pages and on this forum, function often does mean procedure, and I think that this is also the more standard computer science usage outside of Maple.

The above example code is not meant to be a good implementation of Newton's method. It's just intended to show you how to declare f.