Strange result of FourierCosCoefficient on $\cos^2 t$

@user64494 is correct that the result arises because some of Mathematica's solvers return results that are only generically correct (that is, incorrect only for finitely many exceptions or on a manifold of dimension lower than the dimension of the domain of the problem). Clearly (I think it's clear), WRI does not consider all such behavior a bug. However, users may well feel that such behavior is undesirable and may insist on calling it a bug. I can explain where the "genericity problem" arises in this example.

FourierCosCoefficient works by first trying table lookup, which succeeds when the expression is in Fourier-series form (e.g. the OP's second example). When table lookup fails, it tries Integrate. Here the problem comes from integrals like this one, which is only generically true (for n != 0):

Integrate[Cos[n t], t]
(*  Sin[n t]/n  *)

This leads to limits of such functions at t -> Pi, which gives the result 0, which again is incorrect for n == 0. One can get around this specific problem by using the identity $\sin x = x \mathop{\text{sinc}} x$. For example:

Limit[Sin[n t]/n, t -> Pi, Direction -> "FromBelow", Assumptions -> n ∈ Integers]
(*  0  *)

Limit[Sin[n t]/n /. Sin -> (# Sinc[#] &), t -> Pi, Direction -> "FromBelow",
  Assumptions -> n ∈ Integers]
(*  π Sinc[n π]  *)

We can use the Villegas-Gayley trick to catch Limit and apply the identity:

Internal`InheritedBlock[{Limit},
 Unprotect[Limit];
 Limit[f_, args__] /; ! TrueQ[$in] := Block[{$in = True},
   Limit[f /. Sin -> (# Sinc[#] &), args]];
 Protect[Limit];
 FourierCosCoefficient[Cos[t]^2, t, n]
 ]
(*  1/2 (Sinc[(-2 + n) π] + 2 Sinc[n π] + Sinc[(2 + n) π])  *)

FunctionExpand[%, Assumptions -> n ∈ Integers] // InputForm
(*  (KroneckerDelta[2 - n] + 2*KroneckerDelta[n] + KroneckerDelta[2 + n])/2  *)

Table[%, {n, 0, 5}]
(*  {1, 0, 1/2, 0, 0, 0}  *)

I'm pretty sure you don't always want to replace all the sines by sincs, but I think it shows what is going on in the OP's examples.


This is not a bug: Mathematica produces a generic answer by default and most of the coefficients equal zero. The Fourier expansion of Cos[2t] is produced by

FourierCosSeries[Cos[t]^2, t, 5]

$\frac{1}{2} \cos (2 t)+\frac{1}{2} $