Re: question on finding exhaustive search algorithm for a^2+b^2 = N

[Bill Allombert <Bill.Allombert@math.u-bordeaux.fr>]

On Fri, Nov 17, 2023 at 09:19:22PM -0600, Kurt Foster wrote:
> On Nov 17, 2023, at 7:11 PM, American Citizen wrote:
> 
> > My goal is an exhaustive search of 2 squares = N algorithm
> 
> I dashed off some code which may be overkill (it absolutely clobbers the
> problem) but seems to get the job done.  I used your example as a test.
> 
> ? K=bnfinit(x^2+1);
> ? N=5*13*17*29*37;
> ? v=bnfisintnorm(K,N);
> ? w=vecsort(vector(#v,i,abs([polcoeff(v[i],0),polcoeff(v[i],1)])),1)
> %4 = [[64, 1087], [103, 1084], [167, 1076], [191, 1072], [236, 1063], [281,
> 1052], [292, 1049], [359, 1028], [449, 992], [512, 961], [568, 929], [601,
> 908], [607, 904], [664, 863], [673, 856], [743, 796], [796, 743], [856,
> 673], [863, 664], [904, 607], [908, 601], [929, 568], [961, 512], [992,
> 449], [1028, 359], [1049, 292], [1052, 281], [1063, 236], [1072, 191],
> [1076, 167], [1084, 103], [1087, 64]]

Interesting! It seems to work, but according to the documentation
"Computes  a complete system of solutions  (modulo units of positive norm)"
you should have had the same issue as with qfbsolve.
bnfisintnorm returns both 103*x-1084  and 1084*x-103 (with x=I)
but since I is a unit of positive norm, you could get in theory 
103*x-1084 and x*(1084*x-103) = -103+1084*x
which would be identical after taking the absolute values.

Cheers,
Bill.