80 views

Skip to first unread message

Apr 25, 2021, 9:07:21 AMApr 25

to sage-devel

Hi all,

I'm looking for benchmark problems for QQbar/AA arithmetic. Ideally such a problem will:

* Be reducible to a short program that, apart from using QQbar/AA operations, is reasonably self-contained.

* Reflect real-world use, i.e. originate from using Sage to solve an actual mathematical problem, with QQbar/AA arithmetic being a significant ingredient.

I'm interested in such benchmark problems to evaluate how well Calcium (https://fredrikj.net/calcium/) would work to implement the field of algebraic numbers. As such, I'd like to identify any situations where Calcium is slower than Sage's QQbar or where functionality missing.

There are actually two independent types that may serve to implement the field of algebraic numbers: qqbar_t and ca_t. A Calcium qqbar_t uses an "absolute" representation (minimal polynomial), which is nice since results are kept reduced and canonical, but can be slow for field arithmetic. A Calcium ca_t uses a "relative" representation involving symbolic expressions and number fields, somewhat similar to Sage's QQbar or Magma's algebraically closed field. The ca_t type should give the best performance in most cases, but qqbar_t may be better in some circumstances.

Some benchmark problems I've collected so far:

huge_expr (https://ask.sagemath.org/question/52653/)

====================================================

sage qqbar: ~forever (hours)

calcium qqbar: 8.5 s

calcium ca: 8.1 s

heptadecagon (sage/rings/qqbar.py)

====================================================

sage qqbar: ~forever (hours)

calcium qqbar: 2.11 s

calcium ca: 0.011 s

heptadecagon2 (sage/rings/qqbar.py)

====================================================

sage qqbar: 0.060 s

calcium qqbar: 2.11 s

calcium ca: 0.0031 s

symbench_R1 (https://wiki.sagemath.org/symbench)

====================================================

sage qqbar (numerical answer): 0.0041 s

sage qqbar (exactify): 0.049 s

Fredrik

I'm looking for benchmark problems for QQbar/AA arithmetic. Ideally such a problem will:

* Be reducible to a short program that, apart from using QQbar/AA operations, is reasonably self-contained.

* Reflect real-world use, i.e. originate from using Sage to solve an actual mathematical problem, with QQbar/AA arithmetic being a significant ingredient.

I'm interested in such benchmark problems to evaluate how well Calcium (https://fredrikj.net/calcium/) would work to implement the field of algebraic numbers. As such, I'd like to identify any situations where Calcium is slower than Sage's QQbar or where functionality missing.

There are actually two independent types that may serve to implement the field of algebraic numbers: qqbar_t and ca_t. A Calcium qqbar_t uses an "absolute" representation (minimal polynomial), which is nice since results are kept reduced and canonical, but can be slow for field arithmetic. A Calcium ca_t uses a "relative" representation involving symbolic expressions and number fields, somewhat similar to Sage's QQbar or Magma's algebraically closed field. The ca_t type should give the best performance in most cases, but qqbar_t may be better in some circumstances.

Some benchmark problems I've collected so far:

huge_expr (https://ask.sagemath.org/question/52653/)

====================================================

sage qqbar: ~forever (hours)

calcium qqbar: 8.5 s

calcium ca: 8.1 s

heptadecagon (sage/rings/qqbar.py)

====================================================

sage qqbar: ~forever (hours)

calcium qqbar: 2.11 s

calcium ca: 0.011 s

heptadecagon2 (sage/rings/qqbar.py)

====================================================

sage qqbar: 0.060 s

calcium qqbar: 2.11 s

calcium ca: 0.0031 s

symbench_R1 (https://wiki.sagemath.org/symbench)

====================================================

sage qqbar (numerical answer): 0.0041 s

sage qqbar (exactify): 0.049 s

calcium qqbar: 0.095 s

calcium ca (numerical answer): 0.0010 s

calcium ca (qqbar output) 0.057 sFredrik

Apr 25, 2021, 10:26:18 AMApr 25

to sage-...@googlegroups.com

Dear Fredrik,

One technical question: I thought that your ca_t implementation

used multivariate polynomials. This is what Magma does but not

what sage does. The latter uses expression trees and take union

fields anytime there is an exactification to be done. Is it a

misconception of mine?

For "real world" problems

1) algebraic solutions of zero dimensional varieties

defined over QQ

sage: R.<x,y,z> = QQ[]

sage: A = x^4 + y^4 + z^3 - 1

sage: B = x^2*z + 2*y^2*x - x*y*z + 5

sage: C = y*z^2 - 2*(x+y+z) + 2

sage: %time solutions = R.ideal([A,B,C]).variety(QQbar)

CPU times: user 24.3 s, sys: 23.3 ms, total: 24.3 s

Wall time: 24.4 s

Though here is the bottleneck is Groebner basis computation

of multivariate polynomials and not really QQbar (21s over the 24s).

2) Two examples somehow similar to your huge_expr

sage: a = sum(AA(p).sqrt() for p in prime_range(18))

sage: %time a.exactify()

CPU times: user 30.7 s, sys: 13.3 ms, total: 30.8 s

Wall time: 30.8 s

sage: a.degree() # does not terminates

In the above situation, the problem is that sage is representing

the number a as a complicated element in a rather complicated

number field. Computing the minpoly involves apparently too

complicated linear algebra (?). However, if instead of exactifying

in that way but using the minimal polynomial (which looks reasonable

in this range of degrees) the latter would be very quick.

sage: b = sum(QQbar.zeta(n) for n in range(1,8))

sage: b.exactify() # does not terminate

Here sage is stuck at building union fields.

Would be nice to have calcium used in sage!

Best

Vincent

One technical question: I thought that your ca_t implementation

used multivariate polynomials. This is what Magma does but not

what sage does. The latter uses expression trees and take union

fields anytime there is an exactification to be done. Is it a

misconception of mine?

For "real world" problems

1) algebraic solutions of zero dimensional varieties

defined over QQ

sage: R.<x,y,z> = QQ[]

sage: A = x^4 + y^4 + z^3 - 1

sage: B = x^2*z + 2*y^2*x - x*y*z + 5

sage: C = y*z^2 - 2*(x+y+z) + 2

sage: %time solutions = R.ideal([A,B,C]).variety(QQbar)

CPU times: user 24.3 s, sys: 23.3 ms, total: 24.3 s

Wall time: 24.4 s

Though here is the bottleneck is Groebner basis computation

of multivariate polynomials and not really QQbar (21s over the 24s).

2) Two examples somehow similar to your huge_expr

sage: a = sum(AA(p).sqrt() for p in prime_range(18))

sage: %time a.exactify()

CPU times: user 30.7 s, sys: 13.3 ms, total: 30.8 s

Wall time: 30.8 s

sage: a.degree() # does not terminates

In the above situation, the problem is that sage is representing

the number a as a complicated element in a rather complicated

number field. Computing the minpoly involves apparently too

complicated linear algebra (?). However, if instead of exactifying

in that way but using the minimal polynomial (which looks reasonable

in this range of degrees) the latter would be very quick.

sage: b = sum(QQbar.zeta(n) for n in range(1,8))

sage: b.exactify() # does not terminate

Here sage is stuck at building union fields.

Would be nice to have calcium used in sage!

Best

Vincent

Apr 25, 2021, 2:46:02 PMApr 25

to sage-devel

Hi Vincent,

On Sunday, April 25, 2021 at 4:26:18 PM UTC+2 vdelecroix wrote:

Dear Fredrik,

One technical question: I thought that your ca_t implementation

used multivariate polynomials. This is what Magma does but not

what sage does. The latter uses expression trees and take union

fields anytime there is an exactification to be done. Is it a

misconception of mine?

Yes, ca_t uses multivariate polynomials (in fact, multivariate fraction field elements). It uses univariate polynomials in the special case of simple (univariate) number fields.

Symbolic expressions are effectively used to represent non-field operations generating new extension elements (e.g. square roots, when not converted to absolute qqbar elements).

With qqbar_t:

In [2]: %time sum(qqbar(p).sqrt() for p in [2, 3, 5, 7, 11, 13, 17])

CPU times: user 286 ms, sys: 4.23 ms, total: 291 ms

Wall time: 290 ms

Out[2]: 19.0734 (deg 128)

In [3]: %time sum(qqbar(RootOfUnity(n)) for n in range(1,8))

CPU times: user 24.1 ms, sys: 30 µs, total: 24.1 ms

Wall time: 23.5 ms

Out[3]: 0.932507 + 4.46494*I (deg 96)

With ca_t:

In [2]: %time sum(ca(p).sqrt() for p in [2, 3, 5, 7, 11, 13, 17])

CPU times: user 2.77 ms, sys: 265 µs, total: 3.03 ms

Wall time: 2.63 ms

Out[2]: 19.0734 {a+b+c+d+e+f+g where a = 4.12311 [a^2-17=0], b = 3.60555 [b^2-13=0], c = 3.31662 [c^2-11=0], d = 2.64575 [d^2-7=0], e = 2.23607 [e^2-5=0], f = 1.73205 [f^2-3=0], g = 1.41421 [g^2-2=0]}

In [3]: %time sum(ca(qqbar(RootOfUnity(n))) for n in range(1,8))

CPU times: user 3.13 ms, sys: 0 ns, total: 3.13 ms

Wall time: 2.55 ms

Out[3]: 0.932507 + 4.46494*I {a+b+c+d where a = 0.623490 + 0.781831*I [a^6+a^5+a^4+a^3+a^2+a+1=0], b = 0.309017 + 0.951057*I [b^4+b^3+b^2+b+1=0], c = 1.73205*I [c^2+3=0], d = I [d^2+1=0]}

In [2]: %time sum(qqbar(p).sqrt() for p in [2, 3, 5, 7, 11, 13, 17])

CPU times: user 286 ms, sys: 4.23 ms, total: 291 ms

Wall time: 290 ms

Out[2]: 19.0734 (deg 128)

In [3]: %time sum(qqbar(RootOfUnity(n)) for n in range(1,8))

CPU times: user 24.1 ms, sys: 30 µs, total: 24.1 ms

Wall time: 23.5 ms

Out[3]: 0.932507 + 4.46494*I (deg 96)

With ca_t:

In [2]: %time sum(ca(p).sqrt() for p in [2, 3, 5, 7, 11, 13, 17])

CPU times: user 2.77 ms, sys: 265 µs, total: 3.03 ms

Wall time: 2.63 ms

Out[2]: 19.0734 {a+b+c+d+e+f+g where a = 4.12311 [a^2-17=0], b = 3.60555 [b^2-13=0], c = 3.31662 [c^2-11=0], d = 2.64575 [d^2-7=0], e = 2.23607 [e^2-5=0], f = 1.73205 [f^2-3=0], g = 1.41421 [g^2-2=0]}

In [3]: %time sum(ca(qqbar(RootOfUnity(n))) for n in range(1,8))

CPU times: user 3.13 ms, sys: 0 ns, total: 3.13 ms

Wall time: 2.55 ms

Out[3]: 0.932507 + 4.46494*I {a+b+c+d where a = 0.623490 + 0.781831*I [a^6+a^5+a^4+a^3+a^2+a+1=0], b = 0.309017 + 0.951057*I [b^4+b^3+b^2+b+1=0], c = 1.73205*I [c^2+3=0], d = I [d^2+1=0]}

Of course, it is debatable whether the result has really been "computed" here. If you convert the result to a qqbar_t, it will basically take as long as using qqbar_t in the first place. But these examples are trivial; in other cases, the multivariate arithmetic will be helpful.

Fredrik

Apr 25, 2021, 6:35:06 PMApr 25

to sage-...@googlegroups.com

On Sun, 2021-04-25 at 15:07 +0200, Fredrik Johansson wrote:

> Hi all,

>

> I'm looking for benchmark problems for QQbar/AA arithmetic. Ideally such a

> problem will:

>

> * Be reducible to a short program that, apart from using QQbar/AA

> operations, is reasonably self-contained.

> * Reflect real-world use, i.e. originate from using Sage to solve an actual

> mathematical problem, with QQbar/AA arithmetic being a significant

> ingredient.

I'm interested in Gram-Schmidt and the coordinate_vector() methods that
> Hi all,

>

> I'm looking for benchmark problems for QQbar/AA arithmetic. Ideally such a

> problem will:

>

> * Be reducible to a short program that, apart from using QQbar/AA

> operations, is reasonably self-contained.

> * Reflect real-world use, i.e. originate from using Sage to solve an actual

> mathematical problem, with QQbar/AA arithmetic being a significant

> ingredient.

solve linear systems over AA.

For a benchmark, you could conjure a set of long vectors with entries

in AA/QQbar, and orthonormalize them using the kindergarten algorithm

that everyone knows. That will leave a linearly-independent set, and if

you keep the input/output indices matched up, will let you know which

of the input vectors correspond to that linearly independent set. Call

the set of those input vectors "X". Then to benchmark the

coordinate_vector() methods, you can generate random elements in

W=span(X), and solve the system to find their W-coordinates in terms of

the basis X. (Using the de-orthonormalized vectors as the basis makes

the system harder to solve, and more representative of the general

case.)

Unimportant background information:

A Euclidean Jordan algebra (EJA) is a real vector space where every

element "x" has a self-adjoint "left multiplication by x" operation

that I'll write L_x.

It turns out that every EJA is a mishmash of various Hermitian matrix

algebras, which is nice, because it lets you think of their elements as

being some concrete object rather than just an abstract vector space

element. For example, the space of real symmetric matrices with

L_x(y) := (xy + yx)/2

induces an EJA. But when people think of a symmetric matrix, they're

thinking of it in terms of the standard basis, which leads to a

problem; namely that the basis

[1 0] [0 1] [0 0]

[0 0], [1 0], [0 1]

for S^2 is not normalized. This raises a problem: while L_x is always

self-adjoint, its matrix representation won't be symmetric unless your

basis is orthonormal. And you need its matrix to be symmetric if you

want to do anything cool with it in sage, because sage assumes that

when you have a matrix full of numbers, they're with respect to the

standard basis.

In short, this means that a sqrt(2) must become involved, and the "2"

isn't special. Here's where Gram-Schmidt comes in. Afterwards, matrices

aren't vectors in sage, so converting back and forth between

matrix/vector representations involves a lot of to_vector() and

coordinate_vector() over AA.

Things get worse: when constructing a subalgebra, its basis needs to be

orthonormalized as well, and we need to be able to convert back and

forth between the orthonormalized superalgebra coordinates and the

orthonormalized subalgebra coordinates. Thus we wind up

orthonormalizing vectors that already have ugly entries, and solving

systems with uglier and uglier entries in AA. This is a bottleneck when

the algebras get "big," which they do even in simple examples (the

Albert algebra has dimension 27).

Apr 26, 2021, 3:51:18 AMApr 26

to SAGE devel

+1 for Vincent's zero-dimensional variety point finding. I used to

use Magma's AlgebraicallyClosedField a lot for this (for some project

which I did before Sage existed). The idea is to find the points over

Qbar, then construct the absolute field they generate, then find them

again over that field before continuing. I can imagine a nice Sage

function to do that automatically including a field_of_definition_of

points method on such a variety.

John

> --

> You received this message because you are subscribed to the Google Groups "sage-devel" group.

> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/1620f0b7f268a16e477d74f7295d9b1340285604.camel%40orlitzky.com.

use Magma's AlgebraicallyClosedField a lot for this (for some project

which I did before Sage existed). The idea is to find the points over

Qbar, then construct the absolute field they generate, then find them

again over that field before continuing. I can imagine a nice Sage

function to do that automatically including a field_of_definition_of

points method on such a variety.

John

> You received this message because you are subscribed to the Google Groups "sage-devel" group.

> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/1620f0b7f268a16e477d74f7295d9b1340285604.camel%40orlitzky.com.

Apr 26, 2021, 11:38:53 AMApr 26

to sage-devel

On Sunday, April 25, 2021 at 4:26:18 PM UTC+2 vdelecroix wrote:

Dear Fredrik,

One technical question: I thought that your ca_t implementation

used multivariate polynomials. This is what Magma does but not

what sage does. The latter uses expression trees and take union

fields anytime there is an exactification to be done. Is it a

misconception of mine?

For "real world" problems

1) algebraic solutions of zero dimensional varieties

defined over QQ

sage: R.<x,y,z> = QQ[]

sage: A = x^4 + y^4 + z^3 - 1

sage: B = x^2*z + 2*y^2*x - x*y*z + 5

sage: C = y*z^2 - 2*(x+y+z) + 2

sage: %time solutions = R.ideal([A,B,C]).variety(QQbar)

CPU times: user 24.3 s, sys: 23.3 ms, total: 24.3 s

Wall time: 24.4 s

Though here is the bottleneck is Groebner basis computation

of multivariate polynomials and not really QQbar (21s over the 24s).

I managed to isolate the final QQbar part of this computation in a way that I could benchmark outside of Sage. Results:

sage qqbar: 500 ms

calcium ca: 80 ms

calcium ca: 80 ms

calcium qqbar: ~forever (hours)

What takes forever here is evaluating polynomials with coefficients in Q at QQbar arguments, using a generic algorithm (Horner's rule or adding powers). This is slow in the minpoly representation, because you get huge temporary results which each need to be factored. You could throw in some specialized polynomial evaluation code that generates a number field element and then computes the absolute minimal polynomial, but it'd essentially be reinventing the number field arithmetic which Sage's QQbar and Calcium's ca do automatically.

I'd be interested in some way to improve arithmetic in the minpoly representation when the operands lie in the same field (e.g. by quickly recognizing such operands using LLL). My earlier attempts to implement such hacks were unsuccessful (only slowed things down on examples I tried), but perhaps it can make difference here if I get it right.

Fredrik

Apr 26, 2021, 7:10:25 PMApr 26

to sage-devel

Hi Fredrik,

for QQbar/AA arithmetic which originate in actual

math problems, though not in Sage.

Pierre-Vincent Koseleff, Fabrice Rouillier, Cuong Tran

On the Sign of a Trigonometric Expression

ISSAC’15, July 6–9, 2015, Bath, United Kingdom.

http://dx.doi.org/10.1145/2755996.2756664.

Koseleff, P.-V., Pecker, D., and Rouillier, F.

Computing Chebyshev knot diagrams.

http://arxiv.org/abs/1001.5192

Pierre-Vincent Koseleff, Fabrice Rouillier, Cuong Tran

On the Sign of a Trigonometric Expression

ISSAC’15, July 6–9, 2015, Bath, United Kingdom.

http://dx.doi.org/10.1145/2755996.2756664.

Koseleff, P.-V., Pecker, D., and Rouillier, F.

Computing Chebyshev knot diagrams.

http://arxiv.org/abs/1001.5192

Koseleff, P.-V., Pecker, D., and Rouillier, F.

The first rational Chebyshev knots.

Journal of Symbolic Computation 45 (2010), 1341–1358.

See Pierre-Vincent Koseleff's publications page:

There are some interesting algebraic numbers in Theorem 3 of

Shigeki Akiyama

Minimum polyhedron with n vertices

which maybe can also provide benchmark problems.

Happy benchmarking! --Samuel

Reply all

Reply to author

Forward

0 new messages

Search

Clear search

Close search

Google apps

Main menu