Arageli: Library for Doing Exact Computation Segrey S. Lyalin - - PowerPoint PPT Presentation
Arageli: Library for Doing Exact Computation Segrey S. Lyalin Nikolai Yu. Zolotykh University of Nizhni Novgorod, Russia ICMS, 2006 Contents 1 Project 3 2 Classes 4 3 Algorithms 6 4 Examples 8 5 Aragelis (distinctive)
University of Nizhni Novgorod, Russia ICMS, 2006
1 Project 3 2 Classes 4 3 Algorithms 6 4 Examples 8 5 Arageli’s (distinctive) characteristics 29 6 Comparison with other libraries 30 7 Arageli in action 31
Arageli is a C++ library for exact computations in ARithmetic, Algebra, GEometry, Linear and Integer linear programming Sphere of possible applications: exact linear algebra, number theory, linear and integer linear programming, criptography, symbolic computations 100+ files of source code, 50000+ lines
Arageli is a C++ library for exact computations in ARithmetic, Algebra, GEometry, Linear and Integer linear programming Sphere of possible applications: exact linear algebra, number theory, linear and integer linear programming, criptography, symbolic computations 100+ files of source code, 50000+ lines
Eugene Agafonov Max Alekseyev (Uni of San Diego at pres.) Aleksey Bader Sergey S. Lyalin Aleksey Polovinkin Andrey Somsikov Nikolai Yu. Zolotykh University of Nizhni Novgorod
Arbitrary precision integers
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Vectors with entries that belong to T
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Vectors with entries that belong to T
Matrices with entries that belong to T
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Vectors with entries that belong to T
Matrices with entries that belong to T
Sparse polynomials with coefficients that belong to T
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Vectors with entries that belong to T
Matrices with entries that belong to T
Sparse polynomials with coefficients that belong to T
Elements in T modulo certain m ∈ T
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Vectors with entries that belong to T
Matrices with entries that belong to T
Sparse polynomials with coefficients that belong to T
Elements in T modulo certain m ∈ T
Algebraic numbers, i.e. roots of polynomials with coefficients in T = {rational numbers}
Arbitrary precision integers
Rational numbers/functions, i.e. fractions u v , where u, v ∈ T (integers/polynomials)
Vectors with entries that belong to T
Matrices with entries that belong to T
Sparse polynomials with coefficients that belong to T
Elements in T modulo certain m ∈ T
Algebraic numbers, i.e. roots of polynomials with coefficients in T = {rational numbers}
algebraic numbers, multiprecision float point numbers, . . . )
algebraic numbers, multiprecision float point numbers, . . . ) All classes above (except big int) have more than 1 parameter. Supplementary parameters control reducing fractions algorithms, references count methods, etc.
field and integer domains
problems
int a, b; big int c, d; polynomial<rational<> > f, g; . . . . . . . . . . . . . . . . . . . . . . . . . . . . . gcd(a, b); gcd(c, d); gcd(a, c); gcd(f, g);
int a, b; big int c, d; polynomial<rational<> > f, g; . . . . . . . . . . . . . . . . . . . . . . . . . . . . . gcd(a, b); gcd(c, d); gcd(a, c); gcd(f, g);
working, puts out intermediate results of the computation, measures time etc.
Example: big integers
#include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { big int a = “101100111000111100001111100000”; big int b = “-1234567890987654321”; cout << “a = ” << a << endl << “b = ” << b << endl; cout << “2*(a + b) - a%b = ” << 2*(a + b) − a%b << endl; return 0; } a = 101100111000111100001111100000 b = -1234567890987654321 2*(a + b) - a%b = 202200221996782640900764076427
Example: rational numbers
#include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { rational<> a = “22/7”; rational<> b(355, 113); cout << “a = ” << a << endl; cout << “b = ” << b << endl; cout << “Numerator of a = ” << a.numerator() << endl; cout << “Denominator of a = ” << a.denominator() << endl; double af = a; // You can assign rational numbers to float or to double double bf = b; cout << “af = ” << setprecision(5) << af << endl; cout << “bf = ” << setprecision(10) << bf << endl; return 0; }
a = 22/7 b = 355/113 Numerator of a = 22 Denominator of a = 7 af = 3.1429 bf = 3.14159292
Example: polynomials
#include “arageli/arageli.hpp” using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { sparse polynom<rational<> > f = “x^6+x^4-x^2-1”, g = “x^3-2*x^2-x+2”; cout << “f = ” << f << endl; cout << “g = ” << g << endl << endl; cout << “f + g = ” << f + g << endl; cout << “f - g = ” << f − g << endl; cout << “f * g = ” << f * g << endl; cout << “f / g = ” << f / g << endl; cout << “f % g = ” << f % g << endl; cout << “GCD(f, g) = ” << gcd(f, g) << endl; cout << “LCM(f, g) = ” << lcm(f, g) << endl; return 0; }
f = x^6+x^4-x^2-1 g = x^3-2*x^2-x+2 f + g = x^6+x^4+x^3-3*x^2-x+1 f - g = x^6+x^4-x^3+x^2+x-3 f * g = x^9-2*x^8-2*x^5+4*x^4+x-2 f / g = x^3+2*x^2+6*x+12 f % g = 25*x^2-25 GCD(f, g) = x^2-1 LCM(f, g) = x^7-2*x^6+x^5-2*x^4-x^3+2*x^2-x+2
Example: vectors and matrices
#include <arageli/arageli.hpp> using namespace std; using namespace Arageli; typedef rational<> Q; int main(int argc, char *argv[ ]) { matrix<Q> A = “((-1/2, 3/4), (-2/3, 5), (1/7, -5/2))”; matrix<Q> B = “((3/4, 1/6, -7/8), (5/2, 2/5, -9/10))”; vector<Q> c = “(1/4, -4/15, 5)”; vector<Q> d = “(-2/3, -1, 4)”; Q alpha = Q(1, 120), beta = −2; vector<Q> res = ((A*B)*c + d*(A*B) − alpha*d)/beta; // (A*B)*c → c is interpreted as a column // d*(A*B) → d is interpreted as a row cout << “Result:” << endl;
return 0;
} Result: ||968591/50400|| ||174007/15120|| ||-27583/2520 ||
Example: Smith’s normal diagonal form
Integer matrices #include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { matrix<big int> A, B, P, Q; size t rk; sparse polynom<rational<> > d; A = “((0,6,-9), (12,24,9), (30,42,45))”; smith(A, B, P, Q, rk, d); // B = PAQ cout << “A = ” << endl; output aligned(cout, A); cout << “B = ” << endl; output aligned(cout, B); cout << “P = ” << endl; output aligned(cout, P); cout << “Q = ” << endl; output aligned(cout, Q);
cout << “d = ” << d << endl; return 0; } A = ||0 6
||12 24 9 || ||30 42 45|| B = ||3 0 0 || ||0 6 0 || ||0 0 18|| P = ||-1 0 0|| ||11 1 0|| ||18 -1 1|| Q = ||0 -7 15|| ||1 3
||1 2
d = 324
Polynomial matrices #include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { matrix< sparse polynom<rational<> > > A, B, P, Q; size t rk; sparse polynom<rational<> > d; A = “((4-x,6,-15), (1,3-x,-5), (1,2,-4-x))”; smith(A, B, P, Q, rk, d); // B = PAQ cout << “A = ” << endl; output aligned(cout, A); cout << “B = ” << endl; output aligned(cout, B); cout << “P = ” << endl; output aligned(cout, P); cout << “Q = ” << endl; output aligned(cout, Q); cout << “d = ” << d << endl;
return 0; } A = ||-x+4 6
|| 1
|| 1 2
B = ||1 || ||0 x-1 || ||0 x^2-2*x+1|| P = ||-1/15 0 0 || || 1/3
||-x-2
Q = ||0 1 || ||0 1 1/3 || ||1 2/5 -1/15*x+2/5|| d = -x^3+3*x^2-3*x+1
Example: modular arithmetic
#include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { { // Working in Z5 = Z/5Z residue<int> a = “2 (mod 5)”; residue<int> b = “3 (mod 5)”; cout << “a = ” << a << endl; cout << “b = ” << b << endl; cout << “a + b = ” << a + b << endl; cout << “a - b = ” << a − b << endl; cout << “a * b = ” << a * b << endl; cout << “a / b = ” << a / b << endl << endl; }
{ // Working in Q[x]/(x2 + 1)Q[x] typedef residue<sparse polynom<rational<> > > T; T a = “((1+x) (mod (x^2+1)))”; T b = “((1-x) (mod (x^2+1)))”; cout << “a = ” << a << “b = ” << b << endl; cout << “a + b = ” << a + b << endl; cout << “a - b = ” << a − b << endl; cout << “a * b = ” << a * b << endl; cout << “a / b = ” << a / b << endl << endl; } { // Working in Z3[x]/(x2 + 1)Z3[x] typedef residue<sparse polynom<residue<int> > > T; T a = “( ((1(mod 3))*x+(1(mod 3))) (mod ((1(mod 3))*x^2+(1(mod 3)))) )”; T b = “( ((1(mod 3))*x+(2(mod 3))) (mod ((1(mod 3))*x^2+(1(mod 3)))) )”; cout << “a = ” << a << endl << “b = ” << b << endl; cout << “a + b = ” << a + b << endl; cout << “a - b = ” << a − b << endl; cout << “a * b = ” << a * b << endl; cout << “a / b = ” << a / b << endl << endl; }
} a = 2(mod 5) b = 3(mod 5) a + b = 0(mod 5) a - b = 4(mod 5) a * b = 1(mod 5) a / b = 4(mod 5) a = x+1(mod x^2+1)b = -x+1(mod x^2+1) a + b = 2(mod x^2+1) a - b = 2*x(mod x^2+1) a * b = 2(mod x^2+1) a / b = x(mod x^2+1) a = x+1(mod +3)(mod x^2+1(mod +3)) b = x+2(mod +3)(mod x^2+1(mod +3)) a + b = 2(mod 3)*x(mod x^2+1(mod +3)) a - b = 2(mod 3)(mod x^2+1(mod +3)) a * b = 1(mod 3)(mod x^2+1(mod +3)) a / b = 2(mod 3)*x(mod x^2+1(mod +3))
Example: linear system over finite field
#include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { typedef residue<int> T; int mod = 7; // Solving a linear system Ax = b over Z7 matrix<T> A = “((3, 1, 2), (1, 2, 3), (4, 3, 2))”; vector<T> b = “(1, 1, 1)”; vector<T> x; for(matrix<T>::iterator i = A.begin(); i < A.end(); ++i) { i− >module() = mod; i− >normalize(); }
for(matrix<T>::iterator i = b.begin(); i < b.end(); ++i) { i− >module() = mod; i− >normalize(); } cout << “A = ” << endl; output aligned(cout, A); cout << endl; cout << “b = ” << endl; output aligned(cout, b); cout << endl; try { x = solve linsys(A, b); cout << “x = ” << endl;
cout << endl; cout << “Check the result: ” << boolalpha << (A*x == b) << endl; } catch(matrix is singular) { cout << “Error! Matrix is singular!” << endl; } return 0;
} A = ||3(mod 7) 1(mod 7) 2(mod 7)|| ||1(mod 7) 2(mod 7) 3(mod 7)|| ||4(mod 7) 3(mod 7) 2(mod 7)|| b = ||1(mod 7)|| ||1(mod 7)|| ||1(mod 7)|| x = ||2(mod 7)|| ||6(mod 7)|| ||5(mod 7)|| Check the result: true
Example: using iterators
#include <arageli/arageli.hpp> using namespace std; using namespace Arageli; int main(int argc, char *argv[ ]) { typedef sparse polynom<big int> poly; typedef poly::degree iterator degrees; poly S = “213*x^3443+532*x^4432-744*x^44-4235*x^15+292*x+34254”; cout << “S = ” << S << endl << endl; // Replace all degrees by their residues modulo 5 for(degrees di = S.degrees begin(), dj = S.degrees end(); di != dj; ++di) *di %= 5; // Now S can contain monomials this equal degrees // We have to reduce S to its canonical form S.normalize(); cout << “All degrees modulo 5 = ” << S << endl << endl; return 0;
} S = 532*x^4432+213*x^3443-744*x^44-4235*x^15+292*x+34254 All degrees modulo 5 = -744*x^4+213*x^3+532*x^2+292*x+30019
Example: matrix polynomial → polynomial matrix conversion
#include <arageli/arageli.hpp> // 3 −5 7
1 −8
−1 9 13 2
// 3x8 − 1 −5x8 + x3 + 9 13 7x8 − 8x3 + 2
using namespace Arageli; int main(int argc, char *argv[ ]) { sparse polynom<matrix<int> > f = “(((3,-5),(0,7))*x^8+((0,1),(0,-8))*x^3+((-1,9),(13,2)))”; matrix<sparse polynom<int> > x = “((x,0),(0,x))”; cout << “f(x) = ” << f << endl << endl; cout << “Let x = ” << x << endl << endl;
cout << “Then the resulting matrix is” << endl;
return 0; } f(x) = ((3, -5), (0, 7))*x^8+((0, 1), (0, -8))*x^3+((-1, 9), (13, 2)) Let x = ((x, 0), (0, x)) Then the resulting matrix is || 3*x^8-1 -5*x^8+x^3+9 || || 13 7*x^8-8*x^3+2 ||
All functions and operators that have more than one argument can take different types for similar arguments. For example, operator+ for two vectors can perform elementwise addition on vectors with different types of elements. There are several rules (that automatically verifies by a compiler) that define more precisely the meaning of similarity of argument types and allowable type propagation.
All functions and operators that have more than one argument can take different types for similar arguments. For example, operator+ for two vectors can perform elementwise addition on vectors with different types of elements. There are several rules (that automatically verifies by a compiler) that define more precisely the meaning of similarity of argument types and allowable type propagation.
This allows to use our code with user-defined classes.
All functions and operators that have more than one argument can take different types for similar arguments. For example, operator+ for two vectors can perform elementwise addition on vectors with different types of elements. There are several rules (that automatically verifies by a compiler) that define more precisely the meaning of similarity of argument types and allowable type propagation.
This allows to use our code with user-defined classes.
All functions and operators that have more than one argument can take different types for similar arguments. For example, operator+ for two vectors can perform elementwise addition on vectors with different types of elements. There are several rules (that automatically verifies by a compiler) that define more precisely the meaning of similarity of argument types and allowable type propagation.
This allows to use our code with user-defined classes.
– gcc 3.3.3, 4.0.4 – Microsoft Visual Studio C++ 2003, 2005 – Intel C++ 8.0, 9.0
“High-performance, portable C++ library providing data structures and algorithms for arbitrary length integers; for vectors, matrices, and polynomials over the integers and over finite fields; and for arbitrary precision floating point arithmetic”
“High-performance, portable C++ library providing data structures and algorithms for arbitrary length integers; for vectors, matrices, and polynomials over the integers and over finite fields; and for arbitrary precision floating point arithmetic”
at Darmstadt “A library for computational number theory”
“High-performance, portable C++ library providing data structures and algorithms for arbitrary length integers; for vectors, matrices, and polynomials over the integers and over finite fields; and for arbitrary precision floating point arithmetic”
at Darmstadt “A library for computational number theory” – kernel level (multiple precision integers, memory low-level routines) – interfaces level (standardaized functions for multiple precision integers and memory-managment functions) – simple (not parameterized) classes (bigrational, bigfloat,...) – parameterized classes (containers) (vector, matrix, power series, ...) – user inrefaces level (online documentation tool, interpreter)
“High-performance, portable C++ library providing data structures and algorithms for arbitrary length integers; for vectors, matrices, and polynomials over the integers and over finite fields; and for arbitrary precision floating point arithmetic”
at Darmstadt “A library for computational number theory” – kernel level (multiple precision integers, memory low-level routines) – interfaces level (standardaized functions for multiple precision integers and memory-managment functions) – simple (not parameterized) classes (bigrational, bigfloat,...) – parameterized classes (containers) (vector, matrix, power series, ...) – user inrefaces level (online documentation tool, interpreter)
“a Class Library for Numbers”
“High-performance, portable C++ library providing data structures and algorithms for arbitrary length integers; for vectors, matrices, and polynomials over the integers and over finite fields; and for arbitrary precision floating point arithmetic”
at Darmstadt “A library for computational number theory” – kernel level (multiple precision integers, memory low-level routines) – interfaces level (standardaized functions for multiple precision integers and memory-managment functions) – simple (not parameterized) classes (bigrational, bigfloat,...) – parameterized classes (containers) (vector, matrix, power series, ...) – user inrefaces level (online documentation tool, interpreter)
“a Class Library for Numbers”
“An open framework for symbolic computation within the C++ programming language”
New implementation of double description method
There are two possible description of a polyhedron P: P = {x ∈ Rd : Ax ≤ b}, where A ∈ Zm×d (1) P =
s
αiui +
n
βivi, ui ≥ 0, vi ≥ 0,
n
α = 1
(1) → (2) facet enumeration problem (2) → (1) vertex enumeration problem Double description method (DDM) [Motzkin, Raiffa, Thompson, Thrall, 1953] is one of the methods for solving problems (1) ↔ (2) cdd: implementation of DDM by Komei Fukuda [Fukuda, 1996] Skeleton 2: implementation of DDM [Zolotykh, 2006]. It implement a new variation of DDM and uses Arageli.
Experimental results
Problem Input Output Time, sec m d p cdd Skeleton 2 cube18 18 36 262144 61 5 mit729-9 8 729 4862 55 6 prodmT62 24 3461 168 95 9 ccp7 21 64 116754 4760 61 ccc7 21 63 38780 1424 814 The experiment was performed on Intel Pentium D, Clock Speed 2.81GHz, L2cash 1Mb, 2Gb RAM. The program was compiled by C++ MS Visual Studio 2005 compiler with /o2 key. The test problems are from Fukuda’s repository.
www.unn.ru/cs/arageli
Sergey.Lyalin@itlab.unn.ru Zolotykh@vmk.unn.ru