I do not have Maple 7, but I have  Maple 2015 - 2017. But in them (in geom3d package) I did not find a simple way to specify an arbitrary polyhedron, except for a few special cases. Therefore, I propose to use plottools package, in which each face of the polyhedron is given by polygon command. Using  plots:-display command you get a polyhedron bounded by these faces. First, in this way I set the faces of the tetrahedron, at the base of which the regular triangle with side 1, the height I took equal to 3/2. Then using the Pythagorean theorem 2 times, I calculated the height of the pyramid with a square base (also with side 1), which has the same lateral edges as the tetrahedron. Then I set its faces in the same way. plots:-display command allows you to depict everything on a single plot.

with(plots): with(plottools): 
ht:=3/2:
Vt:=[[0,0,0], [1,0,0], [1/2,sqrt(3)/2, 0], [1/2, sqrt(3)/6, ht]]:  # The vertices of the tetrahedron
t1, t2, t3, t4:=polygon(Vt[1..3], color=brown), polygon(subsop(3=NULL,Vt), color=blue), polygon(subsop(2=NULL,Vt), color=yellow), polygon(Vt[2..4], color=green):  # All faces of the tetrahedron
st:=evalf(sqrt((sqrt(3)/3)^2+ht^2));  # Length of lateral edge of the tetrahedron
hs:=evalf(sqrt(st^2-(1/sqrt(2))^2));  # Height of a pyramid with the square base
Vs:=[[0,0,0], [-1,0,0], [-1,1,0], [0,1,0], [-1/2, 1/2, hs]]:  # The vertices of the square pyramid
s1, s2, s3, s4, s5:=polygon(Vs[1..4], color=brown), polygon(subsop(2=NULL,3=NULL,Vs), color=green), polygon(subsop(3=NULL,4=NULL,Vs), color=blue), polygon(subsop(1=NULL,4=NULL,Vs), color=yellow), polygon(subsop(1=NULL,2=NULL,Vs), color=red):  # # All faces of the square pyramid
display(t1, t2, t3, t4, s1, s2, s3, s4, s5, scaling=constrained, axes=normal, view=[-1.35..1.35,0..1.35,0..1.55], orientation=[70,75], lightmodel=light1);

                      


Edit. Fixed errors in the formulas for st and hs.

Your equation has an obvious conststant solution  y=ln(2) . If you want to get some other solution, you can solve it numerically with the specified initial conditions, as in the example below:

ode := diff(y(x), x, x)-exp(y(x))+2 = 0;
ics:=y(0)=0.6, D(y)(0)=0;
Sol:=dsolve({ode,ics}, numeric);
plots:-odeplot(Sol, [x,y(x)], x=0..5,  view=-10..2, color=red);

PS. A solution in the form of a series usually gives a good approximation only in sufficient proximity to the starting point. Compare the two solutions:

ode := diff(y(x), x, x)-exp(y(x))+2 = 0:
ics:=y(0)=0.6, D(y)(0)=0:
Sol1:=dsolve({ode,ics}, numeric):
A:=plots:-odeplot(Sol1, [x,y(x)], x=0..5, color=red):
Sol2:=dsolve({ode,ics}, y(x), 'series');
B:=plot(eval(y(x), Sol2), x=0..5, color=blue):
plots:-display(A, B, view=-5..2);

Edit.
 

If I understand correctly  op(0, e)  returns not the type but a procedure (a command, a function) by which this expression  e  is constructed. Sometimes this coincides with a type, but not always. See:

op(0, n!);
op(0, sin(x));
                                                 factorial
                                                    sin

whattype(n!);
whattype(sin(x));
                                                 function
                                                 function

 Integer  is the command in Maple. See help on this. I do not know why this command is needed. I never used it.                                  

Here are all the options to plot a rectangle:

restart;
A:=plots:-display(plottools:-rectangle([0.2, 0.8], [0.8, 0.2], style=line, color=green, thickness=3)):  # Only the border
B:=plots:-display(plottools:-rectangle([0.2, 0.8], [0.8, 0.2], color=green, linestyle=2)):  # Border in the desired style
C:=plots:-inequal({x>0.2, x<0.8, y>0.2, y<0.8}, x=0..1, y=0..1, color=green, nolines):  # Without any border
E:=plots:-display(plottools:-rectangle([0.2, 0.8], [0.8, 0.2]), style=line,  color=red, thickness=3):  # Border in the desired color
plots:-display(< <A | B >, <C | plots:-display(E,C)> >, view=[0..1,0..1], scaling=constrained, axes=none, size=[250,250]);  # "Array" of all the variants

acer's and vv's ways (S  and  s  procedures) don't work explicitly for  n>50  (don't return an explicit integer):

  S(100);
  s(100);
                                          100*exp(1)*GAMMA(100, 1)
                                          100*exp(1)*GAMMA(100, 1)

Here is another way:

r := rsolve({u(n)=n*(u(n-1)+1), u(0)=0}, u):
SS:=s->simplify(value(eval(r, n=s))):
SS(100);

253686955560127297415270748212280220445147578566298142232775185987449253908386446518940485425152049793267407732328003493609513499849694176709764490323163992000

Of course, it's easy to combine with the previous results, so that the explicit result is returned only at the user's request (for those who are annoyed by long numbers).

If we assume that the cylinder (the rolled carpet) axis coincides with a spatial diagonal of the box, then the same solution can be made substantially shorter as a procedure. The main idea: the cylinder touches the box wall with which the spatial diagonal forms the smallest angle. Formal parameters of the procedure: L - list of the dimensions of the box, r - radius of the cylinder. The procedure returns the maximum length of the rolled carpet that is fitted inside the box.

MaxLength:=proc(L::list, r)
local a, b, c, d, k, s;
a, b, c:=sort(L)[];
d:=sqrt(a^2+b^2+c^2);
k:=a/sqrt(b^2+c^2);
s:=r/k;
evalf(d-2*s);
end proc:

Examples of use:

MaxLength([118, 58, 38],10);
                                                       67.66287620

fsolve(MaxLength([118, 58, 38], x)=120);  # The inverse problem
                                                       2.437052665

The procedure  New_list  helps to do this:

New_list:=(L::list, x::realcons)->[seq(`if`(x<l, 2*l, l), l=L)]:

Example of use: 

L:=[seq(rand(1..20)(), i=1..10)];
New_list(L, 7);
                                  L := [14, 18, 2, 18, 4, 17, 18, 10, 5, 1]

                                     [28, 36, 2, 36, 4, 34, 36, 20, 5, 1]

Addition. If the list items must be doubled, which is less than  x , then simply change the inequality sign:

New_list:=(L::list, x::realcons)->[seq(`if`(l<x, 2*l, l), l=L)]:

Edit.

radnormal(-3*(sqrt(10+4*sqrt(5)))+sqrt(50+20*sqrt(5))+sqrt(10-4*sqrt(5))+sqrt(50-20*sqrt(5)));

                                                                               0

If you want to prove this manually, then here are 2 ways.

First way:

Transfer the term with a coefficient of -3 to the other side of the identity:

sqrt(50+20*sqrt(5))+sqrt(10-4*sqrt(5))+sqrt(50-20*sqrt(5))=3*(sqrt(10+4*sqrt(5)));

The left and right parts are positive numbers. Now you can square both parts. Further simplifications are obvious.

Second way:

From the second and fourth terms we take out the factors, then we group the roots in pairs. We obtain two summands. Then the multipliers before the roots are added under the roots. It remains to simplify expressions under the roots:   

Edit.

All plottings and calculations are performed in polar coordinates.

The first example:

F:=z=x^2+y^2:
solve({eval(F, [y=0,z=8]), x>0});
R:=eval(x, %);
plot3d([r*cos(t),r*sin(t),r^2], r=0..R, t=0..2*Pi, scaling=constrained);
Area:=int(eval(sqrt(1+diff(rhs(F),x)^2+diff(rhs(F),y)^2)*r, [x=r*cos(t),y=r*sin(t)]), [r=0..R, t=0..2*Pi]);

The second example (in order to this body to be clearly visible one quarter of the hemisphere is carved out):

restart;
F:=z = (x^2 + y^2)^(1/2): 
G:=x^2 + y^2 + z^2 = 16:
solve({eval([F, G],y=0)[], x>0});
R:=eval(x, %);
plots:-display(plot3d([r*cos(t),r*sin(t),r], r=0..R, t=0..2*Pi, style=surface, color=khaki), plot3d([4*cos(t)*sin(s),4*sin(t)*sin(s),4*cos(s)], t=Pi/2..2*Pi, s=0..Pi/2, color=khaki),  plot3d([4*cos(t)*sin(s),4*sin(t)*sin(s),4*cos(s)], t=0..2*Pi, s=0..Pi/4, style=surface, color=khaki), scaling=constrained, axes=normal, view=[-4.7..4.7, -4.7..4.7, 0..4.7]);
Volume:=int(eval(sqrt(16-x^2-y^2),[x=r*cos(t),y=r*sin(t)])*r-eval((x^2 + y^2)^(1/2),[x=r*cos(t),y=r*sin(t)])*r, [r=0..R, t=0..2*Pi]);

You wrote "...as the legth of this interval is 1, there is two integers in this interval"  It is wrong. The number of solutions can be from one to two. The function  f  returns the set of solutions for each a :

restart;
solve(abs(a+2*i)<=1, i);
L:=[lhs(%[2,1]),rhs(%[1,1])];
f:=x->eval({floor(L[2]),ceil(L[1])}, a=x);
 

Examples of use:

f(1), f(0.5), f(3);

                                   {-1, 0},  {0},  {-2, -1}

The best and fastest way to enter a matrix is using the angle brackets from the keyboard. Try this and you will forget your palettes.

Example:

A:=<1, 2, 3; 4, 5, 6; 7, 8, 9>;

If you still want to use this expression, then just remove the back quotes around it:

Matrix(2, 3, {(1, 1) = Typesetting:-msub(Typesetting:-mi("n"), Typesetting:-mrow(Typesetting:-mn("1"), Typesetting:-mo(","), Typesetting:-mn("1"))), (1, 2) = Typesetting:-msub(Typesetting:-mi("n"), Typesetting:-mrow(Typesetting:-mn("1"), Typesetting:-mo(","), Typesetting:-mn("2"))), (1, 3) = Typesetting:-msub(Typesetting:-mi("n"), Typesetting:-mrow(Typesetting:-mn("1"), Typesetting:-mo(","), Typesetting:-mn("3"))), (2, 1) = Typesetting:-msub(Typesetting:-mi("n"), Typesetting:-mrow(Typesetting:-mn("2"), Typesetting:-mo(","), Typesetting:-mn("1"))), (2, 2) = Typesetting:-msub(Typesetting:-mi("n"), Typesetting:-mrow(Typesetting:-mn("2"), Typesetting:-mo(","), Typesetting:-mn("2"))), (2, 3) = Typesetting:-msub(Typesetting:-mi("n"), Typesetting:-mrow(Typesetting:-mn("2"), Typesetting:-mo(","), Typesetting:-mn("3")))});

Edit.

Go to the Help page  ?dsolve . There you will find a link to the Help page  dsolve, education , on which there are examples of a number of commands that will give you useful information about your equation. For example, you can find out the type of equation:

 restart;
ode:=diff(y(x),x)=sqrt(1-x^2)/y(x)^2;
DEtools:-odeadvisor(ode);

                                                  ode := diff(y(x), x) = sqrt(-x^2+1)/y(x)^2

                                                                       [_separable]

Using Heron's formula, it is easy to construct a system of equations connecting the lengths of the sides with the lengths of altitudes. Maple easily solves this system and we get formulas expressing the sides of the triangle through its altitudes (a triangle can be an arbitrary triangle - it may not be acute):

p:=(a+b+c)/2: S:=sqrt(p*(p-a)*(p-b)*(p-c)):
Sol:=[solve({hA=2*S/a, hB=2*S/b, hC=2*S/c}, {a,b,c}, explicit)][2];

simplify(eval(Sol, [hA=sqrt(3)/2, hB=sqrt(3)/2, hC=sqrt(3)/2]));  # Examples of use
simplify(eval(Sol, [hA=3, hB=4, hC=12/5]));
simplify(eval(Sol, [hA=2, hB=3, hC=3]));

In the first example, we got an equilateral triangle (I took the equal  lengths of altitudes), in the second example - a rectangular triangle (I took the altitudes in the Pythagorean triangle with sides 3, 4, 5), in the third example we got an obtuse triangle.

To avoid problems with names, the graphical structure for each triangle is created using plots:-display command. The animation was not as fast, on each triangle is given 5 frames:

restart;
with(plots):
with(geometry):
for a from 4 to 6 do
b:=a+1: c:=a+2:
angA:=arccos((a+5)/2/(a+2));
point(A,0,0): point(B,c*cos(angA),c*sin(angA)): point(C,a+1,0):
t1||a:=display(draw(triangle(T1||a,[A,B,C]))):
end do:
 
display(seq(t1||i$5, i=4..6), insequence=true, scaling=constrained);

Addition. When working with graphics, in particular when creating animations, I almost never use  geometry:-draw  command. I prefer  plot  command, as well as the various commands from  plots  and  plottools  packages. Also, instead of a  for loop, it is more convenient to use  seq  command. Here is the same animation with this technique with the addition of labels for the vertices of triangles:

restart;
with(plots):
angA:=arccos((a+5)/2/(a+2)): 
V:=[seq([[0,0],[(a+2)*cos(angA),(a+2)*sin(angA)],[a+1,0],[0,0]], a=4..6)]: NV:=[A,B,C]:
T:=seq(textplot([seq([V[i,j][],NV[j]], j=1..3)], font=[times,bold,18], align={above,left}),i=1..3):
display(seq(display(T[i],plot(V[i], color=red))$5,i=1..3), insequence, scaling=constrained);

                
Edit.

restart;
evalindets(u(i,j)^2*v(i,j)+u(i-1,j), 'specfunc(u)', f->subsop(1 = op(1,f)+1, f));
                                               u(i +1, j)^2*v(i, j) + u(i, j)                                       

or

restart;
evalindets(u[i,j]^2*v[i,j]+u[i-1,j], 'specindex(u)', f->subsop(1 = op(1,f)+1, f));
                                             u[i+1, j]^2*v[i, j]+u[i, j]

Edit.