Before plotting all constants should be specified. Of cause, you can build several graphs at one plot.

Example:

bX:=-1:  aX:=5:

plot([seq(seq((X_0-aX)/a+(X_0-bX)/b, a=1..3), b=1..3)], X_0=bX..aX, thickness=2, color=[red,blue,brown, yellow,pink,green,gold,cyan,grey]);

Another example (a, b are constants, aX, bX are changed):

restart;

a:=-1:  b:=5:

assign(seq(A[n]=plot((X_0-n)/a+(X_0-(n+5))/b, X_0=n..n+5, thickness=2, color=[red,blue,brown, yellow,pink][n]), n=1..5)):

plots[display](seq(A[n], n=1..5));

Just check correctness of these rules for each example.

Solution of (a):

T:=(a,b,c,d)->Matrix([[a,a^2,a^3],[b,c,d]]):

M1:=T(a1+a2,b1+b2,c1+c2,d1+d2);

M2:=T(a1,b1,c1,d1)+T(a2,b2,c2,d2);

LinearAlgebra[Equal](M1,M2);

The plotting depends on the choice of contours and locations for text. If the number of contours is equal to  n, then I divide the range ​​of the function on  n  parts and choose the values ​​of the function at the midpoints of these parts.

restart;

f := (x, y)->x^2+y^2-5;

m, M := minimize(f(x, 0), x = 0 .. 10),  maximize(f(x, 0), x = 0 .. 10);

n, d := 5, (M-m)/n;

C := seq(m+(1/2)*d+d*(i-1), i = 1 .. n);

P := plots[contourplot](f, -10 .. 10, -10 .. 10, color=brown, contours = [C]):

T := plots[textplot]([seq([solve(f(x, x) = C[i])[1]$2, C[i]], i = 1 .. n)], font = [TIMES, ROMAN, 14], align = right):

plots[display](P, T);

 See  Statistics[Regression]  help page. 

I could not find the options for the exact solution. But the exact solution can be easily obtained by the simple formulas:

with(Student[LinearAlgebra]):

a:=<-2|3|2>;  b:=<7|-3|-4>;

Pr:=(a.b)/(Norm(b)^2)*b;   # exact projection

evalf(Pr);   # approximate projection

c:=a-Pr;  # exact orthogonal complement

evalf(c);   # approximate orthogonal complement

Norm(c);  # exact norm of orthogonal complement

evalf(%);   # approximate norm of orthogonal complement

 Addition:  See http://en.wikipedia.org/wiki/Vector_projection

Write as follows:

g := (1/8)*x*sin(13/x)+arcsin(5*x^2/(8*x^2+1));

G := unapply(g, x);

H := unapply(G(G(x)), x);

D(H)(Pi/2);

In calculating the derivative of a function at a point convenient to use the differential operator  D  command instead of  diff  command. Here are 2 equivalent variants with  D  command:

f:=x->(x^12-x*sin(x^11))/x^34+exp(sqrt(x+4))*ln(abs(cos(x)^5-6));

D[1,1](f)(5);

(D@@2)(f)(5);

Compare with the result of work  diff  and  eval  command:

eval(diff(f(x),x$2), x=5);

Example:

Ugen := .1049253920*phi(x)^2+.2490325160*psi(x)^2+0.7218836157e-1*eta(x)^2+0.4163942054e-1*(diff(psi(x), x))^2+0.3590610475e-1*(diff(phi(x), x))^2+0.1916547983e-2*psi(x)*(diff(phi(x), x, x))+0.5733315777e-2*(diff(psi(x), x))*(diff(phi(x), x))-0.3273703041e-1*phi(x)*psi(x)-0.3273703035e-1*psi(x)*eta(x)+0.4980952543e-2*(diff(phi(x), x, x))^2-0.7191876251e-1*phi(x)*eta(x)+0.1153177175e-2*eta(x)*(diff(phi(x), x, x))+0.1073591707e-1*phi(x)*(diff(phi(x), x, x)):

coeff(subs(diff(phi(x),x)*diff(psi(x),x)=t, Ugen),  t);

                                 .5733315777e-2

Maple can not calculate your integral symbolically, but it can do this  numerically for any  n :

f := 2*sin(Pi-Pi*exp(-t)):

A[0] := (1/3)*(int(f, t = 0..3));

A := n->(2/3)*(int(f*cos((2/3)*Pi*n*t), t = 0..3));

seq(evalf(A(n)), n = 1..10);

Addition:  If you want to expand your function into Fourier series  only with cosines, it means that you continue this function from segment [0,3] as even function. I corrected some errors in the formulas. Plotted the original function and series expansion with the first 8 terms:

restart;

f := 2*sin(Pi-Pi*exp(-abs(t))): # by abs the function expanded as even function

A := n->(2/3)*evalf(Int(f*cos((1/3)*Pi*n*t), t = 0..3));

A(0)/2+add(A(n)*cos((1/3)*Pi*n*t), n=1..7);

plot([f, %], t=-3..3, color=[red,blue], thickness=2);

Probably the questioner meant getting real solutions depending on a real parameter m. But Maple does not solve the problem:

restart;

assume(m::realcons);

RealDomain[solve](x^2-(2*(m+1))*x+m^2-2*m+m^2 = 0, x, parametric = full);

                     [{x = m+1+sqrt(-m^2+4*m+1)}, {x = m+1-sqrt(-m^2+4*m+1)}]

Mathematica solves the problem correctly:

Reduce[x^2 - 2*(m + 1)*x + m^2 - 2*m + m^2 == 0, x,
Reals] 

(m == 2 - Sqrt[5] &&
x == 3 - Sqrt[5] - Sqrt[
1 + 4 (2 - Sqrt[5]) - (2 - Sqrt[5])^2]) || (2 - Sqrt[5] < m <
2 + Sqrt[5] && (x == 1 + m - Sqrt[1 + 4 m - m^2] ||
x == 1 + m + Sqrt[1 + 4 m - m^2])) || (m == 2 + Sqrt[5] &&
x == 3 + Sqrt[5] - Sqrt[1 + 4 (2 + Sqrt[5]) - (2 + Sqrt[5])^2])

Explanation:  &&  is logical  and,  ||  is logical  or

S:="110100100011001001010110001000":

L:=[1]:

for i from 2 to length(S) do

if S[i]=S[i-1] then L:=subsop(nops(L)=L[nops(L)]+1, L) else L:=[op(L), 1] fi:

od:

L; 

                               [2, 1, 1, 2, 1, 3, 2, 2, 1, 2, 1, 1, 1, 1, 2, 3, 1, 3]

Addition  -  another variant for long strings:

restart;

S:=StringTools[Random](10^5, 'binary');

L[1]:=1: k:=1:

for i from 2 to length(S) do

if S[i]=S[i-1] then L[k]:=L[k]+1 else k:=k+1: L[k]:=1 fi:

od:

convert(L,list);

fsolve command does not solve the equation with unknown parameters. Use solve command:

eq3 := k*y^2+x^2 = 4:

eq4 := -h*x^2+a*y = 0:

solve({eq3, eq4}, {x, y});

allvalues(%);

StringTools[Random](30, 'binary');

StringTools[CharacterFrequencies](%);

                     "110100100011001001010110001000"

                                  "0" = 18, "1" = 12

for n while ithprime(n)<29 do

end do:

n;

                       10

You can immediately get a decimal answer, if you use the  fsolve command, or solve command but at least one of the coefficients of your equation has float type: