Sinusoidal Generation in Recursive Algorithm

Question

I need to generate sinusoidal values for varying frequencies. I'm making a DTMF tone generator but I must generate my own values of sine using recursive algorithms. The exact wording of how I'm allowed to solve this is "Generate sin while using only real-valued multiplication, addition, and array operations."

If someone can steer me in the right direction you'll help me out a lot and I'd really appreciate it.

Here is what I've come up so far.

frequency = 200;

subplot(2,1,1);

[x, t] = generate_tones(frequency, 0, 8000, 0.065);

plot(t,x);
title('Actual');

subplot(2,1,2);
[x,t] = generate_tones_test(frequency, 0, 8000, 0.065);
plot(t, x);
title('MY Sine function');

% TRIED AND TRUE
function [x,t]=generate_tones(lo,hi,fs,dur)
    t=0:1/fs:dur;
    x=oscil(lo,t)+oscil(hi,t);
    x=[x zeros(size(x))];
    t=[t t+dur];
end

function x=oscil(f,t)
    x=zeros(size(t));
    fprintf('Actual-Values');
    x=sin(2*pi*f*t);
end

% MY CODE
function [x,t]=generate_tones_test(lo,hi,fs,dur)
    t=0:1/fs:dur;
    fprintf('My Values');
    x=sine(lo,t)+sine(hi,t)
    x=[x zeros(size(x))];
    t=[t t+dur];
end

function [output_frequency] = sine(frequency_to_generate, t)
    [m, N] = size(t);

    v = zeros(1,N);
    x = zeros(1,N);

    x(2) = 1;
    x(1) = 0;

    for n = 2:N-1
        T=t(n) * frequency_to_generate;

        x(n+1) = (x(n) + T*v(n))/sqrt(1+(T*T));
        v(n+1) = (v(n) - T*x(n))/(sqrt(1+(T*T)));
    end

    output_frequency = x;
end

Answer 1

To show another way.

The analog generating ODE for sinusoidal wave is $$ \omega ^{\,2} y''(t) + y(t) = 0 $$

Let's see what is the Finite Differences generating equation.

Consider the discrete sinusoidal signal $$ \sin \left( {{{2\pi } \over T}n\,\tau } \right) = \sin \left( {{{2\pi \tau } \over T}n} \right) = \sin \left( {\mu \,n} \right) $$ where the meaning of the symbols is clear.

Then the finite forward difference wrt $n$ is $$ \eqalign{ & \Delta _{\,n} \sin \left( {\mu \,n} \right) = \sin \left( {\mu \left( {n + 1} \right)} \right) - \sin \left( {\mu n} \right) = \cr & = \sin \left( {\mu \left( {n + 1/2 + 1/2} \right)} \right) - \sin \left( {\mu \left( {n + 1/2 - 1/2} \right)} \right) = \cr & = 2\sin \left( {\mu /2} \right)\cos \left( {\mu \left( {n + 1/2} \right)} \right) \cr} $$ and analogously $$ \eqalign{ & \Delta _{\,n} \cos \left( {\mu \,n} \right) = \cos \left( {\mu \left( {n + 1/2 + 1/2} \right)} \right) - \cos \left( {\mu \left( {n + 1/2 - 1/2} \right)} \right) = \cr & = - 2\sin \left( {\mu /2} \right)\sin \left( {\mu \left( {n + 1/2} \right)} \right) \cr} $$

Therefore $$ \eqalign{ & \Delta _{\,n} ^2 \sin \left( {\mu \,n} \right) = - 4\sin ^2 \left( {\mu /2} \right)\sin \left( {\mu \left( {n + 1} \right)} \right) \cr & \quad \quad \Downarrow \cr & \sin \left( {\mu \,\left( {n + 2} \right)} \right) - 2\left( {1 - 2\sin ^2 \left( {\mu /2} \right)} \right)\sin \left( {\mu \,\left( {n + 1} \right)} \right) + \sin \left( {\mu \,\left( n \right)} \right) = 0 \cr & \quad \quad \Downarrow \cr & \sin \left( {\mu \,\left( {n + 2} \right)} \right) = 2\cos \mu \;\sin \left( {\mu \,\left( {n + 1} \right)} \right) - \sin \left( {\mu \,\left( n \right)} \right) \cr} $$ and there are the basis for your generator.

Note that you only need to compute $$ \cos \mu = \cos \left( {2\pi {\tau \over T}} \right) $$ and the seed values $$ \eqalign{ & \sin \left( {\mu \, \cdot 0} \right) = 0 \cr & \sin \left( {\mu \cdot 1} \right) = \sin \left( {2\pi {\tau \over T}} \right) \cr} $$ which you have better to restate at the end of a cycle ($n \tau /T$ integer), to avoid build up of error.
Also, you may introduce them (with appropriate sign) at half of the cycle ...

Answer 2

The sequence $\,a_n\,$ given by the recurrence relation $$a_{n+1} = 2\,\cos(t)\,a_n - a_{n-1} \tag{1}$$ is a sine wave. That is, if $\,a_0 = 0\,$ and $\,a_1 = \sin(t),\,$ then $\,a_n = \sin(n\,t)\,$ for all integer $n.$ In general, the amplitude and phase of the wave will depend on the initial values $\,a_0,\, a_1.$

Your algorithm is similar to the famous Minsky circle algorithm $$ a_{n+1} = a_n + c\, b_n, \,\, b_{n+1} = b_n - c\, a_{n+1} \tag{2}$$ where $\,c = 2\,\sin(t/2)\,$ and $\,a_n=\sin(n\,t),\, b_n=\cos(n\,t+t/2).$

Your code is not wrong, but more complicated than needed. The key idea is to use the addition formula (or theorem) of sine and cosine. For example, you can use two sequences (as you did) instead of just $\,a_n.\,$ That is, $$ a_{n+1} = \sin(t)\,b_n + \cos(t)\,a_n, \,\, b_{n+1} = \cos(t)\,b_n - \sin(t)\,a_n \tag{3}$$ has solution $\,a_n = \sin(n\,t)\,$ as before and also $\,b_n = \cos(n\,t).\,$ The problem here is that roundoff errors will change the amplitude over time which must be adjusted for (as you did). That does not happen using equations $(1)$ or $(2)$.