Thank you for posting such an interesting article. Here is a method of calculating all 159,139 triples at once using 1000x1000 matrices and a lot of memory. the method is based on Theorem 11.1 page 394 from Kenneth Rosen's 1984 book "Elementary Number Theory and its Applications". . a = 1000; /* sqrt of maximum hypotenuse */ m = repeat(t(1:a), 1, a); n = repeat(1:a, a); good = loc( (gcd(m, n)=1) # mod(m - n, 2) # (m > n) # ((m#m + n#n)<=(a#a)) ); m = m[good]; n = n[good]; t = (2#m#n) || (m#m - n#n) || (m#m + n#n); create triples from t; append from t; . If I use this larger set of triples to make a distribution of angles, then the profile is virtually uniform between zero and 45 degrees, so I am thinking that there must be something specific to your algorithm that makes it avoid triangles that are similar in shape to the lower-order ones.