MAC0338  Analise de Algoritmos

Máximo divisor comum

"Analysis of algorithms is the name the author likes
to use to describe investigations such as this. The general idea is
to take a particular algorithm and to determine its quantitative behaviour; occasionally
we also study wheher or not an algorithm is optimal in some sense.

Donald E. Knuth, The Art of Computer Programming: Fundamental Algorithms.

  Está página é baseada na seção 2 do capítulo 31 do CLRS e capítulo 1 do TAOCP volume 1.

 

 

  Suponha que m e n são números inteiros. Dizemos que n divide m se m % n == 0. Note que,

m % n == 0   se e somente se   m == dn para algum inteiro d.
As vezes escrevemos n\m como abreviação de "n divide m.

O máximo divisor comum de dois números inteiros m e n, denotado por mdc(m,n), é o maior número inteiro que divide ambos. Assim,

mdc(m,n)  =  max{ d  |  d divide m e k divide n }.
É claro que se m == n == 0, então o máximo não faz sentido. Logo, estaremos sempre supondo que pelo menos um dos números dados é diferente de zero. Nossa atual preocupação é com o seguinte problema
PROBLEMA:  Dados número inteiros m e n determinar mdc(m,n).

  Por simplicidade, na discussão a seguir estamos supondo que m >= 0 e n >= 0. Isto não é nenhum pecado, já que mdc(m,n) == mdc(abs(m),abs(n))

Vamos examinar dois algoritmos para o problema: um óbvio e lento e outro sofisticado e bem mais rápido.

 

Algoritmo café-com-leite

  Vamos começar com um algoritmo simples e óbvio:

 /*
  * Ao receber dois números inteiros não-negativos m e n, não
  * ambos nulos, está função devolve o máximo divisor de m e n.
  */  
  long
  cafe_com_leite(long m, long n)
  {
    long d;

    if (n == 0) return m;
    if (m == 0) return n;
  
    d =  min(m,n);

    while (/*1*/ m % d != 0 || n % d != 0)
      { /*2*/
        d--;
      }

    /*3*/
    return d;
  }

  Aa descrição do algoritmo utilizou a linguagem C, mas nem sempre será assim. Muitas vezes optaremos por um descrição mais "alto nível" que enfatise os pontos realmente cruciais e não nos percamos com detalhes menos relevantes.

 

Correção

  Passamos agora a verificar a correção do algoritmo.
[Correção da função == a função funciona == a função faz o que promete.]

  A correção de algoritmos iterativos é comumente baseadas em demonstrações da validade de invariantes. Estes invariantes são afirmações ou relações envolvendo os objetos mantidos pelo algoritmo. Eis os nossos invariantes para o algoritmo. Em /*1*/ vale que

(i0) 1 <= d <= min(m,n), e
(i1) m % t != 0 ou n % t != 0 para todo número inteiro t > d,
e em /*2*/ vale que
(i2) em /*2*/ vale que m % d != 0 ou n % d != 0.
Verificaremos que essas afirmações, no momento candidatos a invariantes, são de fato válidas.
[A propósito, eu não sei se é a invariante ou o invariante.]

  Antes de mais nada, note que se de fato (i1) vale então o algoritmo devolve o prometido. É evidente que em /*3*/, antes do algoritmo devover d, vale que

m % d == 0 e n % d == 0.
Como (i1) vale em /*1*/, então (1) também vale em /*3*/. Donde, nenhum número inteiro maior que o número d devolvido divide m e n. Portanto, se o algoritmo devolve algum valor, este é de fato o mdc(m,n). O "se" é um certo excesso de zelo, já que, a rigor, ainda não verificamos que o algoritmo sempre pára quando alimentado com números inteiros satisfazendo as especificações.

  Invariantes são assim mesmo. A validade de alguns principais torna a correção do algoritmo (muitas vezes) evidente. Os invariantes secundários servem para confirmar a validade dos principais.
[Hmmmm, talvez os adjetivos principal e secundário sejam inapropriados.]

  Passemos a verificar a validade dos invariantes. A condição do while faz a validade (i2) em /*2*/ óbvia. Este foi fácil!

  Inicialmente d vale min(m,n) que é pelos menos 1, já que m e n são não nulos. Assim, quando a execução do programa alcança pela primeira vez o while tem-se que todo número inteiro maior que d == min(m,n) não divide m ou não divide n; em particular não divide min(m,n). Logo, (i0) e (i1) valem na primeira vez que a execução do algoritmo alcança o ponto /*1*/.

  Considere agora uma iteração do while em que o comando "d--;" é executado. Suponha que nesta iteração (i0) e (i1) valem em /*1*/, antes da execução do comando. Logo, (i0) e (i1) também valem em /*2*/, imediatamente antes da execução do comando d--;. Em /*2*/ também vale (i2). Portando, em /*2*/ valem (i0) e (i1) com d' == d-1 no papel de d. Como a iteração decrementa o valor de d, então (i0) e (i1) valem na primeira vez que a execução do algoritmo alcança o ponto /*1*/ novamente.

  Resumindo, acabamos de demonstrar que

se (i0) e (i1) valem em /*1*/ antes de uma execução do comando "d--", (i0) e (i1) valem na próxima vez que a execução alcaça o /*1*/
Como (i0) e (i1) valem na primeira vez que a execução passa por /*1*/, então (i0) e (i1) valem sempre que a execução do algoritmo está no ponto /*1*/. Portanto, (i0) e (i1) são invariantes genuinos.

  Invariantes, além de serem uma ferramente útil para demonstra a correção de algoritmos iterativo, eles nos ajudam a compreender o funcionamento do algoritmo. De certa forma, eles "espelham" a maneira que entendemos o algoritmo.

 

Desempenho

  Quantas iterações do while essa função faz? Equivalentemente, quantas vezes o comando "d--" é executado? A resposta é

min(m,n)-1   no pior caso.

Esse número de iterações ocorre quando mdc(m,n) = 1, ou seja, m e n são relativamente primos.

Por exemplo, para a chamada mdc(317811,514229) a função executará 317811 - 1 iterações e para a chamada mdc(2147483647,2147483646) a função faz 2147483646 - 1 iterações, pois mdc(317811,514229) == 1 e mdc(2147483647,2147483646) == 1.

Exercício

  1. Como diz um professor meu, só de farra, use o programa mdc.c para determinar mdc(2147483647,2147483646). Quanto tempo o programa gastou para fazer o serviço?

    Olha o tempo que demorou no computador que eu uso:

        meu_prompt> time mdc 2147483647 2147483646 
        mdc de 2147483647 e 2147483646 e' 1.
    
        real    1m8.292s
        user    1m4.062s
        sys     0m1.101s
       

  Neste caso, costuma-se dizer que a quantidade de tempo consumida pelo algoritmo, no pior caso, é proporcional a min(m,n), ou ainda, que o tempo gasto pelo algoritmo é da ordem de min(m,n). A abreviação de "ordem de blá" é O(blá).
[Em breve, tornaremos o significado de O(blá) mais precisa.]

  Portanto, o número de operações feitas pela função cresce (no ior caso) como uma função linear de min(m,n). Isto significa que se o valor de min(m,n) dobra então espera-se que o tempo gasto pela função também dobre.

  É possível resolver o problema de uma maneira mais eficiente? A resposta é SIM.

 

 

Algoritmo de Euclides

  O máximo divisor comum entre dois número pode ser determinado através de um algoritmo de 2300 anos (cerca de 300 A.C.), o algoritmo de Euclides. Para calcular o mdc(m,n) para 0 <= n < m, o algoritmo de Euclides usa a seguinte recorrência:

mdc(m,0)   ==  m;
mdc(m,n)   ==  mdc(n, m % n),  para n > 0.

Assim, por exemplo,

mdc(12,18) == mdc(18,12) == mdc(12,6) == mdc(6,0) == 6.

 

Correção

A correção da recorrência proposta por Euclides é baseada no seguinte exercício.

Exercício

  1. Seja m um inteiro não-negativo e n um inteiro positivo. Para cada inteiro positivo d tem-se que

    (m % d == 0) e (n % d == 0)
      se e somente se  
    (n % d == 0) e ((m % n) % d == 0).

    Em outras palavras, o pares m,n e m,m%n têm o mesmo conjunto de divisores.

  Eis a implementação recursiva em C do algoritmo de Euclides para determinar mdc(m,n).

 /*
  * Ao receber dois números inteiros não-negativos m e n, não
  * ambos nulos, está função devolve o máximo divisor de m e n.
  */  
  long
  euclides_r(long m, long n)
  {
    if (n==0) return m;
    return euclides_r (n, m % n);
  }

  A seqüência abaixo ilustra as chamadas do algoritmo de Euclides para determinar o mdc de 317811 e 514229:

mdc(317811,514229)
  mdc(514229,317811)
    mdc(317811,196418)
      mdc(196418,121393)
        mdc(121393,75025)
          mdc(75025,46368)
            mdc(46368,28657)
              mdc(28657,17711)
                mdc(17711,10946)
                  mdc(10946,6765)
                    mdc(6765,4181)
                      mdc(4181,2584)
                        mdc(2584,1597)
                          mdc(1597,987)
                            mdc(987,610)
                              mdc(610,377)
                                mdc(377,233)
                                  mdc(233,144)
                                    mdc(144,89)
                                      mdc(89,55)
                                        mdc(55,34)
                                          mdc(34,21)
                                            mdc(21,13)
                                              mdc(13,8)
                                                mdc(8,5)
                                                  mdc(5,3)
                                                    mdc(3,2)
                                                      mdc(2,1)
                                                        mdc(1,0)
mdc de 317811 e 514229 e' 1.

Desempenho

  Quanto tempo o algoritmo leva para parar? Este tempo é certamente proporcional ao número de chamadas recursivas. Suponha que a função euclides_r faz k chamadas recursivas e que no início da 1a. chamada ao algoritmo tem-se que 0 < n <= m. Sejam

(m,n) == (m0,n0), (m1,n1),  (m2,n2),  ... ,  (mk,nk) == (mdc(m,n),0),

os valores dos parâmetros no início de cada uma das chamadas da função. Por exemplo, para m==514229 e n==317811 tem-se

(m0,n0) == (514229,317811),
(m1,n1) == (317811,196418),
(m2,n2) == (196418,121393),
        ... ,
(m27,n27) == (1,0).

Estimaremos o valor de k em função de n == min(m,n). Note que mi+1 == ni e ni+1 == mi % ni para i=1,2,...,k. Note ainda que para inteiros a e b, 0 < b <= a vale que

a % b   <   a/2   (verifique!).

Desta forma tem-se que

n2  ==  m1 % n1  ==  n0 % n1  <  n0/2  ==  n/2  ==  n/21
n4  ==  m3 % n3  ==  n2 % n3  <  n2/2  <  n/4   ==  n/22
n6  ==  m5 % n5  ==  n4 % n5  <  n4/2  <  n/8  ==  n/23
n8  ==  m7 % n7  ==  n6 % n7  <  n6/2  <  n/16  ==  n/24
n10  ==  m9 % n9  ==  n8 % n9  <  n8/2  <  n/32  ==  n/25
...
...
...

  Percebe-se que depois de cada 2 chamadas recursivas o valor do segundo parâmetro é reduzido a menos da sua metade. Seja t o número inteiro tal que

2t <= n < 2t+1

Da primeira desigualdade temos que

t <= lg(n),

onde lg(n) denota o logaritmo de n na base 2, ou seja lg(2) == log2(n). Da desigualde estrita, concluimos que

k <= 2t+1

Logo, o número k de chamadas recursivas é não superior a

2 (t+1) - 1 = 2t + 1 <= 2lg(n)+1.

[Quanto n > m tem-se que k <= 2lg(n)+2.]

  Para o exemplo acima, onde m==514229 e n==317811, temos que

2 lg(n) + 1 == 2 lg(317811) + 1 < 2 x 18.3 + 1 < 37.56

e o número de chamadas recursivas feitas por euclides_r(514229,317811) foram 27.

 

  Resumindo, a quantidade de tempo consumida pelo algoritmo de Euclides é, no pior caso, proporcional a lg(n). Este desempenho é significativamente melhor que o desempenho do algoritmo Café-com-leite, já que a função f(n) = lg(n) cresce muito mais lentamente que a função g(n) = n.

n   (int) lg(n)
4 2
5 2
6 2
10 3
64 6
100 6
128 7
1000 9
1024 10
1000000 19
1000000000 29

Exercícios

  1. Demonstre que para inteiros m, k e n mdc(m,n) == mdc(m + kn, n).
  2.  

  3. É verdade que para todo m >= n > 0 a função mdc faz um número de iterações (do comando while) maior ou igual ao feito (pelo comando do ... while) da função euclides_i abaixo.
       long
       euclides_i(long m, long n)
       {
         long r;
    
         do
           {
             r = m % n;
             m = n;
             n = r;
           }
         while (r != 0);
    
         return m;
       }
    
  4.  

  5. Só de farra, compare, na prática, o consumo de tempo das funções mdc.c, euclides_r e euclides_i.

     

  6. Qual função consome uma quantidade maior de memória, euclides_r ou euclides_i? Faça uma estimativa da quantidade de memória consumida por cada uma dessas funções.

 

 

Números de Fibonacci e o algoritmo de Euclides

  Os 10 primeiros números da seqüência de Fibonacci são:

0 1 2 3 4 5 6 7 8 9 10
0 1 1 2 3 5 8 13 21 34 55

  Os números de Fibonacci são definidos (ou podem ser calculados) pela seguinte função recursiva:

 /*
  * A função receber um número inteiros não-negativos  n, e
  * devolve o n-ésimo número de Fibonacci.
  */  
   long
   fibonacci(long n)
   {
     if (n == 0) return 0;
     if (n == 1) return 1;
     return  fibonacci(n-1) + fibonacci(n-2);
   }

Os números de Fibonacci estão intimamente relacionados com o algoritmo de Euclides, como mostra o exercício abaixo.

Exercício

  1. Se m > n >= 0 e se a chamada euclides_r(m,n) faz k >=1 chamadas recursivas, então m >= fibonacci(k+2) e n >= fibonacci(k+1).
    [Sugestao: demonstre por indução em k.]

  Uma conseqüência imediata deste exercício é

Para todo inteiro t >= 1, se m > n >= 0 e n < fibonacci(t+1), então a chamada euclides_r(m,n) faz menos de t chamadas recursivas.

Por exemplo, fibonacci(28) == 317811 e fibonacci(29) = 514229 e a chamada euclides_r(514229,317811) faz 27 chamadas recursivas. Na verdade, números de Fibonacci consecutivos são os que dão mais trabalho para o algoritmo de Euclides (veja exercício 8)..

  Vejamos agora como utilizando os fatos acima e algo que um passarinho me contou podemos obter uma estimativa melhor para o número de chamadas recursivas de euclides_r(m,n).

  Suponha que m > n >= 4 e que um passarinho me contou que para todo número inteiro não-negativo vale que

phit-2 <= fibonacci(t) <= phit-1,

para t >= 2, onde phi == (1+sqrt(5))/2 que é um número entre 1.618 e 1.619. Por exemplo,

phi26 < 1.61926 < 275689 < fibonacci(28)==317811 < 438954 < 1.61827 < phi27

[Hmmm, parece que o passarinho não mentiu.]

  Seja t o número inteiro tal que [(neste ponto precisamos da hipótese n >= 4]

phit+2 <= n < phit+3 < fibonacci(t+1)

então o número k de chamadas recursivas de euclides_r(m,n) é no máximo

t - 1 < (logphi(n)+2) - 1 == logphi(n) + 1.

  Segundo está nova estimativa temos que euclides_r(514229,317811) faz no máximo

logphi(317811) + 1 > log1.619(317811) + 1 > 26.28 + 1 = 27.45.

e o número de chamadas é 27.
[Uauuuu. Isto foi muito perto! Talvez tenha alguma conta errada.]

Exercícios

  1. Demonstre que a chamada euclides_r(fibonacci(k+1),fibonacci(k)) faz k-1 chamadas recursivas.
    [Sugestao: demonstre por indução em k.]
  2.  

  3. Demonstre que números de Fibonacci consecutivos são relativamente primos.
  4.  

  5. Eis a versão iterativa da função fibonacci.
      /*
      * A função receber um número inteiros não-negativos  n, e
      * devolve o n-ésimo número de Fibonacci.
      */      
      long
      fibonacci(long n)
      {
        long i;
        long fib_ant, fib_atual, fib_prox;
    
        if (n == 0) return 0;
    
        fib_ant = 0; fib_atual = 1; 
        for (i = 1; i < n; i++)
          { 
            fib_prox  = fib_ataul + fib_ant;
            fib_ant   = fib_atual;
            fib_atual = fib_prox;
          }
    
        return  fib_atual;
      }
    
    Qual função você usaria? A recursiva ou a iterativa? Por quê? A propósito, note que, no início da i-ésima iteração vale que
    (i0) fib_ant é o (i-1)-ésimo número de fibonacci,
    (i1) fib_atual é o i-ésimo número de finonacci.
    Quando a função pára tem-se que i == n e portanto, de (i1), conclui-se que a função devolve o n-ésimo número de Fibonacci.

     

  6. Matematicamente falando, quanto é 7 % -3? Qual é o valor que o computador devolve? No computador vale que m == (m/n) * n + (m%n)?

     

  7. Demonstre que 1 + phi == phi2, onde phi == (1+sqrt(5))/2.

     

  8. Demonstre a afirmação do passarinho de que

    phit-2 <= fibonacci(t) <= phit-1,

    para t >= 2.

 

 

Algoritmo de Euclides estendido

  O algoritmode Euclides pode ser facilmente modificado para nos fornecer mais que o mdc(m,n). Podemos estender o algoritmo para que ele compute números inteiros a e b tais que

am + bn == mdc(m,n).
Se n == 0, então podemos tomar a == 1 e b == 0. Caso contrário, aplique o método recursivamente com n e m%n nos papéis de m e n, respectivamente, computando assim inteiros e e f tais que
en + f(m%n) == mdc(n,m%n).
Como m%n == m-(m/n)m  e  mdc(n,m%n) == mdc(m,n), esta equação nos diz que
en + f(m-(m/n)m) == mdc(m,n).
O lado esquerdo da equação pode ser reescrito para mostrar a dependência de m e n
fm + (e-f(m/n))n == mdc(m,n).
Portanto a == f  e  b = (e - f(m/n)) são os inteiros desejados. Por exemplo, para
  mdc(3267,2893)
    mdc(2893,374)
      mdc(374,275)
        mdc(275,99)
          mdc(99,77)
            mdc(77,22)
              mdc(22,11)
                mdc(11,0)
  mdc de 3267 e 2893 e' 11
e temos que
  11 ==    1 x 11   +    0 x 0   ==   0 x 22  +  1 x 11
     ==    1 x 77   +   -3 x 22  ==  -3 x 99  +  4 x 77
     ==    4 x 275  +  -11 x 99  == -11 x 374 + 15 x 275
     =  = 15 x 2893 + -116 x 374
     == -116 x 3267 +  131 x 2893.

  Eis uma versão em C da descrição acima

 /*
  * Ao receber dois números inteiros não-negativos m e n, não
  * ambos nulos, está função devolve o máximo divisor de m e n
  * e nos parâmetros *a e *b devolve números inteiros
  * tais que
  *           (*a)*m + (*b)*n == mdc(m,n).
  */  
  long
  euclides_r_x(long m, long n, long *a, long *b)
  {
     long d;
    
     if (n==0) 
       { 
         *a = 1;
         *b = 0; 
         return m;
       }
     else 
       {
         long e, f;

         d = euclides_r_x (n, m % n, &e, &f);

        /* Neste ponto vale que
              d == e*n + f*(m%n)
                == e*n + f*(m - (m/n)*n)
                == f*m + (e - f*(m/n))*n   */

         *a = f;
         *b = e - f*(m/n);

        return d;
      }
  
    }

  Mas qual a razão de tanta confusão? A principal razão (pelo menos neste ponto da disciplina) é que de uma certa maneira os números a e b são uma prova ou certificado da correção da resposta. Suponha que nosso programa afirme que o máximo divisor comum de dois números gigantescos m e n é d e alám disso nos fornece números inteiros a e b tais que

am + bn == d.
Suponha ainda ques estamos meio deconfiados da correção da resposta. Será que pode existir um inteiro positivo d' maior que d e que divida m e n? A resposta é não! Se existe d' < d que divida m e n, então d' divide am + bn == d, e portanto d' <= d.

  Para verificar que a resposta do nosso programa está correta, basta verificarmos que

am + bn == d,
e que d divide m e n. Esta verificação custa 2 multiplicações, 2 divisões, 1 adição e 3 comparações.

  Algoritmos que devolvem certificados da sua correção sao chamados self-certifying ou simpesmente certifying. Na página pessoal de Kurt Mehlhorn há um link para uma palestra sobre Certifying Algorithms.

 

 

 


Last modified: Thu Mar 11 10:47:28 BRT 2004