Finding polynomial approximations of a sine wave
about a decade ago i did this for an unnamed music synthesizer company who had R&amp;D not too far from my condo in Waltham MA.  (can't imagine who they are.)  i don't have the coefficients.

but try this:

$$\begin{align}
f(x) &amp;\approx \sin\left(\tfrac{\pi}{2} x \right) \qquad \qquad \text{for }-1 \le x \le +1 \\
\\
     &amp;= \tfrac{\pi}{2} x (a_0 + a_1 x^2 + a_2 x^4) \\
\end{align}$$

this guarantees that $f(-x) = -f(x)$. 

To guarantee that $f'(x) \Big|_{x=\pm1} =0$ then

$$ f'(x) = \tfrac{\pi}{2}( a_0 + 3 a_1 x^2 + 5 a_2 x^4) $$

$$ a_0 + 3 a_1 + 5 a_2 = 0 \tag{1}$$

That's the first constraint.  To guarantee that $|f(\pm 1)| = 1$, then

$$ a_0 + a_1 + a_2 = \tfrac{2}{\pi} \tag{2}$$

That's the second constraint.  Eliminating $a_0$ and solving Eqs. (1) and (2) for $a_2$ in terms of $a_1$ (which is left to adjust):

$$ a_0 =  \tfrac{5}{2 \pi} - \tfrac{1}{2} a_1 $$

$$ a_2 = -\tfrac{1}{2 \pi} - \tfrac{1}{2} a_1  $$

Now you have only one coefficient, $a_1$, left to twiddle for best performance:

$$ f(x) = \tfrac{\pi}{2} x \left((\tfrac{5}{2 \pi} - \tfrac{1}{2} a_1) + a_1 x^2 -  (\tfrac{1}{2 \pi} + \tfrac{1}{2} a_1)x^4 \right) $$

This is the way I would twiddle $a_1$ for best performance for a sine wave oscillator.  I would adjust use the above and the symmetry of the sine wave about $x=1$ and place exactly one entire cycle in a buffer with a power of two number of points (say 128, i don't care) and run the FFT on that perfect cycle.

The FFT result bin 1 will be the strength of the sine and should be about $N/2$.  Now you can adjust $a_1$ to bring your 3rd harmonic distortion up and down.  I would start with $a_1 \approx \tfrac{5}{\pi} - 2$ so that $a_0 \approx 1$.  That's in bin 3 of the FFT results  But the 5th harmonic distortion (value in bin 5) will be consequential (it will go up as the 3rd harmonic goes down).  I would adjust $a_1$ so that the strength of the 5th harmonic level is equal to the 3rd harmonic level.  It will be around -70 dB from the 1st harmonic (as I recall).  That will be the nicest-sounding sine wave from a cheap, 3-coefficient, 5th-order, odd-symmetrical polynomial.

Someone else can write the MATLAB code.  How does that sound to you?