Egyiptomi szorzás
Ismert történelem: Nílus áradása után újra ki kell jelölni a határokat.
Papok egy kiválasztott csoportja, a "kötélfeszítők" foglalkoztak evvel.
Kevés írásos emlék maradt fent, az egyik a Rhind matematikai papirusz
(Rhind skót gyüjtő volt), i.e.1650-ből, amit Ahmes írnok jegyzett le.
A tekercs tartalmaz szorzási és osztási algoritmust is.
Az egyiptomiak nem használtak helyiértékeket és nem imerték a nullát.
Ennek eredményeképpen a műveletek, pl. a szorzás nehéz volt és csak
kevesen tudták végrehajtani.
A szorzás: ismételt összeadás.
1a = a
(n+1)a = na + a
In C++:
#include <iostream>
#include <cstdlib>
int mul0( int n, int a)
{
if ( 1 == n ) return a;
return mul0(n-1,a) + a;
}
int main(int argc, char *argv[])
{
if ( argc == 3)
{
int n = atoi(argv[1]);
int a = atoi(argv[2]);
std::cout << mu10(n,a) << std::endl;
}
return 0;
}
Természetesen a>0 és n>0, az egyiptomiak tudomásunk szerint nem használtak
negatív számokat.
Az egyiptomi szorzás (Knuth: orosz paraszt algoritmus) Ahmes példája nyomán:
n = 41 a = 59
1 # 59 a
2 118 (a+a)
4 236 (a+a) + (a+a)
8 # 472
16 944
32 # 1888
41 x 59 = (1x59)+(8x59)+(32x59)
Figyeljük meg a #-jellel jelölt számokat: ez a 41 bináris reprezentációja:
00000000 000101001 (16 biten)
Az algoritmus ellenőrzi, hogy "n" páros vagy páratlan, ezért feltételezhetjük,
bár biztosan nem tudjuk, hogy az egyiptomiak ismerték a páros és páratlan szám
fogalmát. A görögök ismerték.
XXXXXXXXXX -> XXXXX XXXXX páros
XXXXXXXXXXX -> XXXXX X XXXXX páratlan
#include <iostream>
#include <cstdlib>
bool odd(int n) { return n & 0x1; }
int half(int n) { return n >> 1; }
int mul1( int n, int a)
{
if ( 1 == n ) return a;
int result = mul1(half(n),a+a);
if ( odd(n) ) result = result + a;
return result;
}
int main(int argc, char *argv[])
{
if ( argc == 3)
{
int n = atoi(argv[1]);
int a = atoi(argv[2]);
std::cout << mul1(n,a) << std::endl;
}
return 0;
}
Ennek az algoritmusnak a költsége = log(n) + (v(n)-1)
ahol v(n) az 1-esek száma (population count) n bináris ábrázolásában.
- Megj1: 3 x 15 könnyebben számítható, mint 15 x 3, ez könnyű optimizálás
- Megj2: A függvényhívások költségesek
Dolgozzunk kicsit többet, számítsuk ki "na" helyett "r + na" -t.
Gyakran már hardver szinten erre a "multiply-accumulate"-re optimalizálnak.
#include <iostream>
#include <cstdlib>
bool odd(int n) { return n & 0x1; }
int half(int n) { return n >> 1; }
int mulacc0( int r, int n, int a)
{
if ( 1 == n ) return r + a;
if ( odd(n) )
return mulacc0(r+a, half(n), a+a);
else
return mulacc0(r, half(n), a+a);
}
int main(int argc, char *argv[])
{
if ( argc == 3)
{
int n = atoi(argv[1]);
int a = atoi(argv[2]);
std::cout << mulacc0(0,n,a) << std::endl;
}
return 0;
}
Kis átalakítással:
int mulacc1( int r, int n, int a)
{
if ( 1 == n ) return r + a;
if ( odd(n) )
r = r + a;
return mulacc1(r, half(n), a+a);
}
Továbbá,
- n ritkán 1
- ha n páros, akkor nem kell ellenőrizni, hogy 1-e
int mulacc2( int r, int n, int a)
{
if ( odd(n) )
{
r = r + a;
if ( 1 == n ) return r;
}
return mulacc2(r, half(n), a+a);
}
Def: egy algoritmus szigorúan vég-rekurzív, ha a vég-rekurzív hívása
pont ugyanazokkal a formális paraméterekkel történik.
int mulacc2( int r, int n, int a)
{
if ( odd(n) )
{
r = r + a;
if ( 1 == n ) return r;
}
n = half(n);
a = a + a;
return mulacc2(r, n, a);
}
Ez már könnyen alakítható ciklussá:
int mulacc3( int r, int n, int a)
{
while (true )
{
if ( odd(n) )
{
r = r + a;
if ( 1 == n ) return r;
}
n = half(n);
a = a + a;
}
}
Mérjünk!
#include <iostream>
#include <cstdlib>
int eqcounter = 0;
int pluscounter = 0;
int minuscounter = 0;
struct MyInt
{
MyInt( int i ) : value(i) {}
int value;
};
MyInt operator+(MyInt a, MyInt b)
{
++pluscounter;
return MyInt(a.value + b.value);
}
MyInt operator-(MyInt a, int b)
{
++minuscounter;
return MyInt(a.value - b);
}
bool operator==(int a, MyInt b)
{
++eqcounter;
return a == b.value;
}
std::ostream& operator<<( std::ostream& os, MyInt n)
{
os << n.value;
return os;
}
bool odd(MyInt n) { return n.value & 0x1; }
MyInt half(MyInt n) { return MyInt(n.value >> 1); }
MyInt mulacc3( MyInt r, MyInt n, MyInt a)
{
while (true )
{
if ( odd(n) )
{
r = r + a;
if ( 1 == n ) return r;
}
n = half(n);
a = a + a;
}
}
MyInt multiply( MyInt n, MyInt a)
{
while ( !odd(n) )
{
a = a + a;
n = half(n);
}
if ( 1 == n ) return a;
return mulacc3( a, n-1, a);
}
int main(int argc, char *argv[])
{
if ( argc == 3)
{
MyInt n = atoi(argv[1]);
MyInt a = atoi(argv[2]);
std::cout << multiply(n,a) << std::endl;
std::cout << "#plus = " << pluscounter << std::endl;
std::cout << "#minus = " << minuscounter << std::endl;
std::cout << "#eq = " << eqcounter << std::endl;
}
return 0;
}
$ ./a.out 7 8
56
#plus = 4
#minus = 1
#eq = 3
# elso verzio
$ ./a.out 76543 21098
1614904214
#plus = 76542
#minus = 76542
#eq = 76543
# kozbulso verzio
$ ./a.out 76543 21098
1614904214
#plus = 27
#minus = 0
#eq = 17
#utolso verzio
$ ./a.out 76543 21098
1614904214
#plus = 27
#minus = 1
#eq = 12
$