@Kitonum   Continuation of the previous code:

Pic := plots[display](A, Pr, seq((Tr@@k)(F), k = 0 .. 7), C, T, axes = none):

Tr1 := (p, alpha)-> plottools[rotate](p, alpha, [[0, 0, 0], [0, 0, 1]]) :

plots[animate](plots[display], ['Tr1'(Pic, alpha)], alpha = 0 .. 2*Pi, frames = 50);

@lg674  Should be  Pi  instead of pi in your code. Also you have 4 constants  E, A, B, m  instead of 3 costants a, b, c, so the texts of the procedures should be changed.

@Carl Love  Very clear and easy solution, vote up . The same but slightly shorter and with the display of the desired results.

X:= [u,v,w,alpha,beta,gamma]:  G:=LinearAlgebra:-GenerateMatrix: 

M, Sys:= G(convert(Sys, list), diff~(X(t), t$2)):  'M'=M;

C, Sys:= G(-convert(Sys, list), diff~(X(t), t)):  'C'=C;

K, F:= G(-convert(Sys, list), X(t)):  'K'=K; 'F'=F;

@Markiyan Hirnyk  Please explain how it could happen that in the post on January 10 2016, all comments have been removed, including mine, this post has turned into a question with your answer in the form my animation?

You wrote in brackets something like this

restart;

x=solve(y, x);

                                    x=()

@vv  Why to check, we can not instead of

eval(  'evalf'( 'add'('add'(a, i=1..N), j=1..N)),  N=1000);

simply write

add(add(evalf(a, i=1..1000), j=1..1000);

The result is the same, but it requires 10 times less time. 

@Majmaj You can't use  labels = ["x values", "y values"] if axis[2] = [color = white], because everything that refers to  y-axis  will be white.

Markiyan Hirnyk  Why do you never specify the source of the problem? This is very similar to the problems of student competitions.

@tomleslie  for the adjustment. I did not know before about this option.

Unfortunately, some old bugs remain. Here is the example of the definite integral with assuming:

A:=int(abs(cos(n*x)), x=0..Pi);

A assuming n=0;

eval(A) assuming n=0;

value(A) assuming n=0;

evalf(A) assuming n=0;

simplify(A) assuming n=0; 

Of course the integral is easy to calculate by the direct substitution:

eval(A, n=0);

I have put  +1 for OP,  but my vote up inexplicably disappeared.

Addition: Now I see that it appeared again.

Good work, vote up. Just one remark. Your code is too complex for this problem.  Roughly the same can be done much more simple means:

A := plottools[polygon]([[sqrt(3)+1, 0, 0], [1, 1, 0], [1, -1, 0]], color = "Violet"):

B := plottools[polygon]([[sqrt(3)+1, 0, 0.01], [1, 1, 0.01], [1, -1, 0.01]], color = "Green"):

C := plots[display](A, B):

F := t->plottools[rotate](t, alpha, [[0, 0, 0], [0, 0, 1]]):

G := t->plottools[rotate](t, alpha, [[0, 0, 0], [0, 1, 0]]):

plots[animate](plots[display], [(G@F)(C)], alpha = 0 .. 2*Pi, scaling = constrained, axes = none, frames = 180, lightmodel = none);

 The idea of the code - the triangle at the same time rotates around the y-axis and z-axis.

Addition. A slight modification of the code clearly shows that the rotations in 3d around different axes do not commute:

A := plottools[polygon]([[sqrt(3)+1, 0, 0], [1, 1, 0], [1, -1, 0]], color = "Violet"):

B := plottools[polygon]([[sqrt(3)+1, 0, 0.01], [1, 1, 0.01], [1, -1, 0.01]], color = "Green"):

C := plots[display](A, B):

F := t->plottools[rotate](t, alpha, [[0, 0, 0], [0, 0, 1]]):

G := t->plottools[rotate](t, alpha, [[0, 0, 0], [0, 1, 0]]): 

Anim1:=plots[animate](plots[display], [(G@F)(C)], alpha = 0 .. 2*Pi, scaling = constrained, axes = none, frames = 180, lightmodel = none):

Anim2:=plots[animate](plots[display], [(F@G)(C)], alpha = 0 .. 2*Pi, scaling = constrained, axes = none, frames = 180, lightmodel = none):

plots[display](<Anim1 | Anim2>);

@kountouretti  For visualization I used my code in the answer below  "Direct calculation". I took  epsilon = 0.1. We find that for all  n>=287  all the points  (n, x(n))  lie within the grey band. To the point of not merged into one line I took  n  with the step = 10 :

restart;

x := n->n^2/(n^2+31*n+228):

x0 := limit(x(n), n = infinity):

N := epsilon->max(1, floor(op(1, solve(abs(x(n)-x0) <= epsilon, n)[-1]))+1):

 N(0.1);

                                                      287

A := plot([1-.1, 1+.1, [seq([n, x(n)], n = 1 .. 500, 10)]], n = -2 .. 500, -0.1 .. 1.25, style = [line, line, point], color = [black, black, red], linestyle = [3, 3], symbol = solidcircle, symbolsize = 7, labels = [n, "x(n)"], tickmarks = [default, [0.9, 1, 1.1]]):

B := plots[inequal]({y > 0.9, y < 1.1}, x = 0 .. 500, y = 0 .. 1.25, optionsfeasible = [color = "LightGrey"]):

plots[display](A, B);

Edited. 

@AdVen1935   See the corrected file. To upload a file click the thick green up-arrow in the menu of the editor.

Histogram1.mw

@AdVen1935  See attached file. I was a little changed procedure code and your example code.

Histogram.mw