This is done with the identity option to solve.

solve(
   identity(
      cos(t)^6+a*cos(t)^4*sin(t)^2+b*cos(t)^2*sin(t)^4+c*sin(t)^6 = cos(6*t), 
      t
   ), 
   {a,b,c}
);
                   {a = -15, b = 15, c = -1}

These very small imaginary parts are due to round-off error in float computations. I call them spurious imaginary parts, because the true answers are real. They can be removed by

simplify(fnormal([%]), zero);

In some exceptional cases, you may need to specify a tolerance (a number of digits) as a second argument to fnormal. See help page ?fnormal.

The RealDomain package is not appropriate for this. It's better to generate the answers with the imaginary parts and then remove them than it is to try to generate them without the imaginary parts.

The command to turn an expression into a function is unapply.

Mat:= module()
option package;
uses CF= CurveFitting;
export 
     LinModel:= proc(xliste::list, yliste::list, $)
     local a, b, x;
          unapply(
               CF:-LeastSquares(xliste, yliste, x, 'curve'= a*x+b),
               x
          )
     end proc; # LinModel 
end module; #Mat

If you are returning a function, there's no need to pass x to LinModel. Thus, I've made it a local.

What you've called sys is already a set because you put it in braces; then you've put that and the initial condition into more braces. So, your first argument to dsolve isn't a set of equations; rather, it's a set that contains an equation and also contains a set of equations. So, you can make your dsolve command

dsolve({sys[], varphi(0) = 0});

The [] after sys extracts the elements of the set sys. I do not understand why you include the arguments Xa, Xb, Ya, Yb in the command.  

Substitute exact fractions such as beta = 3/10 rather than decimals such as 0.3. If you use a decimal before calling dsolve, it will automatically change it to an exact fraction, which is why your (4) is correct. 

L:= [a,b]:
`/`(L[]);
             a/b

If for some strange reason you want your convert command to work, do

`convert//`:= (L::list)-> `/`(L[]):

Now,
convert([a,b], `/`);
             a/b
works. I say "strange" reason because the first way shown seems so much easier than using convert. But in general, to make convert(..., foo) work, define a procedure named `convert/foo`; or if you're curious how a stock command such as convert(..., list) works, then look at showstat(`convert/list`).

It'll be a lot easier to display summation notation than "..." form, so how about

Sum(q[i], i= 1..n);

?

You haven't supplied any inital conditions. So, is this what you want?

DEtools:-dfieldplot(
   [diff(x1(t),t) = 3*x1(t) + 2*x2(t),
    diff(x2(t),t) = 5*x1(t) + 1*x2(t)
   ],
   [x1,x2](t),
   t= 0..2, x1= -Pi..Pi, x2= -Pi..Pi
);

Increasing the implicitplot option gridrefine from its default 1 is more sophisticated than using a large value of numpoints.

plots:-implicitplot(
     eval(
          [[x^3+y^3-3*a*x*y, x+y+a], x= -3*a..3*a, y= -3*a..3*a],
          a= 1
     )[],
     gridrefine= 2, linestyle= [solid, dash], color= [red, blue]
);

Substitute z= exp(I*theta). Thus dz = I*exp(I*theta)*dtheta. 
z = 1 => theta(2*Pi); z = -1 => theta = Pi, and these limits trace the bottom half of the circle.

int(I*exp(I*t)/sqrt(exp(I*t)), t= 2*Pi..Pi);

The Euler totient is emphatically NOT the number of prime factors; indeed, it's nearly the opposite of that. The totient of n > 0 is the count of k, 0 < k <= n, such that k and n share NO prime factors. The symbol phi is usually used for this function. It's obvious from this definition that for prime p, phi(p) = p-1.

If n = p*q with p and q distinct primes, then phi(n) = (p-1)*(q-1). This fact is a very important aspect of RSA cryptograhy. Usually the only ways to compute the totient of a number all require the number's prime factorization. Since that is currently generally computationally infeasible when p and q are both very large, computing phi(n) is likewise infeasible. However, it hasn't been proven that there's no feasible algorithm. Quantum computation may be able to crack it if and when quantum computers get large enough.

Why don't you want to use Heaviside? A similar alternative is 

f1:= x-> piecewise(x >= 0, 1, 0);

If you want a more-primitive procedure based on if statements, then you need to delay the evaluation of x >= 0 in those cases where it can't be decided (such as when x has no value). Christian showed one way of doing that; here's another:

f1:= proc(x)
local r:= PiecewiseTools:-Is(x >= 0);
   `if`(r::truefalse, `if`(r, 1, 0), 'procname'(x))
end proc:

There's nothing wrong with your Maple: It can do the integral fairly quickly if you tell it, essentially, that the expression under the square root is nonnegative over the interval of integration. The conditions d/L__1 >= 2, d > 0 are sufficient. Use an assuming clause, like this:

simplify(...) assuming d/L__1 >= 2, d > 0;

where ... is exactly what you already have in the parentheses.

If these assumptions are not actually true for your problem, some deeper thought will be required.

It's nearly impossible to accurately compute the eigenvalues of a large matrix by finding the roots of the characteristic polynomial. Instead, use command LinearAlgebra:-Eigenvalues.

I show how to do that in this Post from a year ago: "Numerically solving BVPs that have many parameters" The BVP in that Post is extremely similar to yours; indeed, yours may be a direct simplififation of it. However, I wouldn't recommend that you attempt this if this is your first experience with Maple.