Essays/Gaussian Integers
Introduction
A Gaussian integer is a complex number whose real and imaginary parts are both integers. In the J complex datatype each atom is represented by two 64-bit floating point numbers. Here we use a list of two integers to represent a Gaussian integer. For example:
x=: 2 _5 ] y=: (!21x) , 6^29x 51090942171709440000 36845653286788892983296
The Gaussian integer x is one whose real part is 2 and imaginary part _5 , and the Gaussian integer y is one whose real part is 51090942171709440000 and whose imaginary part is 36845653286788892983296 .
Basic Operations
plus =: +"1 minus =: -"1 times =: (-/@:* , (+/@:* |.)) " 1 conjugate =: 1 _1 *"1 ]
For example:
3 4 plus 2 _5 5 _1 3 4 minus 2 _5 1 9 3 4 times 2 _5 26 _7 conjugate 3 4 3 _4 3 4 times conjugate 3 4 25 0 ] y=: _2 ]\ 3 5 7 _11 _13 17 _19 _23 3 5 7 _11 _13 17 _19 _23 3 5 times y _16 30 76 2 _124 _14 58 _164
Euclidean Ring
The Gaussian integers form an Euclidean ring: There is a function D=: +/@:*:"1 that applies to Gaussian integers to produce non-negative integers, such that for all non-zero Gaussian integers x and y ,
- (D x) <: D x times y
- There exist Gaussian integers q and r such that x = r plus q times y , where (D r)<:D y .
D =: +/@:*:"1 divide =: <. @ (1r2&+) @ ((times conjugate) % D@]) residue=: ] minus [ times divide~
q=:x divide y and r=:y residue x compute q and r that satisfy the above requirements.
x=: 3 4 y=: 2 _5 ] q=: x divide y 0 1 ] r=: y residue x _2 2 D y 29 D r 8 r plus q times y 3 4
GCDs exist in an Euclidean ring and can be computed using the Euclidean algorithm.
gi=: residue/\@|.^:(0 0-.@-:{:)
gi applies to a 2 2 matrix of two Gaussian integers and implements one iteration of the Euclidean algorithm. (It is the identity function if row 1 of the argument is zero.)
y=: 5 7 ,: 11 _13 gi y 11 _13 5 7 gi gi y 5 7 _3 _3 gi gi gi y _3 _3 _1 1 gi^:_ y _1 1 0 0 <"2 gi^:a: y ┌──────┬──────┬─────┬─────┬────┐ │ 5 7│11 _13│ 5 7│_3 _3│_1 1│ │11 _13│ 5 7│_3 _3│_1 1│ 0 0│ └──────┴──────┴─────┴─────┴────┘
GCDs are unique to within multiplication by a unit (one of 1 _1 0j1 0j_1 suitably represented). Of the four unit associates, we favor the one in the first quadrant.
U =: _2]\ 1 0 _1 0 0 1 0 _1
quad0=: ({~ * <./@i. (1 1,:1 0)"_) @ (U times ])
gcd =: quad0 @ {. @ (gi^:_) @ ,: " 1
For example:
5 7 gcd 11 _13
1 1
1 1 residue 5 7
0 0
1 1 residue 11 _13
0 0
] x=: 3 4 times 2 ?.@$ 10^23x
159575762180199862353798 358621571191149837960564
] y=: 3 4 times 2 ?.@$ 10^29x
_150624195482901355474791095802 417094671165568842905828047764
] c=: x gcd y
18 24
c residue x
0 0
c residue y
0 0
t=: gi^:a: x,:y
$t
52 2 2
_5 ,\ ' ' ,. ": ,. 2 %/\ x:^:_1 D {."2 t
7.8347e_13 1.27637e12 7.74701 5.50826 4.2625
3.90546 4.44771 50.3477 2.53547 29.567
5.84651 5.23286 375.061 14.861 6.25772
10.7852 9.0245 18.2838 9.9903 4.23536
8.85927 3.56985 4.39535 4.02444 2.99814
4.80653 9.3591 3.28471 5.87401 26.6918
5.17959 6.88408 6.54967 4.77162 4.44734
4.27171 25.6027 3.75665 5.77035 20.581
4.95218 4.8309 12.683 13.6983 3.40706
5.3091 4.68169 64.2331 5.35669 4.61765
34
The last result is a display of the ratios in D values of the divisors from one iteration to the next, and demonstrates the efficacy of divide and of the Euclidean algorithm.
Collected Definitions
plus =: +"1
minus =: -"1
times =: (-/@:* , (+/@:* |.)) " 1
conjugate =: 1 _1 *"1 ]
D =: +/@:*:"1
divide =: <. @ (1r2&+) @ ((times conjugate) % D@])
residue =: ] minus [ times divide~
gi =: residue/\@|.^:(0 0-.@-:{:)
U =: _2]\ 1 0 _1 0 0 1 0 _1
quad0 =: ({~ * <./@i. (1 1,:1 0)"_) @ (U times ])
gcd =: quad0 @ {. @ (gi^:_) @ ,: " 1
See also
Contributed by Roger Hui.