This section contains functions related to the arithmetic of quadratic number fields. More general functions, more systematic and easier to use, apply to general number fields and will be covered in the next sections. To study the arithmetic of quadratic fields, you may use either sets of functions For instance if you want integral binary quadratic forms (and classes) rather than ideals (and ideal classes). The quadratic setting also treats non-maximal orders transparently, whereas there is no support for non-maximal orders in the general setting.
We create quadratic numbers using quadgen and binary quadratic
forms using Qfb. Quadratic numbers are added, multiplied, etc. using
standard GP operators (+, *, etc.) and quadratic forms are
composed using *.
Reduction of forms modulo SL(2,ℤ) is done via qfbred or
qfbredsl2 when the unimodular transformation matrix is desired.
The functions qfbcomp and qfbpow compute products and powers
in the class group, incorporating reduction at each step. The function
qfbsolve solves the equation q(x,y) = n in integers x,y.
The functions qfbclassno, quadunit, quadregulator
respectively return the ordinary (not narrow) class number,
thefundamental unit, and the regulator. The function quadclassunit is
asymptotically faster as the discriminant tends to infinity and return the
(ordinary) class group structure and the regulator of a quadratic order all
at once, but it is conditional on the GRH. The functions quadhilbert
and quadray compute Hilbert class fields and ray class fields using
complex multiplication or Stark units.
Ordinary class number of the quadratic order of discriminant D, for "small" values of D.
* if D > 0 or flag = 1, use a O(|D|1/2)
algorithm (compute L(1,χD) with the approximate functional equation).
This is slower than quadclassunit as soon as |D| ~ 102 or
so and is not meant to be used for large D.
* if D < 0 and flag = 0 (or omitted), use a O(|D|1/4)
algorithm (Shanks's baby-step/giant-step method). It should
be faster than quadclassunit for small values of D, say
|D| < 1018.
Important warning. In the latter case, this function only
implements part of Shanks's method (which allows to speed it up
considerably). It gives unconditionally correct results for
|D| < 2.1010, but may give incorrect results for larger values
if the class
group has many cyclic factors. We thus recommend to double-check results
using the function quadclassunit, which is about 2 to 3 times slower in
the range |D| ∈ [1010, 1018], assuming GRH. We currently have no
counter-examples but they should exist: we would appreciate a bug report if
you find one.
Warning. Contrary to what its name implies, this routine does not compute the number of classes of binary primitive forms of discriminant D, which is equal to the narrow class number. The two notions are the same when D < 0 or the fundamental unit ϵ has negative norm; when D > 0 and Nϵ > 0, the number of classes of forms is twice the ordinary class number. This is a problem which we cannot fix for backward compatibility reasons. Use the following routine if you are only interested in the number of classes of forms:
? QFBclassno(D) = qfbclassno(D) * if (D > 0 && quadunitnorm(D) > 0, 2, 1) ? QFBclassno(136) %1 = 4 ? qfbclassno(136) %2 = 2 ? quadunitnorm(136) %3 = 1 ? bnfnarrow(bnfinit(x^2 - 136)).cyc %4 = [4] \\ narrow class group is cyclic ~ Z/4Z
Note that the use of bnfnarrow above is only valid
because 136 is a fundamental discriminant: that function is asymptotically
faster (and returns the group structure, not only its order) but only supports
maximal orders.
Here are a few more examples:
? qfbclassno(400000028) \\ D > 0: slow
time = 490 ms.
%1 = 1
? quadclassunit(400000028).no
time = 2 ms. \\ { much faster, assume GRH}
%2 = 1
? qfbclassno(-400000028) \\ D < 0: fast enough
time = 0 ms.
%3 = 7253
? quadclassunit(-400000028).no
time = 0 ms.
%4 = 7253
See also qfbhclassno.
The library syntax is GEN qfbclassno0(GEN D, long flag).
composition of the binary quadratic forms x and y, with reduction of the result.
? x=Qfb(2,3,-10);y=Qfb(5,3,-4); ? qfbcomp(x,y) %2 = Qfb(-2, 9, 1) ? qfbcomp(x,y)==qfbred(qfbcompraw(x,y)) %3 = 1
The library syntax is GEN qfbcomp(GEN x, GEN y).
composition of the binary quadratic forms x and y, without reduction of the result. This is useful e.g. to compute a generating element of an ideal. The result is undefined if x and y do not have the same discriminant.
? x=Qfb(2,3,-10);y=Qfb(5,3,-4);
? qfbcompraw(x,y)
%2 = Qfb(10, 3, -2)
? x=Qfb(2,3,-10);y=Qfb(1,-1,1);
? qfbcompraw(x,y)
*** at top-level: qfbcompraw(x,y)
*** ^ — — — — —
*** qfbcompraw: inconsistent qfbcompraw t_QFB , t_QFB.
The library syntax is GEN qfbcompraw(GEN x, GEN y).
Solves the equation x2 + dy2 = n in integers x and y, where
d > 0 and n is prime. Returns the empty vector [] when no solution
exists. It is also allowed to try n = 4 times a prime but the answer is
then guaranteed only if d is 3 mod 4; more precisely if d != 3 mod
4, the algorithm may fail to find a non-primitive solution.
This function is a special case of qfbsolve applied to the principal
form in the imaginary quadratic order of discriminant -4d (returning the
solution with non-negative x and y). As its name implies,
qfbcornacchia uses Cornacchia's algorithm and runs in time quasi-linear
in log n (using halfgcd); in practical ranges, qfbcornacchia
should be about twice faster than qfbsolve unless we indicate to the
latter that its second argument is prime (see below).
? qfbcornacchia(1, 113) %1 = [8, 7] ? qfbsolve(Qfb(1,0,1), 113) %2 = [8, 7] ? qfbcornacchia(1, 4*113) \\ misses the non-primitive solution 2*[8,7] %3 = [] ? qfbcornacchia(1, 4*109) \\ finds a non-primitive solution %4 = [20, 6] ? p = 122838793181521; isprime(p) %5 = 1 ? qfbcornacchia(24, p) %6 = [10547339, 694995] ? Q = Qfb(1,0,24); qfbsolve(Q,p) %7 = [10547339, 694995] ? for (i=1, 10^5, qfbsolve(Q, p)) time = 345 ms. ? for (i=1, 10^5, qfbcornacchia(24,p)) \\ faster time = 251 ms. ? for (i=1, 10^5, qfbsolve(Q, Mat([p,1]))) \\ just as fast time = 251 ms.
We used Mat([p,1]) to indicate that p1
was the integer factorization of p, i.e., that p is prime. Without it,
qfbsolve attempts to factor p and wastes a little time.
The library syntax is GEN qfbcornacchia(GEN d, GEN n).
Hurwitz class number of x, when
x is nonnegative and congruent to 0 or 3 modulo 4, and 0 for other
values. For x > 5.105, we assume the GRH, and use quadclassunit
with default parameters.
? qfbhclassno(1) \\ not 0 or 3 mod 4 %1 = 0 ? qfbhclassno(3) %2 = 1/3 ? qfbhclassno(4) %3 = 1/2 ? qfbhclassno(23) %4 = 3
The library syntax is GEN hclassno(GEN x).
composition of the primitive positive
definite binary quadratic forms x and y (type t_QFB) using the NUCOMP
and NUDUPL algorithms of Shanks, à la Atkin. L is any positive
constant, but for optimal speed, one should take L = |D/4|1/4, i.e.
sqrtnint(abs(D) >> 2,4), where D is the common discriminant of x and
y. When x and y do not have the same discriminant, the result is
undefined.
The current implementation is slower than the generic routine for small D, and becomes faster when D has about 45 bits.
The library syntax is GEN nucomp(GEN x, GEN y, GEN L).
Also available is GEN nudupl(GEN x, GEN L) when x = y.
n-th power of the primitive positive definite
binary quadratic form x using Shanks's NUCOMP and NUDUPL algorithms;
if set, L should be equal to sqrtnint(abs(D) >> 2,4), where D < 0 is
the discriminant of x.
The current implementation is slower than the generic routine for small discriminant D, and becomes faster for D ~ 245.
The library syntax is GEN nupow(GEN x, GEN n, GEN L = NULL).
n-th power of the binary quadratic form
x, computed with reduction (i.e. using qfbcomp).
The library syntax is GEN qfbpow(GEN x, GEN n).
n-th power of the binary quadratic form
x, computed without doing any reduction (i.e. using qfbcompraw).
Here n must be nonnegative and n < 231.
The library syntax is GEN qfbpowraw(GEN x, long n).
Prime binary quadratic form of discriminant
x whose first coefficient is p, where |p| is a prime number.
By abuse of notation,
p = ± 1 is also valid and returns the unit form. Returns an
error if x is not a quadratic residue mod p, or if x < 0 and p < 0.
(Negative definite t_QFB are not implemented.)
The library syntax is GEN primeform(GEN x, GEN p).
Reduces the binary quadratic form x (updating Shanks's distance function d if x = [q,d] is an extended indefinite form). If flag is 1, the function performs a single reduction step, and a complete reduction otherwise.
The arguments isd, sd, if present, supply the values of floor{sqrt{D}}, and sqrt{D} respectively, where D is the discriminant (this is not checked). If d < 0 these values are useless.
The library syntax is GEN qfbred0(GEN x, long flag, GEN isd = NULL, GEN sd = NULL).
Also available is GEN qfbred(GEN x) (flag is 0, isd
and sd are NULL)
Reduction of the (real or imaginary) binary quadratic form x, returns
[y,g] where y is reduced and g in SL(2,ℤ) is such that
g.x = y; isD, if
present, must be equal to sqrtint(D), where D > 0 is the
discriminant of x.
The action of g on x can be computed using qfeval(x,g)
? q1 = Qfb(33947,-39899,11650); ? [q2,U] = qfbredsl2(q1) %2 = [Qfb(749,2207,-1712),[-1,3;-2,5]] ? qfeval(q1,U) %3 = Qfb(749,2207,-1712)
The library syntax is GEN qfbredsl2(GEN x, GEN isD = NULL).
Solve the equation Q(x,y) = n in coprime integers x and y (primitive solutions), where Q is a binary quadratic form and n an integer, up to the action of the special orthogonal group G = SO(Q,ℤ), which is isomorphic to the group of units of positive norm of the quadratic order of discriminant D = disc Q. If D > 0, G is infinite. If D < -4, G is of order 2, if D = -3, G is of order 6 and if D = -4, G is of order 4.
Binary digits of flag mean: 1: return all solutions if set, else a single solution; return [] if a single solution is wanted (bit unset) but none exist. 2: also include imprimitive solutions.
When flag = 2 (return a single solution, possibly imprimitive), the algorithm returns a solution with minimal content; in particular, a primitive solution exists if and only if one is returned.
The integer n can also be given by its factorization matrix
fa = factor(n) or by the pair [n, fa].
? qfbsolve(Qfb(1,0,2), 603) \\ a single primitive solution
%1 = [5, 17]
? qfbsolve(Qfb(1,0,2), 603, 1) \\ all primitive solutions
%2 = [[5, 17], [-19, -11], [19, -11], [5, -17]]
? qfbsolve(Qfb(1,0,2), 603, 2) \\ a single, possibly imprimitive solution
%3 = [5, 17] \\ actually primitive
? qfbsolve(Qfb(1,0,2), 603, 3) \\ all solutions
%4 = [[5, 17], [-19, -11], [19, -11], [5, -17], [-21, 9], [-21, -9]]
? N = 2^128+1; F = factor(N);
? qfbsolve(Qfb(1,0,1),[N,F],1)
%3 = [[-16382350221535464479,8479443857936402504],
[18446744073709551616,-1],[-18446744073709551616,-1],
[16382350221535464479,8479443857936402504]]
For fixed Q, assuming the factorisation of n is given, the algorithm
runs in probabilistic polynomial time in log p, where p is the largest
prime divisor of n, through the computation of square roots of D modulo
4 p). The dependency on Q is more complicated: polynomial time in log
|D| if Q is imaginary, but exponential time if Q is real (through the
computation of a full cycle of reduced forms). In the latter case, note that
bnfisprincipal provides a solution in heuristic subexponential time
assuming the GRH.
The library syntax is GEN qfbsolve(GEN Q, GEN n, long flag).
Buchmann-McCurley's sub-exponential algorithm for computing the class group of a quadratic order of discriminant D. By default, the results are conditional on the GRH.
This function should be used instead of qfbclassno or
quadregulator
when D < -1025, D > 1010, or when the structure is wanted. It
is a special case of bnfinit, which is slower, but more robust.
The result is a vector v whose components should be accessed using member functions:
* v.no: the class number
* v.cyc: a vector giving the structure of the class group as a
product of cyclic groups;
* v.gen: a vector giving generators of those cyclic groups (as
binary quadratic forms).
* v.reg: the regulator, computed to an accuracy which is the
maximum of an internal accuracy determined by the program and the current
default (note that once the regulator is known to a small accuracy it is
trivial to compute it to very high accuracy, see the tutorial).
* v.normfu (for positive D only) return the norm of the
fundamental unit, either 1 or -1. Note that a result of -1 is
unconditional and no longer depends on the GRH.
The flag is obsolete and should be left alone. In older versions,
it supposedly computed the narrow class group when D > 0, but this did not
work at all; use the general function bnfnarrow.
Optional parameter tech is a row vector of the form [c1, c2], where c1 ≤ c2 are nonnegative real numbers which control the execution time and the stack size, see se:GRHbnf. The parameter is used as a threshold to balance the relation finding phase against the final linear algebra. Increasing the default c1 means that relations are easier to find, but more relations are needed and the linear algebra will be harder. The default value for c1 is 0 and means that it is taken equal to c2. The parameter c2 is mostly obsolete and should not be changed, but we still document it for completeness: we compute a tentative class group by generators and relations using a factorbase of prime ideals ≤ c1 (log |D|)2, then prove that ideals of norm ≤ c2 (log |D|)2 do not generate a larger group. By default an optimal c2 is chosen, so that the result is provably correct under the GRH — a result of Grenié and Molteni states that c2 = 23/6 ~ 3.83 is fine (and even c2 = 15/4 ~ 3.75 for large |D| > 2.41 E8). But it is possible to improve on this algorithmically. You may provide a smaller c2, it will be ignored (we use the provably correct one); you may provide a larger c2 than the default value, which results in longer computing times for equally correct outputs (under GRH).
The library syntax is GEN quadclassunit0(GEN D, long flag, GEN tech = NULL, long prec).
If you really need to experiment with the tech parameter,
it will be more convenient to use
GEN Buchquad(GEN D, double c1, double c2, long prec).
Discriminant of the étale algebra ℚ(sqrt{x}), where x ∈ ℚ*.
This is the same as coredisc(d) where d is the integer
squarefree part of x, so x = d f2 with f ∈ ℚ* and d ∈ ℤ.
This returns 0 for x = 0, 1 for x square and the discriminant of
the quadratic field ℚ(sqrt{x}) otherwise.
? quaddisc(7) %1 = 28 ? quaddisc(-7) %2 = -7
The library syntax is GEN quaddisc(GEN x).
'w})Creates the quadratic number ω = (a+sqrt{D})/2 where a = 0 if D = 0 mod 4, a = 1 if D = 1 mod 4, so that (1,ω) is an integral basis for the quadratic order of discriminant D. D must be an integer congruent to 0 or 1 modulo 4, which is not a square. If v is given, the variable name is used to display g else 'w' is used.
? w = quadgen(5, 'w); w^2 - w - 1 %1 = 0 ? w = quadgen(0, 'w) *** at top-level: w=quadgen(0) *** ^ — — — - *** quadgen: domain error in quadpoly: issquare(disc) = 1
The library syntax is GEN quadgen0(GEN D, long v = -1) where v is a variable number.
When v does not matter, the function
GEN quadgen(GEN D) is also available.
Relative equation defining the Hilbert class field of the quadratic field of discriminant D.
If D < 0, uses complex multiplication (Schertz's variant).
If D > 0 Stark units are used and (in rare cases) a
vector of extensions may be returned whose compositum is the requested class
field. See bnrstark for details.
The library syntax is GEN quadhilbert(GEN D, long prec).
'x})
Creates the "canonical" quadratic
polynomial (in the variable v) corresponding to the discriminant D,
i.e. the minimal polynomial of quadgen(D). D must be an integer
congruent to 0 or 1 modulo 4, which is not a square.
? quadpoly(5,'y) %1 = y^2 - y - 1 ? quadpoly(0,'y) *** at top-level: quadpoly(0,'y) *** ^ — — — — -- *** quadpoly: domain error in quadpoly: issquare(disc) = 1
The library syntax is GEN quadpoly0(GEN D, long v = -1) where v is a variable number.
Relative equation for the ray
class field of conductor f for the quadratic field of discriminant D
using analytic methods. A bnf for x2 - D is also accepted in place
of D.
For D < 0, uses the σ function and Schertz's method.
For D > 0, uses Stark's conjecture, and a vector of relative equations may be
returned. See bnrstark for more details.
The library syntax is GEN quadray(GEN D, GEN f, long prec).
Regulator of the quadratic order of positive discriminant D in time
Õ(D1/2) using the continued fraction algorithm. Raise
an error if D is not a discriminant (fundamental or not) or if D is a
square. The function quadclassunit is asymptotically faster (and also
in practice for D > 1010 or so) but depends on the GRH.
The library syntax is GEN quadregulator(GEN D, long prec).
'w})
A fundamental unit u of the real quadratic order
of discriminant D. The integer D must be congruent to 0 or 1 modulo 4
and not a square; the result is a quadratic number (see Section se:quadgen).
If D is not a fundamental discriminant, the algorithm is wasteful: if D =
df2 with d fundamental, it will be faster to compute quadunit(d)
then raise it to the power quadunitindex(d,f); or keep it in
factored form.
If v is given, the variable name is used to display u
else 'w' is used. The algorithm computes the continued fraction
of (1 + sqrt{D}) / 2 or sqrt{D}/2 (see GTM 138, algorithm 5.7.2).
Although the continued fraction length is only O(sqrt{D}),
the function still runs in time Õ(D), in part because the
output size is not polynomially bounded in terms of log D.
See bnfinit and bnfunits for a better alternative for large
D, running in time subexponential in log D and returning the
fundamental units in compact form (as a short list of S-units of size
O(log D)3 raised to possibly large exponents).
The library syntax is GEN quadunit0(GEN D, long v = -1) where v is a variable number.
When v does not matter, the function
GEN quadunit(GEN D) is also available.
Given a fundamental discriminant D, returns the index of the unit group
of the order of conductor f in the units of ℚ(sqrt{D}). This function
uses the continued fraction algorithm and has O(D1/2 + ϵ
fϵ) complexity; quadclassunit is asymptotically faster but
depends on the GRH.
? quadunitindex(-3, 2) %1 = 3 ? quadunitindex(5, 2^32) \\ instantaneous %2 = 3221225472 ? quadregulator(5 * 2^64) / quadregulator(5) time = 3min, 1,488 ms. %3 = 3221225472.0000000000000000000000000000
The conductor f can be given in factored form or as
[f, factor(f)]:
? quadunitindex(5, [100, [2,2;5,2]]) %4 = 150 ? quadunitindex(5, 100) %5 = 150 ? quadunitindex(5, [2,2;5,2]) %6 = 150
If D is not fundamental, the result is undefined; you may use the following script instead:
index(d, f) =
{ my([D,F] = coredisc(d, 1));
quadunitindex(D, f * F) / quadunitindex(D, F)
}
? index(5 * 10^2, 10)
%7 = 10
The library syntax is GEN quadunitindex(GEN D, GEN f).
Returns the norm (1 or -1) of the fundamental unit of the quadratic
order of discriminant D. The integer D must be congruent to 0 or 1
modulo 4 and not a square. This is of course equal to norm(quadunit(D))
but faster.
? quadunitnorm(-3) \\ the result is always 1 in the imaginary case %1 = 1 ? quadunitnorm(5) %2 = -1 ? quadunitnorm(17345) %3 = -1 ? u = quadunit(17345) %4 = 299685042291 + 4585831442*w ? norm(u) %5 = -1
This function computes the parity of the continued fraction
expansion and runs in time Õ(D1/2). If D is fundamental,
the function bnfinit is asymptotically faster but depends of the GRH.
If D = df2 is not fundamental, it will usually be faster to first compute
quadunitindex(d, f). If it is even, the result is 1, else the result
is quadunitnorm(d). The narrow class number of the order of
discriminant D is equal to the class number if the unit norm is 1 and to
twice the class number otherwise.
Important remark. Assuming GRH, using bnfinit is much
faster, running in time subexponential in log D (instead of exponential
for quadunitnorm). We give examples for the maximal order:
? GRHunitnorm(bnf) = vecprod(bnfsignunit(bnf)[,1]) ? bnf = bnfinit(x^2 - 17345, 1); GRHunitnorm(bnf) %2 = -1 ? bnf = bnfinit(x^2 - nextprime(2^60), 1); GRHunitnorm(bnf) time = 119 ms. %3 = -1 ? quadunitnorm(nextprime(2^60)) time = 24,086 ms. %4 = -1
Note that if the result is -1, it is unconditional because (if GRH is false) it could happen that our tentative fundamental unit in bnf is actually a power uk of the true fundamental unit, but we would still have Norm(u) = -1 (and k odd). We can also remove the GRH assumption when the result is 1 with a little more work:
? v = bnfunits(bnf)[1][1] \\ a unit in factored form ? v[,2] %= 2; ? nfeltissquare(bnf, nffactorback(bnf, v)) %7 = 0
Under GRH, we know that v is the fundamental unit, but as
above it can be a power uk of the true fundamental unit u. But the
final two lines prove that v is not a square, hence k is odd and
Norm(u) must also be 1. We modified the factorization matrix
giving v by reducing all exponents modulo 2: this allows to computed
nffactorback even when the factorization involves huge exponents.
And of course the new v is a square if and only if the original one was.
The library syntax is long quadunitnorm(GEN D).