Consideramos o problema de interpolação no intervalo [0,2p[
com funções trigonométricas (senos e co-senos, ou
mais simplesmente através da exponencial complexa). Dado um conjunto
de n nós de interpolação t0,
... , tn-1 em [0,2p[ pretende-se
interpolar os pontos (t0, f0), ...
, (tn-1, fn-1).
Através da mudança de variável x = eit,
este problema está relacionado com a interpolação
polinomial
q(t) = a0+ a1x + .... + an-1xn-1 = a0+ a1eit + .... + an-1ei(n-1)t
e atribuindo xk = exp(itk), obtemos um sistema através da matriz de Vandermonde:
é
ê ê ê ê ë |
|
ù
ú ú ú ú û |
= | é
ê ê ê ê ë |
|
ù
ú ú ú ú û |
Interessa-nos especialmente o caso de nós igualmente espaçados
tk=2pk/n,
com k=0, ..., n-1.
Nessa situação obtemos como elementos da matriz xkm
= exp(imtk) = exp(2pikm/n)
para k,m=0, ..., n-1. Em particular, a matriz de Vandermonde é
simétrica, com a primeira linha e coluna igual a 1 (caso k=0
ou m=0). Obtemos assim um sistema
é
ê ê ê ê ë |
|
ù
ú ú ú ú û |
é
ê ê ê ê ë |
|
ù
ú ú ú ú û |
= | é
ê ê ê ê ë |
|
ù
ú ú ú ú û. |
que abrevidamente escrevemos na forma Wa = f.
[W*W]km = | n-1
S r =0 |
e -2p ikr/n e2p irm/n = | n-1
S r =0 |
e2p ir(m-k)/n |
[W*W]km = | n-1
S r =0 |
e2p ir(m-k)/n = | n-1
S r =0 |
Zr = | Zn-1
Z -1 |
= 0 |
Definição: Dado um vector v = (v0,
..., vn-1) designamos por Fv =W*v,
a Transformada de Fourier Discreta de v, ou seja, é o
vector cujas componentes são dadas por
[Fv]m = | n-1
S r = 0 |
vr e -2pi rm/n, (m = 0, ..., n-1) |
Observação: A Inversa da Transformação
de Fourier Discreta pode ser obtida de forma semelhante. Basta reparar
que WFv = WW*v = nv.
Portanto considerando a inversa é dada por F-1u
= (1/n) Wu, o que dá uma expressão
semelhante à anterior, excepto que o sinal na exponencial é
positivo e há uma divisão por n. Alguns autores podem
considerar a raiz de 1/n em ambas as expressões ou colocar
a divisão por n na definição de Fv.
A noção de Transformação de Fourier Discreta
(DFT) tem várias aplicações, e o caso da interpolação
é apenas um aspecto.
O cálculo da DFT envolve O(n2) operações,
mas através de uma sucessiva decomposição é
possível definir um algoritmo que reduz esse número de operações
para O(n log(n)). Trata-se de um algoritmo de Cooley e Tukey
(de 1965, apesar de ter sido usado 160 anos antes, num trabalho de Gauss
sobre trajectórias de asteróides). Esse algoritmo, ou variantes,
são denominados Transformação de Fourier Rápida
(FFT).
Valores Reais - Interpolação com senos e co-senos:
Iremos simplificar os resultados anteriores, para o caso em que fk
são reais e q(t) também.
No caso em que os valores de vr são reais
verifica-se que [Fv]n-m são os conjugados
de [Fv]m porque
[Fv]n-m = | n-1
S r = 0 |
vr exp(-2pi r (n-m)/n ) = | n-1
S r = 0 |
vr exp(+2pi r m/n ) = [Fv]m* |
exp(i (n-m) tk) = exp(i(n-m) 2pk/n) = exp(-im 2pk/n) = exp(-i mtk).
Número ímpar de pontos:
No caso em que temos n=2p+1 pontos, o valor nos nós tk
, é dado por
fk = q(tk) = a0 + a1exp(itk) + .... + apexp(iptk)+ ap+1exp(i(p+1)tk) + .... + a2pexp(2iptk)
Em termos do vector solução, vimos que a
= Fv
é da forma a = ( a0
,
a1,
..., ap , ap*, ..., a1*
).
Também vimos que exp(i(2p+1-m)tk)
= exp(-imtk), logo
exp(i(p+1)tk) = exp(-iptk), ..., exp(2iptk) = exp(-itk) .
A função Q(t) = a0+ 2 Re[ a1eit + .... + apeipt] é real e, como vimos, é interpoladora nos 2p+1 pontos, trata-se do "polinómio" trigonométrico de grau p que se escreve através de senos e co-senos da seguinte forma:
Q(t) = a0+ 2 [ Re(a1)cos(t) + .... + Re(ap)cos(pt) ] - 2 [ Im(a1)sin(t)+ .... + Im(ap)sin(pt) ]
porque Re[ak exp(ikt) ] = Re(ak)cos(kt) - Im(ak)sin(kt) .
Assim, no caso real, é possível obter uma função
real interpoladora em 2p+1 pontos usando a parte real obtida pela interpolação
com p+1 funções base complexas 1, eit,
..., eipt, ou de forma equivalente, usando 2p+1
funções base reais, 1, cos(t), ..., cos(pt), ..., sin(t),
..., sin(pt).
Em termos de cálculo, basta avaliar a primeiras p+1 componentes
dada pela Transformação de Fourier Discreta (a0
,
a1,
..., ap).