Je dispose de
1- X entier codé sur n bits
2 - Y entier codé sur n bits

je sais de plus que la norme du vecteur (X,Y) peut être codée sur n bits
(par contre X^2 tout comme Y^2 peut dépasser n bits!)

quelqu'un aurrait-il/elle un algorithme éfficace pour calculer la valeur entière de
sqrt(X^2 + Y^2) Sans passer par des floatting point

pour info n=40 dans mon cas, le nombre de bits effectivement utiliés va de 20 a 38 ou 39.

Actuellement j'ai créé des buffers de 80 bits pour stocker X^2 et Y^2
un buffer de 81 bits pour la somme et je cacule le sqrt avec la méthode
"de la potence " adaptée à la base 2.

Cela me semble un peu lourd!!


Merci de vos réponses

Voici mon code qui fonctionne mais me semble non optimal

Code C : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
 
#define Nbits   40
unsigned long long L1, L2;
unsigned short  W1[Nbits], W2[Nbits],
                W1_sqr[2*Nbits], W2_sqr[2*Nbits],
                W_sm[2*(Nbits+1)]; // en principe 2*Nbits+1 suffisent mais un nombre paire est + facile a programmer dans 'pendule'
void Bit40( unsigned long long L, unsigned short * W)
   {
   short i;
   for (i=0; i < Nbits; i++)
      {
      W[i] = L & 1;
      L = L>> 1;
      }
   }
void LSQR( unsigned short * W,unsigned short * W2)
   {
   short i;
   memset(W2,0,2*Nbits*sizeof(short));
   for (i=0; i< Nbits; i++) if ( W[i] )
      {
      short retenue =0;
      short k;
      for (k=0; k < Nbits; k++)
         {
         short u= W2[i+k] + W[k] + retenue;
         if ( u <= 1)
            {
            W2[i+k] = u;
            retenue = 0;
            }
         else
            {
            W2[i+k] = u - 2;
            retenue = 1;
            }
         }
      }
   }
void Add2( unsigned short * W1_sqr, unsigned short * W2_sqr,unsigned short * W_sm)
   {
   short k,   retenue = 0;
   memset(W_sm,0,2*(Nbits+1)*sizeof(short));
   for (k=0; k < 2*Nbits; k++)
      {
      short m = retenue+W1_sqr[k] + W2_sqr[k] ;
      if ( m <=1)
         {
         W_sm[k] = m;
         retenue = 0;
         }
      else
         {
         W_sm[k] = m-2;
         retenue = 1;
         }
      }
   W_sm[2*Nbits]=retenue; //W_sm[2*Nbits+1] est toujours nul
   }
void potence (unsigned short *W_sm, unsigned long long * square, unsigned long long * reste)
   {
   *square=0;
   * reste =0;
   short i = Nbits;
   while (1)
      {
      long long L;
      short j=2*i+1;
      short k=2*i;
      * reste = (* reste) * 4  + (W_sm[j] << 1) +  W_sm[k];
      L = 4* (*square) + 1;
      if ( L <= (* reste) )
         {
         (* reste) -= L;
         (*square) = 2*(*square)  +1;
         }
      else
         (* square) =  2 * ( * square);
      i--;
      if ( i < 0) break;
      }
   }
void  Norme(unsigned long long L1, unsigned long long L2, unsigned long long * L, unsigned long long * R)
   {
   Bit40(L1,W1);
   Bit40(L2,W2);
   LSQR(W1,W1_sqr);
   LSQR(W2,W2_sqr);
   Add2(W1_sqr,W2_sqr,W_sm);
   potence (W_sm,  L,R);
   }

Note: J'utilise les unsigned short car au final le SW doit aller sur un DSP de type C5510 qui ne connait pas d'entité de 8 bits.