If you mean the polar equation for any circle, then your equation  r^2=(r*cos(t)-a)^2+(r*sin(t)-b)^2  is wrong, because at the right  r  is the polar radius for a point on the circle and at the right there is the radius of the circle. Here is the correct calculation with a plotting:

solve((r*cos(t)-a)^2+(r*sin(t)-b)^2=R^2, r);
plot(eval([%], {a=-1,b=0, R=2}), t=0..2*Pi, color=red, coords=polar);

The obvious solution for the general matrix - write it down as a procedure. As the input of this procedure can be submitted any specific  nxn matrix:

ModalMatrix:=proc(A)
uses LinearAlgebra;
simplify(fnormal~(evalf~(Eigenvectors(A)[2])), zero);
end proc:

Examples of use:

ModalMatrix(<-4,-2; 3,1>);
ModalMatrix(<-4,-2,5; -2,1,7; 5,7,9>);
                     

Since starting with n = 3 symbolic solution can be very cumbersome, simplifying the commands used.

Simple  Rand  procedure solves your problem. It returns a random number in a range  -r .. r . r  can be specified in any format: integer, fraction, float, symbolic.

Rand:=proc(r::positive)
local rf, p;
rf:=convert(r, float);
p:=rand(-rf .. rf);
p();
end proc:

Example of use ( 10 random numbers in the range  -Pi .. Pi ):

seq(Rand(Pi), i=1..10);

0.065963526, 0.743846604, 0.980407970, 1.155281063, 0.793998054, -3.075735010, 0.546438190, -1.067155770, 2.994200916, 1.581718683

The procedure was written in Maple 2016. In older versions you may need  the call of  RandomTools[Generate]  command.

A bit shorter:

A := simplify([a=log[12](27), log[36](24)]);
B := simplify(subs(ln(3)=k*ln(2), A));
simplify(eval(B[2], k=solve(B[1], k))); 

plot3d([r*cos(theta), r*sin(theta), 1-r^2*cos(theta)^2], r = 0 .. 1, theta = 0 .. 2*Pi, color = khaki, filled = true, style = surface, numpoints = 40000);

The formal parameters of the procedure  P  are the inital values of the velocity  v1  and the coordinate of the point  x1  and the number of the steps  n :

P:=proc(v1, x1, n)
local F, v, x, m, Δt, i;
global V, X;
F:=rand(-0.5..0.5);
v[1]:=v1; x[1]:=x1; m:=1; Δt:=1;
for i from 2 to n do
v[i] := v[i - 1] +(1/m)*F()*Δt;
x[i] := x[i - 1] + 0.5*(v[i - 1] + v[i]);
od;
V:=convert(v,list);
X:=convert(x,list);
V, X;
end proc: 
 

Example of use:

P(1, 0, 10);

[1, .9794824628, .7573440398, 1.053436739, 1.478762216, .9900980562, .4951472162, .6204418912, 1.056397926, .6999520291],   [0, .9897412315, 1.858154483, 2.763544872, 4.029644350, 5.264074486, 6.006697122, 6.564491676, 7.402911584, 8.281086562]

Use  solve  command instead of  fsolve  for your example. Here a simplier example of the matrix  <x,1; 1,x> , which leads to the solution of the system:

solve({-c1+c2=0, c1-c2=0, c1^2+c2^2=1}, explicit);

     {c1 = (1/2)*sqrt(2), c2 = (1/2)*sqrt(2)}, {c1 = -(1/2)*sqrt(2), c2 = -(1/2)*sqrt(2)}

a:=convert(2.1345, string);
b:=[seq(a[i], i=1..length(a))];
remove(t->t=".", b);
parse~(%);

Another way (shorter):

ListTools:-Reverse(convert(op(1, 2.1345), base, 10));
                                              

Error reason is that the syntax of  VectorCalculus:-int  is different from the  int  command  from  Maple kernel, which you have used. Therefore, do not download the whole package  VectorCalculus, and use the long form of command call as  VectorCalculus:-command

I corrected only the first error in the calculation of the volume. The remaining errors you can easily fix yourself.

See the corrected file.

error_integral1.mw

Slightly reduce the right end of  x - range:

f:=x->piecewise(x>-Pi and x<=Pi, 3*x):
plot(f(x), x=-Pi..Pi, discont=true):
a:=-Pi:  b:=Pi:
p:=b-a:
fp:=x->f(x-floor((x-a)/p)*p):
plot(fp(x), x=-6*Pi..7*Pi-0.001, discont=true, tickmarks=[piticks,piticks]);

Example:

A:=plot3d([cos(phi)*sin(theta), sin(phi)*sin(theta), cos(theta)], phi=-arccos(0.9)-0.01 .. Pi/2-arccos(0.9), theta=arccos(0.9) .. `if`(phi>=-arccos(0.9) and phi<=arccos(0.9), arcsin(0.9/cos(phi)), Pi/2), style=surface, color=khaki, scaling=constrained, numpoints=25000):
B:=plots[display](seq(plottools[rotate](plots[display](A), Pi/2*k, [[0,0,0], [0,0,1]]), k=0..3)):
C:=plottools[reflect](B, [[0,0,0], [1,0,0], [0,1,0]]):
Cube:=plottools[cuboid]([-0.9,-0.9,-0.9], [0.9,0.9,0.9], color=pink, transparency=0.8):
plots[display](B, C, Cube, axes=none);

Edit.
 

For example,  if  i=4  then

`or`(seq(R[i]=R[k], k=1..i-1))   means

R[4]=R[1]  or  R[4]=R[2]  or  R[4]=R[3]

In order to more clearly show that your function is unbounded in the neighborhood of the origin (and the top and bottom), it is helpful to reduce the ranges for x and y axes, and to increase grid and to equate the scales along the axes.

plot3d(2*x/(x^2+y^2), x = -5 .. 5, y = -5 .. 5, style=surface, color=khaki, view=-10..10, grid=[500,500], scaling=constrained, axes=normal);

If I understand the problem, the walk is over, if we find ourselves at a point where we already were in this walk. I completely rewrote your procedure  Randwalk3 . Procedure  RandWalk  describes only one random walk. To specify several walks, we simply repeat this procedure.

RandWalk:=proc(N)
local P, i, R3, r;
P[1]:=[0,0]; P[2]:=[1,0];
R3:=rand(1..3);    # 1: go straight on; 2: turn right; 3: turn left
for i from 3 to N do
r:=R3();
      if r=1 then P[i]:=[P[i-1][1], P[i-1][2]+1]      # go straight on
      elif r=2 then P[i]:=[P[i-1][1]+1, P[i-1][2]]   # turn right
      else P[i]:=[P[i-1][1]-1, P[i-1][2]]                 # turn left
      end if;
      if `or`(seq(P[i]=P[k], k=1..i-1))  then  break end if; 
end do; 
convert(P, list);
end proc:

Example of use:

for j from 1 to 15 do
RandWalk(1500);  L[j]:=nops(%);
print(%%, L[j]);
od:
L:=convert(L, list);  # The list of the lengths of these random walks

Edit.

First you can specify every piece of your space curve by plots[spacecurve]  command then all together by plots[display]  command.

Example:
L1:=plots[spacecurve]([cos(t), sin(t), t], t=0..2*Pi, color=red):
L2:=plots[spacecurve]([cos(t), t-2*Pi, sin(t)+2*Pi], t=2*Pi..4*Pi, color=red):
plots[display](L1, L2, axes=normal, scaling=constrained);

Of course quite easy to write a special procedure that automates this process, if necessary.

Addition.  Here is the procedure which automates the plotting of a piecewise curve in 3D:

PiecewisePlot:=proc(P::piecewise,Opt::list)
local n, m, t, i, L;
uses plots;
n:=nops(P);
if n::odd then error `nops(P) should be even` fi;
t:=op(indets(op(1,P)));
for i from 1 to n/2 do
L[i]:=spacecurve(op(2*i,P),t=op([2*i-1,1,1],P)..op([2*i-1,2,2],P), op(Opt));
od;
display(seq(L[i], i=1..n/2));
end proc:

Example of use:

L:=piecewise(0 <= t and t < 2*Pi,[cos(t), sin(t), t], 2*Pi <= t and t <= 4*Pi, [cos(t), t-2*Pi, sin(t)+2*Pi], t>4*Pi and t<=5*Pi, [1, t-2*Pi, 2*Pi], t>5*Pi and t<=6*Pi, [cos(t),3*Pi,sin(t-5*Pi)+2*Pi]);
PiecewisePlot(L, [color=red, thickness=3, scaling=constrained, axes=normal]);

Edited.