An $n \times n$ circulant matrix $M$ is diagonalized by the discrete Fourier basis. This gives us a way to solve $Mx = b$ very efficiently, as follows. Let $F$ be an $n \times n$ matrix whose columns are the discrete Fourier basis vectors. Then $MF = F \Lambda $, where $\Lambda$ is a diagonal matrix whose diagonal entries are the eigenvalues of $M$. So $F^{-1} M F = \Lambda$ is a diagonal matrix. Now, notice that
\begin{align}
& Mx = b \\
\iff &F^{-1} M F F^{-1} x = F^{-1} b \\
\iff& \Lambda \hat x = \hat b
\end{align}
where $\hat x = F^{-1} x$ and $\hat b = F^{-1} b$.
It is now effortless to compute $\hat x = \Lambda^{-1} \hat b$, because we are merely inverting a diagonal matrix. We then obtain
$$
x = F \hat x
$$
and so we have solved $Mx = b$.
Normally multiplying a matrix by a vector takes $O(n^2)$ floating point operations. So already, we are able to solve $Mx = b$ using $O(n^2)$ rather than the usual $O(n^3)$ floating point operations. But, there is something even more tremendous that happens here: computing $F^{-1} x$ can in fact be accomplished with just $O(n \log(n))$ floating point operations, thanks to the wonderful Fast Fourier Transform algorithm. We can also compute $F^{-1} b$ and $F \hat x$ using the blazingly fast FFT algorithm. So, we are able to solve $Mx = b$ with just $O(n \log(n))$ floating point operations. This is a radical improvement over the usual $O(n^3)$ cost of solving a linear system of equations.
Finally, in your problem, if $A$ represents a convolution operator on $\mathbb R^n$ or on $\mathbb C^n$ (using periodic boundary conditions) then $A$ is circulant, so $M = A^T A + \lambda I$ is also circulant, so the above technique can be applied.
If $A$ represents a convolution operator on images (that is, on $p \times q$ rectangular arrays of real or complex numbers), then the above argument still works, but now the "discrete Fourier basis" is a basis for $\mathbb C^{p \times q}$.
