Programming Praxis

5 Responses to “Pythagorean Quadruple”

1. Zack said

   February 25, 2020 at 10:20 AM

   Interesting problem. Here is my take on it using Julia: https://pastebin.com/eYYR51sm

   For N = 1000 I found 85490 quadruples (allowing for multiple instances of the same number in the quadruple set).

   Cheers
2. matthew said

   February 26, 2020 at 12:07 AM

   Here’s a C++ solution that runs in O(n²) time, using O(n²) space (assuming hash table access time is constant):

   #include <vector>
   #include <utility>
   #include <iostream>
   #include <cstdio>
   
   // Find all "pythagorean quadruples" (a,b,c,d) with
   // 1 <= a <= b <= c <= N and
   // a^2 + b^2 + c^2 = d^2
   
   // Meet-in-the-middle: build table of all possible ways
   // of making a*a+b*b, then search all possible values
   // of d*d-c*c
   
   int main(int argc, char *argv[]) {
     int n = strtoul(argv[1],0,0);
     std::unordered_map<int,std::vector<std::pair<int,int>>> m;
     for (int a = 1; a <= n; a++) {
       for (int b = a; b <= n; b++) {
         m[a*a+b*b].push_back(std::make_pair(a,b));
       }
     }
     // sqrt(3c^2) = sqrt(3)*c = 2c, near enough
     for (int c = 1; c <= n; c++) {
       for (int d = c+1; d <= 2*c; d++) {
         for (auto p: m[d*d-c*c]) {
           auto a = p.first, b = p.second;
           if (c >= p.second) {
             std::cout << a << b << c << d << "\n";
           }
         }
       }
     }
   }
   

   $ g++ -O3 -Wall pquads.cpp -o pquads
   $ time ./pquads 1000 | wc -l
   85490
   
   real	0m0.468s
   user	0m0.443s
   sys	0m0.032s
   
3. matthew said

   February 26, 2020 at 10:36 AM

   Just in case overflow in all those squares is a problem, I thought I’d see what gcc can do to check things at runtime:

   $ g++ -O3 -Wall -fsanitize=signed-integer-overflow pquads.cpp -o pquads
   $ time ./pquads 1000 | wc -l
   85490
   
   real 0m0.447s
   user 0m0.426s
   sys  0m0.029s
   $ g++ -O3 -Wall -fsanitize=undefined pquads.cpp -o pquads
   $ time ./pquads 1000 | wc -l
   85490
   
   real	0m0.651s
   user 0m0.627s
   sys  0m0.032s
   
4. matthew said

   February 26, 2020 at 1:48 PM

   Here’s a slightly improved version (that prints spaces between the values for a,b,c,d):

   int main(int argc, char *argv[]) {
     int n = strtoul(argv[1],0,0);
     std::unordered_map<int,std::vector<std::pair<int,int>>> m;
     for (int a = 1; a <= n; a++) {
       for (int b = a; b <= n; b++) {
         m[a*a+b*b].emplace_back(a,b);
       }
     }
     for (int c = 1; c <= n; c++) {
       for (int d = c+1; d <= 2*c; d++) {
         for (auto p: m[d*d-c*c]) {
           auto a = p.first, b = p.second;
           if (b <= c) std::cout << a << " " << b << " " << c << " " << d << "\n";
         }
       }
     }
   }
   
5. Daniel said

   March 18, 2020 at 1:13 AM

   Here’s a solution in Python.

   from collections import defaultdict
   from itertools import count
   
   def pythagorean_quadruples(N):
       bc_lookup = defaultdict(list)
       for b in range(1, N + 1):
           for c in range(b, N + 1):
               bc_lookup[b * b + c * c].append((b, c))
       output = 0
       for a in range(1, N + 1):
           for d in count(a):
               if d * d - a * a > 2 * N * N:
                   break
               for b, c in bc_lookup.get(d * d - a * a, []):
                   if a <= b:
                       output += 1
       return output
   
   print(pythagorean_quadruples(1000))
   

   Output:

   85490
   

Leave a Reply

Enter your comment here...

Fill in your details below or click an icon to log in:

Email (required) (Address never made public)

Name (required)

Website

You are commenting using your WordPress.com account. ( Log Out /  Change )

You are commenting using your Facebook account. ( Log Out /  Change )

Cancel

Connecting to %s

Notify me of new comments via email.

Notify me of new posts via email.

Δ