1.3.3 Raízes de Equações — Raízes de Polinômios — O Método de Bairstow
1. Introdução
O Método de Bairstow consiste em dividir um polinômio em fatores isolados que permitam aproximar suas raízes.
Possui a grande vantagem de conseguir extrair todas as raízes de um polinômio arbitrário sem a necessidade de percorrer todo o espaço de busca,como fariam o métodos de Muller ou Newton Raphson, que sempre convergem à aproximação mais próxima.
Além disso, ainda possui o benefício de fornecer as raízes complexas à uma taxa de convergência quadrática utilizando-se apenas de aritmética real.
2. Desenvolvimento do Método
Seja a função um polinômio geral:
(eq1)
$$ f_{n} = a_{0} + a_{1}x + a_{2}x^2 + \cdots + a_{n}x^n $$
Se dividirmos o polinômio 1 pelo fator $$x-t $$ , obtemos outro polinômio de um grau a menos:
(eq2)
$$ f_{n-1} = b_{1} + b_{2}x + b_{3}x^2 + \cdots + b_{n}x^{n-1} $$
Caso esta divisão não gere resto, então $$t $$ é uma raiz do polinômio $$f $$. Mas, tomando o pior caso, onde há um resto $$R $$ . Nesse caso $$R = b_0 $$ e os coeficientes poderão ser calculados pela relação:
(eq3)
$$ b_{n} = a_{n} $$
e
(eq4)
$$ b_{i} = a_{i} + b_{i+1}t $$
onde $$i=n-1, n-2, \ldots, 0 $$ .
O Método de Bairstow segue um procedimento semelhante ao descrito acima, contudo, ao invés de dividir $$f $$ por um fator $$x-t $$ , ele divide o polinômio pelo fator $$x*x – r*x – s $$ . A introdução deste fator quadrático serve para permitir a determinação de raízes complexas, enquanto que o fator linear permitiria apenas aproximações reais.
Ao dividirmos o polinômio pelo fator de segunda ordem, obtemos:
(eq5)
$$f_{n-2} = b_{2} + b_{3} x + b_{4} x^2 + \cdots + b_{n-1} x^{n-3} + b_{n} x^{n-2}$$
isso gera uma expressão para o resto que é:
(eq6)
$$ R = b_{1} \left(x -r \right) + b_{0} $$
Aplicando o algoritmo de Briot-Ruffine com $$x*x – r*x – s $$ como divisor, obtemos a seguinte relação de recorrência:
(eq7)
$$ b_{n} = a_{n} $$
$$\Downarrow $$
(eq8)
$$ b_{n-1} = a_{n-1} + r b_{n} $$
$$\Downarrow $$
(eq9)
$$ b_j = a_{j} + r b_{j+1} + s b_{j+2} $$
onde $$j = n-2, n-1, \ldots, 0 $$ .
Com isso, o Método de Bairstow se resume a determinar os valores de $$r $$ e $$s $$ que tornam o fator quadrático $$(x*x – r*x – s) $$ exato. Ou seja, ele utiliza uma estratégia para encontrar os valores de $$r $$ e $$s $$ para qual $$R $$ é zero. Quando $$R=0 $$ , $$r $$ será uma raiz do polinômio.
Da equação do resto (Equação 6) é possível notar que quando $$R $$ é nulo, $$b_{0} $$ e $$b_{1} $$ são nulos. A técnica usada por Bairstow segue o desenvolvimento do Método de Newton, expandindo $$b_{0} $$ e $$b_{1} $$ em expressões que os aproximem de zero. Como $$b_{0} $$ e $$b_{1} $$ são funções de $$r $$ e $$s $$ , é possível expandir $$R $$ por Série de Taylor para duas variáveis. Bairstow utiliza uma expansão de até a primeira ordem, gerando:
$$ b_{1}(r+\Delta r, s+\Delta s) = b_{1} + \left( \dfrac{\partial}{\partial r} b_{1} \right) \Delta r + \left( \dfrac{\partial}{\partial s} b_{1} \right) \Delta s $$
$$ b_{0}(r+\Delta r, s+\Delta s) = b_{0} + \left( \dfrac{\partial}{\partial r} b_{0} \right) \Delta r + \left( \dfrac{\partial}{\partial s} b_{0} \right) \Delta s $$
Como $$b_{0} $$ e $$b_{1} $$ serão nulos quando $$R=0 $$ , então é possível reescrever as equações acimas como:
(eq12)
$$ -b_{1} = \left( \dfrac{\partial}{\partial r} b_{1} \right) \Delta r + \left( \dfrac{\partial}{\partial s} b_{1} \right) \Delta s $$
$$ -b_{0} = \left( \dfrac{\partial}{\partial r} b_{0} \right) \Delta r + \left( \dfrac{\partial}{\partial s} b_{0} \right) \Delta s $$
A partir das duas últimas equações, podemos obter os valores de $$b_{0} $$ e $$b_{1} $$ . O problema agora está nas derivadas parciais de $$b $$ em relação a $$r $$ e $$s $$ . Se elas forem determinadas o problema é reduzido à um sistema de duas equações lineares com apenas duas incógnitas: $$\Delta r$$ e $$\Delta s $$.
Tomando as derivadas como constantes, renomeamos elas em termos de um vetor $$c $$ , tal que:
(eq14)
$$ \dfrac{\partial}{\partial r} b_{0} = c_{1} $$
$$ \dfrac{\partial}{\partial s} b_{0} = c_{2} $$
$$ \dfrac{\partial}{\partial r} b_{1} = c_{2} $$
$$ \dfrac{\partial}{\partial s} b_{1} = c_{3} $$
A partir das suposições acimas, o sistema de equações pode ser reescrito com suas derivadas parciais substituídas pela variável $$c $$ , tal que
(eq18)
$$ c_{2} \Delta r + c_{3} \Delta s = -b_{1} $$
O Método de Bairstow utiliza este sistema de equações lineares para ajustar variáveis complexas da mesma forma que as reais. As derivadas parciais, obtidas pela divisão sintética de $$b $$, são:
(eq20)
portanto, generalizando para todos termos,
(eq22)
onde $$i = n-1, n-2, \ldots, 1 $$ .
Com este conjunto de equações é possível determinar $$\Delta r $$ e $$\Delta s $$, necessários para as aproximações de $$r $$ e $$s $$ , sendo o processo de cálculo do erro relativo de cada iteração é obtido pelas expressões:
(eq23)
e
(eq24)
Quando ambos erros estiverem abaixo do máximo tolerado, $$r $$ e $$s $$ estarão aproximados de forma que as raízes $$x $$ serão determinadas por
- Pela expressão abaixo, caso o polinômio que deseja-se encontrar a raiz seja quadrático:
(eq25)
$$ x = \dfrac{1}{2}r + \dfrac{1}{2} \sqrt{r^2 + 4s} $$
- Pela expressão abaixo, caso o polinômio seja de grau um:
(eq26)
$$ x = – \dfrac{s}{r} $$
- Caso o polinômio seja de grau maior ou igual a três, o Método de Bairstow será aplicado ao quociente para o cálculo de novos valores $$r $$ e $$s $$ usando a função abaixo. Os valores de $$r $$ e $$s $$ calculados anteriormente servem como parâmetro inicial para a próxima iteração. </p> </li> </ol>
</p>
3. Implementação
Implementação em Python:
def bairstow(a, rr=1, ss=1, errto=0.001, imax=500):
"""
Calculate the roots of a function (using the Bairstow Method)
bairstow(a, polDegree, errto, rr, ss, imax)
* a[]: list with the coeficients
* rr: initial aprox. of "r" param
* ss: initial aprox. of "s" param
* errto: tolerated error
* imax: max iterCountation allowed (not inf-loop)
return: a list with all roots of polynom
Author: Pedro Garcia [sawp@sawp.com.br]
see: http://www.sawp.com.br
License: Creative Commons
http://creativecommons.org/licenses/by-nc-nd/2.5/br/
Jan 2010
"""
def filt(l):
r = l.real
i = l.imag
if abs(r) < errto:
r = 0.0
if abs(i) < errto:
i = 0.0
return complex(r, i)
r = rr
s = ss
n = len(a) - 1
ea1 = errto + 1
ea2 = errto + 1
iterCount = 0
roots = [ complex(0,0) for i in xrange(len(a)) ]
b = [ complex(0,0) for i in xrange(len(a)) ]
c = [ complex(0,0) for i in xrange(len(a)) ]
while n >= 3 or iterCount < imax:
iterCount = 0
while ea1 >= errto and ea2 >= errto or iterCount < imax:
iterCount += 1
b[n] = a[n]
b[n-1] = a[n-1] + r * b[n]
c[n] = b[n]
c[n-1] = b[n-1] + r * c[n]
for i in xrange(n-2, -1, -1):
b[i] = a[i] + r * b[i+1] + s * b[i+2]
c[i] = b[i] + r * c[i+1] + s * c[i+2]
det = c[2]*c[2] - c[3]*c[1]
if det != 0.0:
dr = (-b[1]*c[2] + b[0]*c[3])/det
ds = (-b[0]*c[2] + b[1]*c[1])/det
r = r + dr
s = s + ds
if r != 0.0:
ea1 = abs(dr/r)
if s != 0.0:
ea2 = abs(ds/s)
else:
r += errto
s += errto
iterCount = 0
(r1, i1, r2, i2) = quadratic_solve(r, s)
roots[n] = complex(r1, i1)
roots[n-1] = complex(r2, i2)
n = n - 2
for i in xrange(n):
a[i] = b[i+2]
if n == 2:
r = -a[1]/a[2]
s = -a[0]/a[2]
(r1, i1, r2, i2) = quadratic_solve(r, s)
roots[n] = complex(r1, i1)
roots[n-1] = complex(r2, i2)
else:
roots[n] = complex(-a[0]/a[1], 0.0)
print iterCount
roots[1:len(roots)].sort(key=lambda r: r.real)
return map(filt, roots[1:len(roots)])
Note que a função quadratic\_solve foi incluída para deixar explícita a resolução da fórmula quadrática. A implementação desta função é:
def quadratic_solve(r, s):
"""
Calculate the roots of a parabolic equation: a*x^2 + b*x + c = 0 (*)
quadratic_solve(r, s)
* r: b coeficient of (*)
* s: a*c, coeficients of (*)
return: a tuple with roots of (*)
Author: Pedro Garcia [sawp@sawp.com.br]
see: http:!www.sawp.com.br
License: Creative Commons
http://creativecommons.org/licenses/by-nc-nd/2.5/br/
Jan 2009
"""
delta = r*r + 4.0*s
if delta > 0:
r1 = (r + sqrt(delta))/2.0
r2 = (r - sqrt(delta))/2.0
i1 = 0.0
i2 = 0.0
else:
r1 = r/2.0
r2 = r1
i1 = sqrt(abs(delta))/2.0
i2 = -i1
return (r1, i1, r2, i2)
Observe que o programa altera os valores de $$r $$ e $$s $$ com algum incremento aleatório e repete o procedimento da resolução das equações quando a determinante é nula, forçando o programa a convergir à uma aproximação melhor.
Um exemplo da utilização desta sub-rotina pode ser encontrado em: http://www.sawp.com.br/code/rootfind/bairstow.py
A utilização desta função não é tão simples, pois exige como parâmetro de entrada um vetor com os coeficientes do polinômio. Embora a utilização de listas para representar um polinômio seja algo comum para muitos programadores, certamente é uma abordagem menos intuitiva que passar uma função como referência.
Também há disponível uma outra implementação deste algoritmo em Fortran 90 no endereço http://www.sawp.com.br/code/rootfind/bairstow.f90
</p> </p>4. Copyright
Este documento é disponível sob a licença Creative Commons. As regras dos direitos de cópia deste conteúdo estão acessíveis em http://creativecommons.org/licenses/by-nc-nd/2.5/br/.