Karim Belabas on Sat, 06 Dec 2008 18:28:05 +0100

[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: nfsnf problems

* David Loeffler [2008-12-03 13:06]:
> I'm having some trouble with the "nfsnf" function.
> Firstly, in quite a lot of cases -- particularly where the matrix is
> close to being in Smith form already -- I get the error message "bug2 in
> nfsmith". This happens in both 2.3.3 and 2.4.2alpha. Secondly, when I do
> get an answer I'm not convinced it's the right one.
> Here's an example:
> E = nfinit(x^2 - x + 2);
> M = [1, 0, x; 0, x, 0; 0,0,2+x];
> MM = nfalgtobasis(E, M);
> N = [[1, 0; 0, 1], [1, 0; 0, 1], [1, 0; 0, 1]];
> nfsnf(E, [MM, N, N])
> If I've understood everything correctly, that should calculate the SNF
> of the matrix [1,0,w; 0,w,0; 0,0,2+w] where w^2 - w + 2 = 0. But it returns
> [[8, 3; 0, 1], [1, 0; 0, 1], [1, 0; 0, 1]]
> Converting these back to polmods, that's saying that the matrix has
> elementary factors 1, 1, and the ideal generated by 8 and w + 2 (since
> E.zk[1] = 1, E.zk[2] = x-1.) But it's clear just by looking at the
> matrix that its elementary factors are 1, w and 2 + w.
> If I change the [1,3] entry of M to 0 instead of x, it gives me the bug2
> in nfsmith error message.
> Any suggestions what's going on here?

1) The nfsnf function was not really specified in the PARI/GP documentation
so I am not sure what was intended in the implementation.

2) Taking as a reference Section 1.7 of Henri Cohen's GTM 193 (Advanced
Topics in Computational Number Theory), in particular Algorithm 1.7.3,
and reverse-engineering the code to try and link the two together,
I wrote a new specification (section 3.6.4 and 3.6.108 of the User's

3) Based on that specification, I fixed 3 major problems in the code:

-- one related to both the "bug2" and the "wrong result" problem [ rescale
the column so that diagonal term is 1, as the rest of the algorithm
apparently assumed ],

-- one which did not show up in your endeavours, and apparently replaced
b_i <--> b_i^{-1}. There's a comment in GTM 193 saying that one should
manipulate b_i^{-1} and not b_i to avoid ideal inversions; the code
assumed it had been done from the start (i.e the input had been suitably

-- another obscure one, which I introduced a long time ago while
cleaning this code without checking what it meant. (Basically: if some
test fail, go on computing garbage instead of restarting the main loop.
The test hardly ever failed in simple examples...)

One final problem, which I have not dealt with: the algorithm only makes
sense for "integral" pseudo-matrices ( a_{i,j} in B_i A_j^(-1) for all i,j ).
Neither the code, nor the regression tests seem to be interested in that
condition. (On the contrary, it is violated in every single regression test
involving nfsnf.)

I am not too sure what to do with this; for the time being, I do not check it
(rather costly), and let the algorithm go on even if I detect it is not
satisfied (would be very cheap but prevent existing code from "working"). In
that case, the elementary divisors (D_i) will not be integral ideals, and
statements like

  M ~ Z_K / D_1 + ... + Z_K/D_n

become rather pointless.

All the above was committed to the "testing / unstable" svn repository.

Thanks for your report !


P.S: in current svn, your example can be simplified to 

  E = nfinit(x^2 - x + 2);
  M = [1, 0, x; 0, x, 0; 0,0,2+x];
  N = [1, 1, 1];
  nfsnf(E, [M, N, N])

... and returns the expected

%4 = [[8, 2; 0, 1], [2, 0; 0, 1], [1, 0; 0, 1]]

In fact, disregarding the possible simplification, it *must* be changed since
nfalgtobasis no longer applies recursively on t_VEC / t_COL / t_MAT, so

  MM = nfalgtobasis(E, M);

would no longer work. The function

  MM = matalgtobasis(E, M);

should be used instead. But in fact, this step has become completely
unnecessary since number field elements in arbitrary form are now accepted in
relative matrices. Just like "ideal lists" can now accept ideals in any form, 
hence the "1" above instead of matid(2).

P.S2: The reason why nfalgtobasis was changed to no longer be recursive is that
now it can convert arbitrary number field elements to "basis form". (Just
like nfbasistoalg converted arbitrary elements to "algebraic form".)

Before the change, it converted elements *NOT* already in basis form to basis
form, and returned "junk" (vector of vectors) when applied to elements in
"basis form".

Karim Belabas, IMB (UMR 5251)  Tel: (+33) (0)5 40 00 26 17
Universite Bordeaux 1          Fax: (+33) (0)5 40 00 69 50
351, cours de la Liberation    http://www.math.u-bordeaux.fr/~belabas/
F-33405 Talence (France)       http://pari.math.u-bordeaux.fr/  [PARI/GP]