Maybe you want it:

A:=Matrix(2, 2, {(1, 1) = (1/6)*sqrt(3)+(1/2)*I, (1, 2) = (1/6)*sqrt(3)-(1/2)*I, (2, 1) = (1/6)*sqrt(3)-(1/2)*I, (2, 2) = (1/6)*sqrt(3)+(1/2)*I});
A1:=(sqrt(3)/6)%*(Matrix(2,(i,j)->Re(A[i,j]))/(sqrt(3)/6))%+((1/2)%*(Matrix(2,(i,j)->A[i,j]-Re(A[i,j]))/(1/2)));   # Or
A2:=(sqrt(3)/6)%*(Matrix(2,(i,j)->Re(A[i,j]))/(sqrt(3)/6))%+((1/2)%*(Matrix(2,(i,j)->Im(A[i,j]))/(1/2))%*I);

If you want to calculate the values of an expression for different values of a variable, it is convenient to immediately specify it as a function, and not as an expression. Then in the future it will be very easy to calculate the value of this function symbolically or numerically.

Your example:
f := x->2*sin(x)-x^2/10:
f(5),  f(5.);
                                    2*sin(5)-5/2,  -4.417848549                    

Here are 3 ways to solve the problem by an example.

The problem: find the smallest root of equation  sin(100*x)+0.1*x=0  in the range  1 .. 2 . In this range there are several dozen roots.

The first way by  fsolve  command. Using  plot  command, we reduce the range to a single root:

restart;
eq:=sin(100*x)+0.1*x=0:
fsolve(sin(100*x)+0.1*x=0, x = 1 .. 2);
plot(lhs(eq), x = 1 .. 2);
plot(lhs(eq), x = 1 .. 1.1);
fsolve(lhs(eq), x = 1 .. 1.02);

The second way by  RootFinding:-Analytic  command, which provides all the real and complex roots in the range  re=1..2, im=-1..1  and we select the smallest real one:

[RootFinding:-Analytic(lhs(eq), re=1..2, im=-1..1)];
min(%);

The third way by  Student:-Calculus1:-Roots  command, which provides all the real roots in this range and we select the smallest one:

evalf(Student:-Calculus1:-Roots(lhs(eq), x = 1 .. 2));
min(%);



 

One solution we can get by  fsolve  command:

fsolve({x^2 + y^2 + z^2 = 4, x + y + z = 0, x*sin(y*z) = -1});
                                 {x = -1.630958259, y = 0.8860579420, z =0 .7449003167}

The plots shows that there are 2 real solutions (I eliminated z-variable from the system):
plots:-implicitplot([x^2 + y^2 + (-x-y)^2 = 4, x*sin(y*(-x-y)) = -1], x=-3..3, y=-3..3, color=[red,blue], gridrefine=5):
plots:-implicitplot([x^2 + y^2 + (-x-y)^2 = 4, x*sin(y*(-x-y)) = -1], x=-1.7..-1.2, y=0.5..1, color=[red,blue], gridrefine=5):
plots:-display(<%% | %>);

The second solution:
fsolve({x^2 + y^2 + z^2 = 4, x + y + z = 0, x*sin(y*z) = -1}, {x=-1.7..-1.6,y=0.7..0.8, z=0..1});
                    {x = -1.630958259, y = 0.7449003167, z =0.8860579420}

                
Edit.

Use  RealDomain:-solve  instead of  solve . See the corrected file:

mathias1.mw 

You did not write how many values of  a  we should consider, so I took  a=1  and  a=2  only:

seq([fsolve(x^3 - a*x + 1 = 0)], a=1..2);  # All the real roots
plot([seq(x^3 - a*x + 1, a=1..2)], x=-2..2, -2..3, color=[red,blue]);  # The corresponding plots 

Because of your 2d-mode, we can not determine which commands you use. 1d mode explains the difference:

sum(beta__i, i = 0 .. 2); 
sum(beta[i], i = 0 .. 2);

de := 100*diff(v(t),t)=100*9.81-k*v(t);
inc:=v(0)=0;
sol:=dsolve({de, inc}, v(t));
k=solve(limit(rhs(sol), t=infinity)=5, k) assuming k>0;
evalf(%);
plot([5, eval(rhs(sol),%)], t=0..5, linestyle=[3,1], thickness=[0,2], color=[black,red], labels=[t,v(t)], labelfont=[times,bold, 16], size=[800,300]);

        
     
We see that for  t>3  the speed  v(t)  is practically equal to  5 . 

Edit.                   

Q:=[]:
for i from 1 by 1 to 3 do 
simplfloat:=rand(-1.0..1.0): 
a:=simplfloat():     
eq:=a+i:
Q:=[op(Q),[i,a,eq]]:
end do:
writedata("all.dat", Q):
 

We can make a simple visual qualitative analysis of this equation if we rewrite it in the form  c=f(x)  (vv's nice idea) . Using the graphical interpretation, it is easy to see how many real roots this equation has, depending on the value of the parameter c . The numbers  x0  and  y0  are the coordinates of the red point.

c=expand(solve(10*c*x^7-6*x^3+6=0, c));
f:=rhs(%):
solve({diff(f,x)=0, x>0}, explicit):
x0:=eval(x, %);
y0:=eval(f, x=x0);
plots:-display(plot(f, x=-3..3, -0.1..0.2, color=blue, discont, size=[500,350], labels=[x,c]), plot([[x0,y0]], style=point, color=red, symbol=solidcircle, symbolsize=12));

Based on the properties of the function  x->3/(5*x^4)-3/(5*x^7)  and its plot, the following conclusions can be drawn:

1.If  c<0, there is the unique positive solution that tends to  1  if c tends to 0 .
2. If  c=0, there is the unique solution  x=1 .
3. If  0<c<36*14^(2/3)*(1/1715), there are 3 solutions. If  c  tends to  0, then one of them tends to  -infinity, the second to 1, and the third to  + infinity. 
4. If  с=36*14^(2/3)*(1/1715), there are 2 solutions: x=(1/2)*14^(1/3)  and the unique solution <0 .
5. If  c<36*14^(2/3)*(1/1715)  there are  the unique  solution  <0 .

The following simple procedure numerically finds a sorted list of all real roots for any given value c .

RealRoots:=proc(c)
uses RealDomain;
sort([evalf(solve(10*c*x^7-6*x^3+6=0))]);
end proc:

Examples of use:

 seq(RealRoots(c), c=-0.2..0.2, 0.01);  # All the real roots in the range -0.2..0.2 when  c  changes with the step 0.01
map(s->select(`>`, s, 0), [%]);  # Only the positive roots

See corrected file   command_syntax_new.mw

Edit.

If your code is left as it is, then for large values  j  huge symbolic expressions will accumulate, which very slows down the work. So I instead of  rand()  wrote  rand(0...evalf(2*Pi)) , which gives a random angle in the range  0..2*Pi. I also took  c=1. All this greatly accelerated the work of the code. For visualization, I took the first 2 points and 100 steps. With  randomize()  command each time we run this code, we can see a different picture.

 restart;
M:=10: N:=1000: c:=1: r:=rand(0...evalf(2*Pi)):
randomize():
for i to M do
X[i, 0] := 0;
Y[i, 0] := 0;
X[i, 1] := 1;
Y[i, 1] := 0;
for j from 2 to N do
Vinkel := r();
X[i, j] := X[i, j-1]+cos(Vinkel);
Y[i, j]:= Y[i, j-1]+sin(Vinkel)
end do:
end do:
plot([[seq([X[1,j],Y[1,j]], j=0..100)], [seq([X[2,j],Y[2,j]], j=0..100)]], color=[red,blue], scaling=constrained, axes=box);  # Plotting for two points
j:='j':
plots:-animate(plot,[[['seq'([X[1,j],Y[1,j]], j=0..round(s))], ['seq'([X[2,j],Y[2,j]], j=0..round(s))]], color=[red,blue], scaling=constrained, axes=box], s=0..100, frames=200, paraminfo=false);  # Animation for two points 

The plot for  j=0..1000 :

 plot([[seq([X[1,j],Y[1,j]], j=0..1000)], [seq([X[2,j],Y[2,j]], j=0..1000)]], color=[red,blue], scaling=constrained, axes=box);

Edit.

Examples:

sol1 := solve({x+y+z=1, x-y+3*z=7});
sol2 := solve({x+s*y+z=1, x-y+3*z=7*t}, {x,y,z});
Params1:=indets(map(rhs, sol1), name);
Params2:=indets(map(rhs, sol2), name);

Here's the first option as I first understood the problem:

restart;
M0:=100: R3:=rand(1..3):
X[1] := 1: Y[1] := 0:
randomize():
for i from 2 to M0 do
r := R3();
if r = 1 then X[i] := X[i-1]; Y[i] := Y[i-1]+1;
elif r = 2 then X[i] := X[i-1]+1; Y[i] := Y[i-1];
else X[i] := X[i-1]-1; Y[i] := Y[i-1];
end if;
end do:

A:=plots:-animate(plot,[['seq'([X[j], Y[j]],  j = 1 .. round(n))], linestyle=2, color=blue, scaling = CONSTRAINED], n=1..M0, frames=180,paraminfo=false):
B:=plots:-animate(plots:-display,['plottools:-disk'([X[round(n)], Y[round(n)]],0.2, color=red)], n=1..M0, frames=180, paraminfo=false):
plots:-display(A, B, axes=box);

        

Here is another option (probably more correct, see tomleslie's answer and my comment to it):

restart;
M0:=100: R3:=rand(1..3):
X[1] := 0: Y[1] := 0: X[2] := 1: Y[2] := 0:
randomize():
for i from 3 to M0 do
r := R3();
if r=1 then X[i]:=2*X[i-1]-X[i-2]; Y[i]:=2*Y[i-1]-Y[i-2];
elif r=2 then X[i]:=X[i-1]+Y[i-1]-Y[i-2]; Y[i]:=Y[i-1]-X[i-1]+X[i-2];
else X[i]:=X[i-1]-Y[i-1]+Y[i-2]; Y[i]:=Y[i-1]+X[i-1]-X[i-2];
end if;
end do:

A:=plots:-animate(plot,[['seq'([X[j], Y[j]],  j = 1 .. round(n))], linestyle=2, color=blue], n=1..M0, frames=180,paraminfo=false):
B:=plots:-animate(plots:-display,['plottools:-disk'([X[round(n)], Y[round(n)]],0.2, color=red)], n=1..M0, frames=180, paraminfo=false):
plots:-display(A, B, axes=box, view=[min(seq(X[i],i=1..M0))-1..max(seq(X[i],i=1..M0))+1,min(seq(Y[i],i=1..M0))-1..max(seq(Y[i],i=1..M0))+1], scaling = constrained);

     

Edit.

Random_walk.mw

It's true for any real  n :

expr:=exp(-n*ln(2*Pi));
simplify(expr, exp) assuming n::real;

# Or
RealDomain:-simplify(expr, exp);