I do not think that there is any programatic trace left of the menu choice: I see that the 'Typesetting' object is identical in either case. This is in line with the design of other features in the "new Maple" which are menu driven, without a visible programatic counterpart.

If this is the result you want:

M := matrix(2,2,[(x^4+x^2+1)^2,z*(f[x]-z)/(f(x)-y),
x*(f[x]-z)/(f(x)-y)^2,(f[x]-z)^2/(f(x)-y)]):
subs(f[x]-z=lambda*(f(x)-y),eval(M));
                [  4    2     2                      ]
                [(x  + x  + 1)          z lambda     ]
                [                                    ]
                [   x lambda             2           ]
                [   --------       lambda  (f(x) - y)]
                [   f(x) - y                         ]

No, it is not possible. Not only kernel but library, Java version and code, other dlls,  and presumably GUI, change from version to version

If this is the result you want:

M := matrix(2,2,[(x^4+x^2+1)^2,z*(f[x]-z)/(f(x)-y),
x*(f[x]-z)/(f(x)-y)^2,(f[x]-z)^2/(f(x)-y)]):
subs(f[x]-z=lambda*(f(x)-y),eval(M));
                [  4    2     2                      ]
                [(x  + x  + 1)          z lambda     ]
                [                                    ]
                [   x lambda             2           ]
                [   --------       lambda  (f(x) - y)]
                [   f(x) - y                         ]

No, it is not possible. Not only kernel but library, Java version and code, other dlls,  and presumably GUI, change from version to version

p:=sin(n*Pi/2):
simplify(subs(n=2*m,p)) assuming m::integer;
                                  0

simplify(subs(n=2*m+1,p)) assuming m::integer;

                                    m
                                (-1)

Yes, under the transformation (x -> -y, y -> -x) , it holds s -> -s, t -> t and the equation remains invariant:

subs([x=-y,y=-x],eq):
subsindets(%,`^`,x->trigsubs(x)[2]);

                                            2         2
              cos(x + y) sin(x - y) = cos(x)  - sin(y)

I wonder whether there is an algorithm that detects these symmetries.

eq:=cos(x+y)*sin(x-y)=cos(x)^2-sin(y)^2:
solve(eq,x);

1/2*Pi-y, -1/2*Pi-y, 
arctan((2*sin(y)^2-1)/(2*sin(y)*cos(y)-1)*(2-4*sin(y)*cos(y))^(1/2),
(2-4*sin(y)*cos(y))^(1/2)), 
arctan(-1/2*(2*sin(y)^2-1)/(2*sin(y)*cos(y)-1)*(2-4*sin(y)*cos(y))^(1/2),
-1/2*(2-4*sin(y)*cos(y))^(1/2))

Sounds like a user friendly way to trace what Mathematica does as here:

integrate 1/(2+cos(x)) dx

By the way, it is interesting that functions work  with round brackets.

And this is a too common term for a google search. So, in this case, I find it simpler to edit directly the URL to:  http://www.maplesoft.com/support/help/view.aspx?path=use

For reference, Maxima 5.18 quickly gives an error here:

residue(1/sqrt(1+x^4),x,x=(1+I)/sqrt(2));
`taylor': attempt to expand 1/sqrt(x^4+1)
at a point depending on x
 -- an error.  To debug this try debugmode(true);

I wonder whether an algorithm exist implementing integration by residues, for integrals like this one. And if so, whether such algorithm has been implemented in some other CAS.

In regards to the transformation to the cosine form, this one liner is a possibility:

J:=Pi*sum((-1)^n*cos(a*Pi*(2*n+1)/2),n=0..infinity):
`/`(op((simplify@(x->x/exp(1/2*I*a*Pi)))~([numer,denom](J))));
                                   Pi
                            1/2 ---------
                                    a Pi
                                cos(----)
                                     2

I wonder whether an algorithm exist implementing integration by residues, for integrals like this one. And if so, whether such algorithm has been implemented in some other CAS.

In regards to the transformation to the cosine form, this one liner is a possibility:

J:=Pi*sum((-1)^n*cos(a*Pi*(2*n+1)/2),n=0..infinity):
`/`(op((simplify@(x->x/exp(1/2*I*a*Pi)))~([numer,denom](J))));
                                   Pi
                            1/2 ---------
                                    a Pi
                                cos(----)
                                     2

It may depend on the case whether splitting the expression into numerator and denominator may be advantageous. E.g. for the expression of above you could generate a list of pairs [numerator,denominator] like:

vfX3 := diff(Rp(X,Y),X)/Rp(X,Y)*(2*Rp(X,Y)-Z)/(Rp(X,Y)-Z)
-diff(Rm(X,Y),X)/Rm(X,Y)*(2*Rm(X,Y)-Z)/(Rm(X,Y)-Z):

map([numer,denom],[op(vfX3)]);
    /d          \
  [[|-- Rp(X, Y)| (2 Rp(X, Y) - Z), Rp(X, Y) (Rp(X, Y) - Z)],
    \dX         /

          /d          \
        [-|-- Rm(X, Y)| (2 Rm(X, Y) - Z), Rm(X, Y) (Rm(X, Y) - Z)]]
          \dX         /

manipulate each term, then recombine each numerator - denominator pair, convert back into sum, etc. From my point of view, the main issue here is again to keep under control the automatic simplification/normalization that may undo or spoil your manipulations at the stage of recombination. The effect of these automatic operations is difficult to predict, in part as they are not well documented. So, I am not very optimistic about the implementation of a procedure of generic usefulness. It may work though if you have to repeat similar calculations with similar expressions.

It may depend on the case whether splitting the expression into numerator and denominator may be advantageous. E.g. for the expression of above you could generate a list of pairs [numerator,denominator] like:

vfX3 := diff(Rp(X,Y),X)/Rp(X,Y)*(2*Rp(X,Y)-Z)/(Rp(X,Y)-Z)
-diff(Rm(X,Y),X)/Rm(X,Y)*(2*Rm(X,Y)-Z)/(Rm(X,Y)-Z):

map([numer,denom],[op(vfX3)]);
    /d          \
  [[|-- Rp(X, Y)| (2 Rp(X, Y) - Z), Rp(X, Y) (Rp(X, Y) - Z)],
    \dX         /

          /d          \
        [-|-- Rm(X, Y)| (2 Rm(X, Y) - Z), Rm(X, Y) (Rm(X, Y) - Z)]]
          \dX         /

manipulate each term, then recombine each numerator - denominator pair, convert back into sum, etc. From my point of view, the main issue here is again to keep under control the automatic simplification/normalization that may undo or spoil your manipulations at the stage of recombination. The effect of these automatic operations is difficult to predict, in part as they are not well documented. So, I am not very optimistic about the implementation of a procedure of generic usefulness. It may work though if you have to repeat similar calculations with similar expressions.

There was a typo in 's' above, where an Rp should have been Rm. Once corrected, vfX above is recovered:

sp1:=ln(Rp/Rm*(Rp-Z)/(Rm-Z)):
sp2:=subs(Rp=Rp(X,Y),Rm=Rm(X,Y),sp1):
s:=diff(Rp(X,Y),X)=(X-1)/Rp(X,Y),diff(Rm(X,Y),X)=
(X+1)/Rm(X,Y), diff(Rp(X,Y),Y)=Y/Rp(X,Y), 
diff(Rm(X,Y),Y)=Y/Rm(X,Y):
vf2X:=diff(sp2,X):
vf3X:=subs(s,vf2X):
fX:=op(indets(vf3X,specop(`*`,`+`))):
vf5X:=map(x->x*vf3X/fX,fX):
subs(X-1=Xp,X+1=Xm,vf5X):
collect(%,[Xp,Xm],simplify):
subs(Xp=X-1,Xm=X+1,Rp(X,Y)=Rp,Rm(X,Y)=Rm,%);

               (2 Rp - Z) (X - 1)   (2 Rm - Z) (X + 1)
               ------------------ - ------------------
                    2                    2
                  Rp  (Rp - Z)         Rm  (Rm - Z)

There was a typo in 's' above, where an Rp should have been Rm. Once corrected, vfX above is recovered:

sp1:=ln(Rp/Rm*(Rp-Z)/(Rm-Z)):
sp2:=subs(Rp=Rp(X,Y),Rm=Rm(X,Y),sp1):
s:=diff(Rp(X,Y),X)=(X-1)/Rp(X,Y),diff(Rm(X,Y),X)=
(X+1)/Rm(X,Y), diff(Rp(X,Y),Y)=Y/Rp(X,Y), 
diff(Rm(X,Y),Y)=Y/Rm(X,Y):
vf2X:=diff(sp2,X):
vf3X:=subs(s,vf2X):
fX:=op(indets(vf3X,specop(`*`,`+`))):
vf5X:=map(x->x*vf3X/fX,fX):
subs(X-1=Xp,X+1=Xm,vf5X):
collect(%,[Xp,Xm],simplify):
subs(Xp=X-1,Xm=X+1,Rp(X,Y)=Rp,Rm(X,Y)=Rm,%);

               (2 Rp - Z) (X - 1)   (2 Rm - Z) (X + 1)
               ------------------ - ------------------
                    2                    2
                  Rp  (Rp - Z)         Rm  (Rm - Z)