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 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734 735 736 737 738 739 740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 825 826 827 828 829 830 831 832 833 834 835 836 837 838 839 840 841 842 843 844 845 846 847 848 849 850 851 852 853 854 855 856 857 858 859 860 861 862 863 864 865 866 867 868 869 870 871 872 873 874 875 876 877 878 879 880 881 882 883 884 885 886 887 888 889 890 891 892 893 894 895 896 897 898 899 900 901 902 903 904 905 906 907 908 909 910 911 912 913 914 915 916 917 918 919 920 921 922 923 924 925 926 927 928 929 930 931 932 933 934 935 936 937 938 939 940 941 942 943 944 945 946 947 948 949 950 951 952 953 954 955 956 957 958 959 960 961 962 963 964 965 966 967 968 969 970 971 972 973 974 975 976 977 978 979 980 981 982 983 984 985 986 987 988 989 990 991 992 993 994 995 996 997 998 999 1000 1001 1002 1003 1004 1005 1006 1007 1008 1009 1010 1011 1012 1013 1014 1015 1016 1017 1018 1019 1020 1021 1022 1023 1024 1025 1026 1027 1028 1029 1030 1031 1032 1033 1034 1035 1036 1037 1038 1039 1040 1041 1042 1043 1044 1045 1046 1047 1048 1049 1050 1051 1052 1053
|
.. _qqbar:
**qqbar.h** -- algebraic numbers represented by minimal polynomials
===============================================================================
A :type:`qqbar_t` represents a real or complex algebraic number
(an element of `\overline{\mathbb{Q}}`) by its unique reduced
minimal polynomial in `\mathbb{Z}[x]` and an isolating complex interval.
The precision of isolating intervals is maintained automatically to
ensure that all operations on :type:`qqbar_t` instances are exact.
This representation is useful for working with
individual algebraic numbers of moderate degree (up to 100, say).
Arithmetic in this representation is expensive: an arithmetic operation
on numbers of degrees *m* and *n* involves computing and then factoring an
annihilating polynomial of degree *mn* and potentially also performing
numerical root-finding. For doing repeated arithmetic, it is generally
more efficient to work with the :type:`ca_t` type in a fixed
number field.
The :type:`qqbar_t` type is used internally by the :type:`ca_t` type
to represent the embedding of number fields in `\mathbb{R}` or
`\mathbb{C}` and to decide predicates for algebraic numbers.
Types and macros
-------------------------------------------------------------------------------
.. type:: qqbar_struct
.. type:: qqbar_t
A *qqbar_struct* consists of an *fmpz_poly_struct* and an *acb_struct*.
A *qqbar_t* is defined as an array of length one of type
*qqbar_struct*, permitting a *qqbar_t* to be passed by
reference.
.. type:: qqbar_ptr
Alias for ``qqbar_struct *``, used for *qqbar* vectors.
.. type:: qqbar_srcptr
Alias for ``const qqbar_struct *``, used for *qqbar* vectors
when passed as readonly input to functions.
.. macro:: QQBAR_POLY(x)
Macro returning a pointer to the minimal polynomial of *x* which can be used as an *fmpz_poly_t*.
.. macro:: QQBAR_COEFFS(x)
Macro returning a pointer to the array of *fmpz* coefficients of the
minimal polynomial of *x*.
.. macro:: QQBAR_ENCLOSURE(x)
Macro returning a pointer to the enclosure of *x* which can be used as an *acb_t*.
Memory management
-------------------------------------------------------------------------------
.. function:: void qqbar_init(qqbar_t res)
Initializes the variable *res* for use, and sets its value to zero.
.. function:: void qqbar_clear(qqbar_t res)
Clears the variable *res*, freeing or recycling its allocated memory.
.. function:: qqbar_ptr _qqbar_vec_init(slong len)
Returns a pointer to an array of *len* initialized *qqbar_struct*:s.
.. function:: void _qqbar_vec_clear(qqbar_ptr vec, slong len)
Clears all *len* entries in the vector *vec* and frees the
vector itself.
Assignment
-------------------------------------------------------------------------------
.. function:: void qqbar_swap(qqbar_t x, qqbar_t y)
Swaps the values of *x* and *y* efficiently.
.. function:: void qqbar_set(qqbar_t res, const qqbar_t x)
void qqbar_set_si(qqbar_t res, slong x)
void qqbar_set_ui(qqbar_t res, ulong x)
void qqbar_set_fmpz(qqbar_t res, const fmpz_t x)
void qqbar_set_fmpq(qqbar_t res, const fmpq_t x)
Sets *res* to the value *x*.
.. function:: void qqbar_set_re_im(qqbar_t res, const qqbar_t x, const qqbar_t y)
Sets *res* to the value `x + yi`.
.. function:: int qqbar_set_d(qqbar_t res, double x)
int qqbar_set_re_im_d(qqbar_t res, double x, double y)
Sets *res* to the value *x* or `x + yi` respectively. These functions
performs error handling: if *x* and *y* are finite, the conversion succeeds
and the return flag is 1. If *x* or *y* is non-finite (infinity or NaN),
the conversion fails and the return flag is 0.
Properties
-------------------------------------------------------------------------------
.. function:: slong qqbar_degree(const qqbar_t x)
Returns the degree of *x*, i.e. the degree of the minimal polynomial.
.. function:: int qqbar_is_rational(const qqbar_t x)
Returns whether *x* is a rational number.
.. function:: int qqbar_is_integer(const qqbar_t x)
Returns whether *x* is an integer (an element of `\mathbb{Z}`).
.. function:: int qqbar_is_algebraic_integer(const qqbar_t x)
Returns whether *x* is an algebraic integer, i.e. whether its minimal
polynomial has leading coefficient 1.
.. function:: int qqbar_is_zero(const qqbar_t x)
int qqbar_is_one(const qqbar_t x)
int qqbar_is_neg_one(const qqbar_t x)
Returns whether *x* is the number `0`, `1`, `-1`.
.. function:: int qqbar_is_i(const qqbar_t x)
int qqbar_is_neg_i(const qqbar_t x)
Returns whether *x* is the imaginary unit `i` (respectively `-i`).
.. function:: int qqbar_is_real(const qqbar_t x)
Returns whether *x* is a real number.
.. function:: void qqbar_height(fmpz_t res, const qqbar_t x)
Sets *res* to the height of *x* (the largest absolute value of the
coefficients of the minimal polynomial of *x*).
.. function:: slong qqbar_height_bits(const qqbar_t x)
Returns the height of *x* (the largest absolute value of the
coefficients of the minimal polynomial of *x*) measured in bits.
.. function:: int qqbar_within_limits(const qqbar_t x, slong deg_limit, slong bits_limit)
Checks if *x* has degree bounded by *deg_limit* and height
bounded by *bits_limit* bits, returning 0 (false) or 1 (true).
If *deg_limit* is set to 0, the degree check is skipped,
and similarly for *bits_limit*.
.. function:: int qqbar_binop_within_limits(const qqbar_t x, const qqbar_t y, slong deg_limit, slong bits_limit)
Checks if `x + y`, `x - y`, `x \cdot y` and `x / y` certainly have
degree bounded by *deg_limit* (by multiplying the degrees for *x* and *y*
to obtain a trivial bound). For *bits_limits*, the sum of the bit heights
of *x* and *y* is checked against the bound (this is only a heuristic).
If *deg_limit* is set to 0, the degree check is skipped,
and similarly for *bits_limit*.
Conversions
-------------------------------------------------------------------------------
.. function:: void _qqbar_get_fmpq(fmpz_t num, fmpz_t den, const qqbar_t x)
Sets *num* and *den* to the numerator and denominator of *x*.
Aborts if *x* is not a rational number.
.. function:: void qqbar_get_fmpq(fmpq_t res, const qqbar_t x)
Sets *res* to *x*. Aborts if *x* is not a rational number.
.. function:: void qqbar_get_fmpz(fmpz_t res, const qqbar_t x)
Sets *res* to *x*. Aborts if *x* is not an integer.
Special values
-------------------------------------------------------------------------------
.. function:: void qqbar_zero(qqbar_t res)
Sets *res* to the number 0.
.. function:: void qqbar_one(qqbar_t res)
Sets *res* to the number 1.
.. function:: void qqbar_i(qqbar_t res)
Sets *res* to the imaginary unit `i`.
.. function:: void qqbar_phi(qqbar_t res)
Sets *res* to the golden ratio `\varphi = \tfrac{1}{2}(\sqrt{5} + 1)`.
Input and output
-------------------------------------------------------------------------------
.. function:: void qqbar_print(const qqbar_t x)
Prints *res* to standard output. The output shows the degree
and the list of coefficients
of the minimal polynomial followed by a decimal representation of
the enclosing interval. This function is mainly intended for debugging.
.. function:: void qqbar_printn(const qqbar_t x, slong n)
Prints *res* to standard output. The output shows a decimal
approximation to *n* digits.
.. function:: void qqbar_printnd(const qqbar_t x, slong n)
Prints *res* to standard output. The output shows a decimal
approximation to *n* digits, followed by the degree of the number.
For example, *print*, *printn* and *printnd* with `n = 6` give
the following output for the numbers 0, 1, `i`, `\varphi`, `\sqrt{2}-\sqrt{3} i`:
.. code::
deg 1 [0, 1] 0
deg 1 [-1, 1] 1.00000
deg 2 [1, 0, 1] 1.00000*I
deg 2 [-1, -1, 1] [1.61803398874989484820458683436563811772 +/- 6.00e-39]
deg 4 [25, 0, 2, 0, 1] [1.4142135623730950488016887242096980786 +/- 8.67e-38] + [-1.732050807568877293527446341505872367 +/- 1.10e-37]*I
0
1.00000
1.00000*I
1.61803
1.41421 - 1.73205*I
0 (deg 1)
1.00000 (deg 1)
1.00000*I (deg 2)
1.61803 (deg 2)
1.41421 - 1.73205*I (deg 4)
Random generation
-------------------------------------------------------------------------------
.. function:: void qqbar_randtest(qqbar_t res, flint_rand_t state, slong deg, slong bits)
Sets *res* to a random algebraic number with degree up to *deg* and
with height (measured in bits) up to *bits*.
.. function:: void qqbar_randtest_real(qqbar_t res, flint_rand_t state, slong deg, slong bits)
Sets *res* to a random real algebraic number with degree up to *deg* and
with height (measured in bits) up to *bits*.
.. function:: void qqbar_randtest_nonreal(qqbar_t res, flint_rand_t state, slong deg, slong bits)
Sets *res* to a random nonreal algebraic number with degree up to *deg* and
with height (measured in bits) up to *bits*. Since all algebraic numbers
of degree 1 are real, *deg* must be at least 2.
Comparisons
-------------------------------------------------------------------------------
.. function:: int qqbar_equal(const qqbar_t x, const qqbar_t y)
Returns whether *x* and *y* are equal.
.. function:: int qqbar_equal_fmpq_poly_val(const qqbar_t x, const fmpq_poly_t f, const qqbar_t y)
Returns whether *x* is equal to `f(y)`. This function is more efficient
than evaluating `f(y)` and comparing the results.
.. function:: int qqbar_cmp_re(const qqbar_t x, const qqbar_t y)
Compares the real parts of *x* and *y*, returning -1, 0 or +1.
.. function:: int qqbar_cmp_im(const qqbar_t x, const qqbar_t y)
Compares the imaginary parts of *x* and *y*, returning -1, 0 or +1.
.. function:: int qqbar_cmpabs_re(const qqbar_t x, const qqbar_t y)
Compares the absolute values of the real parts of *x* and *y*, returning -1, 0 or +1.
.. function:: int qqbar_cmpabs_im(const qqbar_t x, const qqbar_t y)
Compares the absolute values of the imaginary parts of *x* and *y*, returning -1, 0 or +1.
.. function:: int qqbar_cmpabs(const qqbar_t x, const qqbar_t y)
Compares the absolute values of *x* and *y*, returning -1, 0 or +1.
.. function:: int qqbar_cmp_root_order(const qqbar_t x, const qqbar_t y)
Compares *x* and *y* using an arbitrary but convenient ordering
defined on the complex numbers. This is useful for sorting the
roots of a polynomial in a canonical order.
We define the root order as follows: real roots come first, in
descending order. Nonreal roots are subsequently ordered first by
real part in descending order, then in ascending order by the
absolute value of the imaginary part, and then in descending
order of the sign. This implies that complex conjugate roots
are adjacent, with the root in the upper half plane first.
.. function:: ulong qqbar_hash(const qqbar_t x)
Returns a hash of *x*. As currently implemented, this function
only hashes the minimal polynomial of *x*. The user should
mix in some bits based on the numerical value if it is critical
to distinguish between conjugates of the same minimal polynomial.
This function is also likely to produce serial runs of values for
lexicographically close minimal polynomials. This is not
necessarily a problem for use in hash tables, but if it is
important that all bits in the output are random,
the user should apply an integer hash function to the output.
Complex parts
-------------------------------------------------------------------------------
.. function:: void qqbar_conj(qqbar_t res, const qqbar_t x)
Sets *res* to the complex conjugate of *x*.
.. function:: void qqbar_re(qqbar_t res, const qqbar_t x)
Sets *res* to the real part of *x*.
.. function:: void qqbar_im(qqbar_t res, const qqbar_t x)
Sets *res* to the imaginary part of *x*.
.. function:: void qqbar_re_im(qqbar_t res1, qqbar_t res2, const qqbar_t x)
Sets *res1* to the real part of *x* and *res2* to the imaginary part of *x*.
.. function:: void qqbar_abs(qqbar_t res, const qqbar_t x)
Sets *res* to the absolute value of *x*:
.. function:: void qqbar_abs2(qqbar_t res, const qqbar_t x)
Sets *res* to the square of the absolute value of *x*.
.. function:: void qqbar_sgn(qqbar_t res, const qqbar_t x)
Sets *res* to the complex sign of *x*, defined as 0 if *x* is zero
and as `x / |x|` otherwise.
.. function:: int qqbar_sgn_re(const qqbar_t x)
Returns the sign of the real part of *x* (-1, 0 or +1).
.. function:: int qqbar_sgn_im(const qqbar_t x)
Returns the sign of the imaginary part of *x* (-1, 0 or +1).
.. function:: int qqbar_csgn(const qqbar_t x)
Returns the extension of the real sign function taking the
value 1 for *x* strictly in the right half plane, -1 for *x* strictly
in the left half plane, and the sign of the imaginary part when *x* is on
the imaginary axis. Equivalently, `\operatorname{csgn}(x) = x / \sqrt{x^2}`
except that the value is 0 when *x* is zero.
Integer parts
-------------------------------------------------------------------------------
.. function:: void qqbar_floor(fmpz_t res, const qqbar_t x)
Sets *res* to the floor function of *x*. If *x* is not real, the
value is defined as the floor function of the real part of *x*.
.. function:: void qqbar_ceil(fmpz_t res, const qqbar_t x)
Sets *res* to the ceiling function of *x*. If *x* is not real, the
value is defined as the ceiling function of the real part of *x*.
Arithmetic
-------------------------------------------------------------------------------
.. function:: void qqbar_neg(qqbar_t res, const qqbar_t x)
Sets *res* to the negation of *x*.
.. function:: void qqbar_add(qqbar_t res, const qqbar_t x, const qqbar_t y)
void qqbar_add_fmpq(qqbar_t res, const qqbar_t x, const fmpq_t y)
void qqbar_add_fmpz(qqbar_t res, const qqbar_t x, const fmpz_t y)
void qqbar_add_ui(qqbar_t res, const qqbar_t x, ulong y)
void qqbar_add_si(qqbar_t res, const qqbar_t x, slong y)
Sets *res* to the sum of *x* and *y*.
.. function:: void qqbar_sub(qqbar_t res, const qqbar_t x, const qqbar_t y)
void qqbar_sub_fmpq(qqbar_t res, const qqbar_t x, const fmpq_t y)
void qqbar_sub_fmpz(qqbar_t res, const qqbar_t x, const fmpz_t y)
void qqbar_sub_ui(qqbar_t res, const qqbar_t x, ulong y)
void qqbar_sub_si(qqbar_t res, const qqbar_t x, slong y)
void qqbar_fmpq_sub(qqbar_t res, const fmpq_t x, const qqbar_t y)
void qqbar_fmpz_sub(qqbar_t res, const fmpz_t x, const qqbar_t y)
void qqbar_ui_sub(qqbar_t res, ulong x, const qqbar_t y)
void qqbar_si_sub(qqbar_t res, slong x, const qqbar_t y)
Sets *res* to the difference of *x* and *y*.
.. function:: void qqbar_mul(qqbar_t res, const qqbar_t x, const qqbar_t y)
void qqbar_mul_fmpq(qqbar_t res, const qqbar_t x, const fmpq_t y)
void qqbar_mul_fmpz(qqbar_t res, const qqbar_t x, const fmpz_t y)
void qqbar_mul_ui(qqbar_t res, const qqbar_t x, ulong y)
void qqbar_mul_si(qqbar_t res, const qqbar_t x, slong y)
Sets *res* to the product of *x* and *y*.
.. function:: void qqbar_mul_2exp_si(qqbar_t res, const qqbar_t x, slong e)
Sets *res* to *x* multiplied by `2^e`.
.. function:: void qqbar_sqr(qqbar_t res, const qqbar_t x)
Sets *res* to the square of *x*.
.. function:: void qqbar_inv(qqbar_t res, const qqbar_t x)
Sets *res* to the multiplicative inverse of *y*.
Division by zero calls *flint_abort*.
.. function:: void qqbar_div(qqbar_t res, const qqbar_t x, const qqbar_t y)
void qqbar_div_fmpq(qqbar_t res, const qqbar_t x, const fmpq_t y)
void qqbar_div_fmpz(qqbar_t res, const qqbar_t x, const fmpz_t y)
void qqbar_div_ui(qqbar_t res, const qqbar_t x, ulong y)
void qqbar_div_si(qqbar_t res, const qqbar_t x, slong y)
void qqbar_fmpq_div(qqbar_t res, const fmpq_t x, const qqbar_t y)
void qqbar_fmpz_div(qqbar_t res, const fmpz_t x, const qqbar_t y)
void qqbar_ui_div(qqbar_t res, ulong x, const qqbar_t y)
void qqbar_si_div(qqbar_t res, slong x, const qqbar_t y)
Sets *res* to the quotient of *x* and *y*.
Division by zero calls *flint_abort*.
.. function:: void qqbar_scalar_op(qqbar_t res, const qqbar_t x, const fmpz_t a, const fmpz_t b, const fmpz_t c)
Sets *res* to the rational affine transformation `(ax+b)/c`, performed as
a single operation. There are no restrictions on *a*, *b* and *c*
except that *c* must be nonzero. Division by zero calls *flint_abort*.
Powers and roots
-------------------------------------------------------------------------------
.. function:: void qqbar_sqrt(qqbar_t res, const qqbar_t x)
void qqbar_sqrt_ui(qqbar_t res, ulong x)
Sets *res* to the principal square root of *x*.
.. function:: void qqbar_rsqrt(qqbar_t res, const qqbar_t x)
Sets *res* to the reciprocal of the principal square root of *x*.
Division by zero calls *flint_abort*.
.. function:: void qqbar_pow_ui(qqbar_t res, const qqbar_t x, ulong n)
void qqbar_pow_si(qqbar_t res, const qqbar_t x, slong n)
void qqbar_pow_fmpz(qqbar_t res, const qqbar_t x, const fmpz_t n)
void qqbar_pow_fmpq(qqbar_t res, const qqbar_t x, const fmpq_t n)
Sets *res* to *x* raised to the *n*-th power.
Raising zero to a negative power aborts.
.. function:: void qqbar_root_ui(qqbar_t res, const qqbar_t x, ulong n)
void qqbar_fmpq_root_ui(qqbar_t res, const fmpq_t x, ulong n)
Sets *res* to the principal *n*-th root of *x*. The order *n*
must be positive.
.. function:: void qqbar_fmpq_pow_si_ui(qqbar_t res, const fmpq_t x, slong m, ulong n)
Sets *res* to the principal branch of `x^{m/n}`. The order *n*
must be positive. Division by zero calls *flint_abort*.
.. function:: int qqbar_pow(qqbar_t res, const qqbar_t x, const qqbar_t y)
General exponentiation: if `x^y` is an algebraic number, sets *res*
to this value and returns 1. If `x^y` is transcendental or
undefined, returns 0. Note that this function returns 0 instead of
aborting on division zero.
Numerical enclosures
-------------------------------------------------------------------------------
The following functions guarantee a polished output in which both the real
and imaginary parts are accurate to *prec* bits and exact when exactly
representable (that is, when a real or imaginary part is a sufficiently
small dyadic number).
In some cases, the computations needed to polish the output may be
expensive. When polish is unnecessary, :func:`qqbar_enclosure_raw`
may be used instead. Alternatively, :func:`qqbar_cache_enclosure`
can be used to avoid recomputations.
.. function:: void qqbar_get_acb(acb_t res, const qqbar_t x, slong prec)
Sets *res* to an enclosure of *x* rounded to *prec* bits.
.. function:: void qqbar_get_arb(arb_t res, const qqbar_t x, slong prec)
Sets *res* to an enclosure of *x* rounded to *prec* bits, assuming that
*x* is a real number. If *x* is not real, *res* is set to
`[\operatorname{NaN} \pm \infty]`.
.. function:: void qqbar_get_arb_re(arb_t res, const qqbar_t x, slong prec)
Sets *res* to an enclosure of the real part of *x* rounded to *prec* bits.
.. function:: void qqbar_get_arb_im(arb_t res, const qqbar_t x, slong prec)
Sets *res* to an enclosure of the imaginary part of *x* rounded to *prec* bits.
.. function:: void qqbar_cache_enclosure(qqbar_t res, slong prec)
Polishes the internal enclosure of *res* to at least *prec* bits
of precision in-place. Normally, *qqbar* operations that need
high-precision enclosures compute them on the fly without caching the results;
if *res* will be used as an invariant operand for many operations,
calling this function as a precomputation step can improve performance.
Numerator and denominator
-------------------------------------------------------------------------------
.. function:: void qqbar_denominator(fmpz_t res, const qqbar_t y)
Sets *res* to the denominator of *y*, i.e. the leading coefficient
of the minimal polynomial of *y*.
.. function:: void qqbar_numerator(qqbar_t res, const qqbar_t y)
Sets *res* to the numerator of *y*, i.e. *y* multiplied by
its denominator.
Conjugates
-------------------------------------------------------------------------------
.. function:: void qqbar_conjugates(qqbar_ptr res, const qqbar_t x)
Sets the entries of the vector *res* to the *d* algebraic conjugates of
*x*, including *x* itself, where *d* is the degree of *x*. The output
is sorted in a canonical order (as defined by :func:`qqbar_cmp_root_order`).
Polynomial evaluation
-------------------------------------------------------------------------------
.. function:: void _qqbar_evaluate_fmpq_poly(qqbar_t res, const fmpz * poly, const fmpz_t den, slong len, const qqbar_t x)
void qqbar_evaluate_fmpq_poly(qqbar_t res, const fmpq_poly_t poly, const qqbar_t x)
void _qqbar_evaluate_fmpz_poly(qqbar_t res, const fmpz * poly, slong len, const qqbar_t x)
void qqbar_evaluate_fmpz_poly(qqbar_t res, const fmpz_poly_t poly, const qqbar_t x)
Sets *res* to the value of the given polynomial *poly* evaluated at
the algebraic number *x*. These methods detect simple special cases and
automatically reduce *poly* if its degree is greater or equal
to that of the minimal polynomial of *x*. In the generic case, evaluation
is done by computing minimal polynomials of representation matrices.
.. function:: int qqbar_evaluate_fmpz_mpoly_iter(qqbar_t res, const fmpz_mpoly_t poly, qqbar_srcptr x, slong deg_limit, slong bits_limit, const fmpz_mpoly_ctx_t ctx)
int qqbar_evaluate_fmpz_mpoly_horner(qqbar_t res, const fmpz_mpoly_t poly, qqbar_srcptr x, slong deg_limit, slong bits_limit, const fmpz_mpoly_ctx_t ctx)
int qqbar_evaluate_fmpz_mpoly(qqbar_t res, const fmpz_mpoly_t poly, qqbar_srcptr x, slong deg_limit, slong bits_limit, const fmpz_mpoly_ctx_t ctx)
Sets *res* to the value of *poly* evaluated at the algebraic numbers
given in the vector *x*. The number of variables is defined by
the context object *ctx*.
The parameters *deg_limit* and *bits_limit*
define evaluation limits: if any temporary result exceeds these limits
(not necessarily the final value, in case of cancellation), the
evaluation is aborted and 0 (failure) is returned. If evaluation
succeeds, 1 is returned.
The *iter* version iterates over all terms in succession and computes
the powers that appear. The *horner* version uses a multivariate
implementation of the Horner scheme. The default algorithm currently
uses the Horner scheme.
Polynomial roots
-------------------------------------------------------------------------------
.. function:: void qqbar_roots_fmpz_poly(qqbar_ptr res, const fmpz_poly_t poly, int flags)
void qqbar_roots_fmpq_poly(qqbar_ptr res, const fmpq_poly_t poly, int flags)
Sets the entries of the vector *res* to the *d* roots of the polynomial
*poly*. Roots with multiplicity appear with repetition in the
output array. By default, the roots will be sorted in a
convenient canonical order (as defined by :func:`qqbar_cmp_root_order`).
Instances of a repeated root always appear consecutively.
The following *flags* are supported:
- QQBAR_ROOTS_IRREDUCIBLE - if set, *poly* is assumed to be
irreducible (it may still have constant content), and no polynomial
factorization is performed internally.
- QQBAR_ROOTS_UNSORTED - if set, the roots will not be guaranteed
to be sorted (except for repeated roots being listed
consecutively).
.. function:: void qqbar_eigenvalues_fmpz_mat(qqbar_ptr res, const fmpz_mat_t mat, int flags)
void qqbar_eigenvalues_fmpq_mat(qqbar_ptr res, const fmpq_mat_t mat, int flags)
Sets the entries of the vector *res* to the eigenvalues of the
square matrix *mat*. These functions compute the characteristic polynomial
of *mat* and then call :func:`qqbar_roots_fmpz_poly` with the same
flags.
.. function:: int _qqbar_roots_poly_squarefree(qqbar_ptr roots, qqbar_srcptr coeffs, slong len, slong deg_limit, slong bits_limit)
Writes to the vector *roots* the `d` roots of the polynomial
with algebraic number coefficients *coeffs* of length *len* (`d = len - 1`).
Given the polynomial `f = a_0 + \ldots + a_d x^d` with coefficients in
`\overline{\mathbb{Q}}`, we construct an annihilating polynomial with
coefficients in `\mathbb{Q}` as
`g = \prod (\tilde a_0 + \ldots + \tilde a_d x^d)`
where the product is taken over all
combinations of algebraic conjugates `\tilde a_k` of the input
coefficients.
The polynomial `g` is subsequently factored to find candidate
roots.
The leading coefficient `a_d` must be nonzero and the polynomial `f`
polynomial must be squarefree.
To compute roots of a general polynomial which may have repeated roots,
it is necessary to perform a squarefree factorization before calling
this function.
An option is to call :func:`gr_poly_roots` with a ``qqbar`` context object,
which wraps this function
and takes care of the initial squarefree factorization.
Since the product `g` can explode in size very quickly, the *deg_limit*
and *bits_limit* parameters allow bounding the degree and working precision.
The function returns 1 on success and 0 on failure indicating that
such a limit has been exceeded.
Setting nonpositive values for the limits removes the restrictions;
however, the function can still fail and return 0 in that case if `g`
exceeds machine size.
Note: to compute algebraic number roots of polynomials of various other
types, use :func:`gr_poly_roots_other`.
Roots of unity and trigonometric functions
-------------------------------------------------------------------------------
The following functions use word-size integers *p* and *q*
instead of *fmpq_t* instances to express rational numbers.
This is to emphasize that
the computations are feasible only with small *q* in this representation
of algebraic numbers since the
associated minimal polynomials have degree `O(q)`.
The input *p* and *q* do not need to be reduced *a priori*,
but should not be close to the word boundaries (they may be added
and subtracted internally).
.. function:: void qqbar_root_of_unity(qqbar_t res, slong p, ulong q)
Sets *res* to the root of unity `e^{2 \pi i p / q}`.
.. function:: int qqbar_is_root_of_unity(slong * p, ulong * q, const qqbar_t x)
If *x* is not a root of unity, returns 0.
If *x* is a root of unity, returns 1.
If *p* and *q* are not *NULL* and *x* is a root of unity,
this also sets *p* and *q* to the minimal integers with `0 \le p < q`
such that `x = e^{2 \pi i p / q}`.
.. function:: void qqbar_exp_pi_i(qqbar_t res, slong p, ulong q)
Sets *res* to the root of unity `e^{\pi i p / q}`.
.. function:: void qqbar_cos_pi(qqbar_t res, slong p, ulong q)
void qqbar_sin_pi(qqbar_t res, slong p, ulong q)
int qqbar_tan_pi(qqbar_t res, slong p, ulong q)
int qqbar_cot_pi(qqbar_t res, slong p, ulong q)
int qqbar_sec_pi(qqbar_t res, slong p, ulong q)
int qqbar_csc_pi(qqbar_t res, slong p, ulong q)
Sets *res* to the trigonometric function `\cos(\pi x)`,
`\sin(\pi x)`, etc., with `x = \tfrac{p}{q}`.
The functions tan, cot, sec and csc return the flag 1 if the value exists,
and return 0 if the evaluation point is a pole of the function.
.. function:: int qqbar_log_pi_i(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{log}(x) / (\pi i)` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `-1 < y \le 1` and returns 1.
If *y* is not algebraic, returns 0.
.. function:: int qqbar_atan_pi(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{atan}(x) / \pi` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `|y| < \tfrac{1}{2}` and returns 1.
If *y* is not algebraic, returns 0.
.. function:: int qqbar_asin_pi(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{asin}(x) / \pi` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `|y| \le \tfrac{1}{2}` and returns 1.
If *y* is not algebraic, returns 0.
.. function:: int qqbar_acos_pi(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{acos}(x) / \pi` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `0 \le y \le 1` and returns 1.
If *y* is not algebraic, returns 0.
.. function:: int qqbar_acot_pi(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{acot}(x) / \pi` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `-\tfrac{1}{2} < y \le \tfrac{1}{2}` and returns 1.
If *y* is not algebraic, returns 0.
.. function:: int qqbar_asec_pi(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{asec}(x) / \pi` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `0 \le y \le 1` and returns 1.
If *y* is not algebraic, returns 0.
.. function:: int qqbar_acsc_pi(slong * p, ulong * q, const qqbar_t x)
If `y = \operatorname{acsc}(x) / \pi` is algebraic, and hence
necessarily rational, sets `y = p / q` to the reduced such
fraction with `-\tfrac{1}{2} \le y \le \tfrac{1}{2}` and returns 1.
If *y* is not algebraic, returns 0.
Guessing and simplification
-------------------------------------------------------------------------------
.. function:: int qqbar_guess(qqbar_t res, const acb_t z, slong max_deg, slong max_bits, int flags, slong prec)
Attempts to find an algebraic number *res* of degree at most *max_deg* and
height at most *max_bits* bits matching the numerical enclosure *z*.
The return flag indicates success.
This is only a heuristic method, and the return flag neither implies a
rigorous proof that *res* is the correct result, nor a rigorous proof
that no suitable algebraic number with the given *max_deg* and *max_bits*
exists. (Proof of nonexistence could in principle be computed,
but this is not yet implemented.)
The working precision *prec* should normally be the same as the precision
used to compute *z*. It does not make much sense to run this algorithm
with precision smaller than O(*max_deg* ยท *max_bits*).
This function does a single iteration at the target *max_deg*, *max_bits*,
and *prec*. For best performance, one should invoke this function
repeatedly with successively larger parameters when the size of the
intended solution is unknown or may be much smaller than a worst-case bound.
.. function:: int qqbar_express_in_field(fmpq_poly_t res, const qqbar_t alpha, const qqbar_t x, slong max_bits, int flags, slong prec)
Attempts to express *x* in the number field generated by *alpha*, returning
success (0 or 1). On success, *res* is set to a polynomial *f* of degree
less than the degree of *alpha* and with height (counting both the numerator
and the denominator, when the coefficients of *g* are
put on a common denominator) bounded by *max_bits* bits, such that
`f(\alpha) = x`.
(Exception: the *max_bits* parameter is currently ignored if *x* is
rational, in which case *res* is just set to the value of *x*.)
This function looks for a linear relation heuristically using a working
precision of *prec* bits. If *x* is expressible in terms of *alpha*,
then this function is guaranteed to succeed when *prec* is taken
large enough. The identity `f(\alpha) = x` is checked
rigorously, i.e. a return value of 1 implies a proof of correctness.
In principle, choosing a sufficiently large *prec* can be used to
prove that *x* does not lie in the field generated by *alpha*,
but the present implementation does not support doing so automatically.
This function does a single iteration at the target *max_bits* and
and *prec*. For best performance, one should invoke this function
repeatedly with successively larger parameters when the size of the
intended solution is unknown or may be much smaller than a worst-case bound.
Symbolic expressions and conversion to radicals
-------------------------------------------------------------------------------
.. function:: void qqbar_get_quadratic(fmpz_t a, fmpz_t b, fmpz_t c, fmpz_t q, const qqbar_t x, int factoring)
Assuming that *x* has degree 1 or 2, computes integers *a*, *b*, *c*
and *q* such that
.. math::
x = \frac{a + b \sqrt{c}}{q}
and such that *c* is not a perfect square, *q* is positive, and
*q* has no content in common with both *a* and *b*. In other words,
this determines a quadratic field `\mathbb{Q}(\sqrt{c})` containing
*x*, and then finds the canonical reduced coefficients *a*, *b* and
*q* expressing *x* in this field.
For convenience, this function supports rational *x*,
for which *b* and *c* will both be set to zero.
The following remarks apply to irrationals.
The radicand *c* will not be a perfect square, but will not
automatically be squarefree since this would require factoring the
discriminant. As a special case, *c* will be set to `-1` if *x*
is a Gaussian rational number. Otherwise, behavior is controlled
by the *factoring* parameter.
* If *factoring* is 0, no factorization is performed apart from
removing powers of two.
* If *factoring* is 1, a complete factorization is performed (*c*
will be minimal). This can be very expensive if the discriminant
is large.
* If *factoring* is 2, a smooth factorization is performed to remove
small factors from *c*. This is a tradeoff that provides pretty
output in most cases while avoiding extreme worst-case slowdown.
The smooth factorization guarantees finding all small factors
(up to some trial division limit determined internally by FLINT),
but large factors are only found heuristically.
.. function:: int qqbar_set_fexpr(qqbar_t res, const fexpr_t expr)
Sets *res* to the algebraic number represented by the symbolic
expression *expr*, returning 1 on success and 0 on failure.
This function performs a "static" evaluation using *qqbar* arithmetic,
supporting only closed-form expressions with explicitly algebraic
subexpressions. It can be used to recover values generated by
:func:`qqbar_get_expr_formula` and variants.
For evaluating more complex expressions involving other
types of values or requiring symbolic simplifications,
the user should preprocess *expr* so that it is in a form
which can be parsed by :func:`qqbar_set_fexpr`.
The following expressions are supported:
* Integer constants
* Arithmetic operations with algebraic operands
* Square roots of algebraic numbers
* Powers with algebraic base and exponent an explicit rational number
* NumberI, GoldenRatio, RootOfUnity
* Floor, Ceil, Abs, Sign, Csgn, Conjugate, Re, Im, Max, Min
* Trigonometric functions with argument an explicit rational number times Pi
* Exponentials with argument an explicit rational number times Pi * NumberI
* The Decimal() constructor
* AlgebraicNumberSerialized() (assuming valid data, which is not checked)
* PolynomialRootIndexed()
* PolynomialRootNearest()
Examples of formulas that are not supported, despite the value being
an algebraic number:
* ``Pi - Pi`` (general transcendental simplifications are not performed)
* ``1 / Infinity`` (only numbers are handled)
* ``Sum(n, For(n, 1, 10))`` (only static evaluation is performed)
.. function:: void qqbar_get_fexpr_repr(fexpr_t res, const qqbar_t x)
Sets *res* to a symbolic expression reflecting the exact internal
representation of *x*. The output will have the form
``AlgebraicNumberSerialized(List(coeffs), enclosure)``.
The output can be converted back to a ``qqbar_t`` value using
:func:`qqbar_set_fexpr`. This is the recommended format for
serializing algebraic numbers as it requires minimal computation,
but it has the disadvantage of not being human-readable.
.. function:: void qqbar_get_fexpr_root_nearest(fexpr_t res, const qqbar_t x)
Sets *res* to a symbolic expression unambiguously describing *x*
in the form ``PolynomialRootNearest(List(coeffs), point)``
where *point* is an approximation of *x* guaranteed to be closer
to *x* than any conjugate root.
The output can be converted back to a ``qqbar_t`` value using
:func:`qqbar_set_fexpr`. This is a useful format for human-readable
presentation, but serialization and deserialization can be expensive.
.. function:: void qqbar_get_fexpr_root_indexed(fexpr_t res, const qqbar_t x)
Sets *res* to a symbolic expression unambiguously describing *x*
in the form ``PolynomialRootIndexed(List(coeffs), index)``
where *index* is the index of *x* among its conjugate roots
in the builtin root sort order.
The output can be converted back to a ``qqbar_t`` value using
:func:`qqbar_set_fexpr`. This is a useful format for human-readable
presentation when the numerical value is important, but serialization
and deserialization can be expensive.
.. function:: int qqbar_get_fexpr_formula(fexpr_t res, const qqbar_t x, ulong flags)
Attempts to express the algebraic number *x* as a closed-form expression
using arithmetic operations, radicals, and possibly exponentials
or trigonometric functions, but without using ``PolynomialRootNearest``
or ``PolynomialRootIndexed``.
Returns 0 on failure and 1 on success.
The *flags* parameter toggles different methods for generating formulas.
It can be set to any combination of the following. If *flags* is 0,
only rational numbers will be handled.
.. macro:: QQBAR_FORMULA_ALL
Toggles all methods (potentially expensive).
.. macro:: QQBAR_FORMULA_GAUSSIANS
Detect Gaussian rational numbers `a + bi`.
.. macro:: QQBAR_FORMULA_QUADRATICS
Solve quadratics in the form `a + b \sqrt{d}`.
.. macro:: QQBAR_FORMULA_CYCLOTOMICS
Detect elements of cyclotomic fields. This works by trying plausible
cyclotomic fields (based on the degree of the input), using LLL
to find candidate number field elements, and certifying candidates
through an exact computation. Detection is heuristic and
is not guaranteed to find all cyclotomic numbers.
.. macro:: QQBAR_FORMULA_CUBICS
QQBAR_FORMULA_QUARTICS
QQBAR_FORMULA_QUINTICS
Solve polynomials of degree 3, 4 and (where applicable) 5 using
cubic, quartic and quintic formulas (not yet implemented).
.. macro:: QQBAR_FORMULA_DEPRESSION
Use depression to try to generate simpler numbers.
.. macro:: QQBAR_FORMULA_DEFLATION
Use deflation to try to generate simpler numbers.
This allows handling number of the form `a^{1/n}` where *a* can
be represented in closed form.
.. macro:: QQBAR_FORMULA_SEPARATION
Try separating real and imaginary parts or sign and magnitude of
complex numbers. This allows handling numbers of the form `a + bi`
or `m \cdot s` (with `m > 0, |s| = 1`) where *a* and *b* or *m* and *s* can be
represented in closed form. This is only attempted as a fallback after
other methods fail: if an explicit Cartesian or magnitude-sign
represented is desired, the user should manually separate the number
into complex parts before calling :func:`qqbar_get_fexpr_formula`.
.. macro:: QQBAR_FORMULA_EXP_FORM
QQBAR_FORMULA_TRIG_FORM
QQBAR_FORMULA_RADICAL_FORM
QQBAR_FORMULA_AUTO_FORM
Select output form for cyclotomic numbers. The *auto* form (equivalent
to no flags being set) results in radicals for numbers of low degree,
trigonometric functions for real numbers, and complex exponentials
for nonreal numbers. The other flags (not fully implemented) can be
used to force exponential form, trigonometric form, or radical form.
Internal functions
-------------------------------------------------------------------------------
.. function:: void qqbar_fmpz_poly_composed_op(fmpz_poly_t res, const fmpz_poly_t A, const fmpz_poly_t B, int op)
Given nonconstant polynomials *A* and *B*, sets *res* to a polynomial
whose roots are `a+b`, `a-b`, `ab` or `a/b` for all roots *a* of *A*
and all roots *b* of *B*. The parameter *op* selects the arithmetic
operation: 0 for addition, 1 for subtraction, 2 for multiplication
and 3 for division. If *op* is 3, *B* must not have zero as a root.
.. function:: void qqbar_binary_op(qqbar_t res, const qqbar_t x, const qqbar_t y, int op)
Performs a binary operation using a generic algorithm. This does not
check for special cases.
.. function:: int _qqbar_validate_uniqueness(acb_t res, const fmpz_poly_t poly, const acb_t z, slong max_prec)
Given *z* known to be an enclosure of at least one root of *poly*,
certifies that the enclosure contains a unique root, and in that
case sets *res* to a new (possibly improved) enclosure for the same
root, returning 1. Returns 0 if uniqueness cannot be certified.
The enclosure is validated by performing a single step with the
interval Newton method. The working precision is determined from the
accuracy of *z*, but limited by *max_prec* bits.
This method slightly inflates the enclosure *z* to improve the chances
that the interval Newton step will succeed. Uniqueness on this larger
interval implies uniqueness of the original interval, but not
existence; when existence has not been ensured a priori,
:func:`_qqbar_validate_existence_uniqueness` should be used instead.
.. function:: int _qqbar_validate_existence_uniqueness(acb_t res, const fmpz_poly_t poly, const acb_t z, slong max_prec)
Given any complex interval *z*, certifies that the enclosure contains a
unique root of *poly*, and in that case sets *res* to a new (possibly
improved) enclosure for the same root, returning 1. Returns 0 if
existence and uniqueness cannot be certified.
The enclosure is validated by performing a single step with the
interval Newton method. The working precision is determined from the
accuracy of *z*, but limited by *max_prec* bits.
.. function:: void _qqbar_enclosure_raw(acb_t res, const fmpz_poly_t poly, const acb_t z, slong prec)
void qqbar_enclosure_raw(acb_t res, const qqbar_t x, slong prec)
Sets *res* to an enclosure of *x* accurate to about *prec* bits
(the actual accuracy can be slightly lower, or higher).
This function uses repeated interval Newton steps to polish the initial
enclosure *z*, doubling the working precision each time. If any step
fails to improve the accuracy significantly, the root is recomputed
from scratch to higher precision.
If the initial enclosure is accurate enough, *res* is set to this value
without rounding and without further computation.
.. function:: int _qqbar_acb_lindep(fmpz * rel, acb_srcptr vec, slong len, int check, slong prec)
Attempts to find an integer vector *rel* giving a linear relation between
the elements of the real or complex vector *vec*, using the LLL algorithm.
The working precision is set to the minimum of *prec* and the relative
accuracy of *vec* (that is, the difference between the largest magnitude
and the largest error magnitude within *vec*). 95% of the bits within the
working precision are used for the LLL matrix, and the remaining 5% bits
are used to validate the linear relation by evaluating the linear
combination and checking that the resulting interval contains zero.
This validation does not prove the existence or nonexistence
of a linear relation, but it provides a quick heuristic way to eliminate
spurious relations.
If *check* is set, the return value indicates whether the validation
was successful; otherwise, the return value simply indicates whether
the algorithm was executed normally (failure may occur, for example,
if the input vector is non-finite).
In principle, this method can be used to produce a proof that no linear
relation exists with coefficients up to a specified bit size, but this has
not yet been implemented.
.. raw:: latex
\newpage
|