What's wrong to fit periodic data with polynomials? Suppose we have toy daily temperate data and we want to fit a model.

A reasonable thing to do is fitting a periodic model with Fourier basis
$$
f(x)=\beta_0+\beta_1 \cos(2\pi x/24)+\beta_2 \sin(2\pi x/24)
$$
So the Fourier basis expansion of data matrix $\mathbf X$ is
$$
\begin{bmatrix} 1&\cos 0 & \sin 0 \\ 1&\cos \frac \pi 4 & \sin \frac \pi 4 \\ \cdots & \cdots & \cdots\\ 1&\cos \frac {7\pi} 4 &\sin \frac {7\pi} 4  \end{bmatrix}
$$
On the other hand, suppose I do not know Fourier expansion and only know polynomial fit. So I fit data with a third order polynomial, where
$$
f(x)=\beta_0+\beta_1 x+\beta_2 x^2 +\beta_3 x^3 
$$
The polynomial basis expansion of data matrix $\mathbf X$ is (for demo purpose, I am not using orthogonal polynomial which will be ill conditioned in real world problems.)
$$
\begin{bmatrix} 1&\ 0 & 0 & 0 \\ 1& 3 & 3^2 & 3^3\\  & \cdots\\ 1&\ 21 &21^2  & 21^3\end{bmatrix}
$$
The two fits are shown below, and they are "similar". My question is, what's wrong with the polynomial fit on periodical data? In this case we do not have the extrapolation and the time should always be $[0,23]$
    d = data.frame(t=c(0,3,6,9,12,15,18,21),
    temp=c(-2.2,-2.8,-6.1,-3.9,0,1.1,-0.6,-1.1))
    X=cbind(1,cos(2*pi*d$t/24),sin(2*pi*d$t/24))
    coeff = solve(t(X) %*% X, t(X) %*% d$temp)
    X2=cbind(1,d$t, d$t^2, d$t^3)
    coeff_2 = solve(t(X2) %*% X2, t(X2) %*% d$temp)
    plot(d$t,d$temp,type='b')

    d_new = seq(0,24,0.1)
    X=cbind(1,cos(2*pi*d_new/24),sin(2*pi*d_new/24))
    X2=cbind(1,d_new, d_new^2, d_new^3)

    lines(d_new,X %*% coeff, type='l',col='red')
    lines(d_new,X2 %*% coeff_2, type='l', col='blue')


 A: In just the dataset you've provided, the only real downside to using polynomials over the Fourier basis is the issue of discontinuity at $T = 0$ and $T = 24$. As you stated, you can add constraints to fix this up if you really wished to. 
But more typically for this type of data, we observe several cycles. In this case, it would be the number of days of data. The whole point is to take advantage of the fact that 3pm on Monday has very similar features to 3pm on Tuesday. This relation would not show up at all in the "vanilla" polynomial expansion, and so you would not be borrowing at all from different cycles for estimation. For similar reasons, you would have almost no hope of getting a good extrapolation, even just 1 day out, where as even from a very basic Fourier expansion, you could say "I think at 3pm tomorrow, it will probably be the same heat as it usually is at 3pm". 
A: The wrong is that to exactly capture the simplest periodic process such as a monochrome sine wave you need infinite number of polynomial terms. Look at Taylor expansion formula.
Intuitively you want to fit function that (in some sense) looks like  your underlying process. This way you'll have the fewest number of parameters to estimate. Say you have a round hole, and need to fit a cork into it. If your cork is square it's harder to fit it well than if the cork were round.
A: Discontinuity at $T=0$ and $T=24$ is problem. In fact, the plot is misleading because it only plots $T$ up to $21$. If we change the plot code as:
plot(d$t,d$temp,type='b',xlim=c(0,24),ylim=c(-7.5,1.5))
We can see 3rd order polynomial is not a good fit: 
At time $0$, the temperature is $-1.7$, but next day at time $0$ the temperature at $-7.04$ !:

In addition, it is very nature to have function input $T$ as any real number, instead of limited to 0 to 23. 
For example, when $T=25$ it just means 1:00 in next day and $T=-1$ means 23:00 in previous day. Using polynomial basis we need to make to inside 0 to 23 to generate output. 
But with Fourier basis expansion, everything is build in.
A: If You fit data from a limited timeinterval, say one day, using splines, this does not take into account, the values of the preceding and following intervals. You find this effect even with fitting non periodic data with a polynom: the fitted data are "over reacting " to the last and first span interval. One  way to smoothen this, is to repeat the data of the interval to be fitted three times, make the polynomial fit over the long interval and use as a "better" fit only the fitted data of the middle interval.
But certainly, the use of a periodic as a "basis function" is the best approach if You know, that the effect considered is periodic.  Polynoms with a limited number of coefficients can not fit a periodic signal.As already said by Aksakal
