Pythagorean triplet bruteforce algorithm breaks at 969

Question

I have this code written in C to find all possible Pythagorean triplets within a certain number range. The original algorithm I wrote, just nested for loops and if(pow(a, 2) + pow(b, 2) == pow(c, 2)), worked just fine. However, my new, more optimized algorithm, which only loops through a and b, setting c to sqrt(pow(a, 2) + pow(b, 2)) and checking if ceil(c) == c, immediately begins accepting all whole numbers the moment b hits 969.

Furthermore, when I ran the second algorithm with a smaller amount and checked the results, the amount of triplets output is exactly 0.969x the limit of b (the top loop).

This is an extremely weird and interesting phenomenon, and I am unsure what makes 969 such a special number.

My new algorithm:

#include <stdio.h>
#include <math.h>
#include <string.h>

int main(void) {
    unsigned int limit;

    printf("Max value to bruteforce: ");
    scanf("%d", &limit);

    unsigned int triples[limit][3]; // `0.969 * limit` is exact
    unsigned int i = 0;

    printf("Bruteforcing...\n");

    for(unsigned int b = 1; b < limit; b++) {
        for(unsigned int a = 1; a <= b; a++) {
            double c = sqrt(pow(a, 2) + pow(b, 2));
            
            if(ceil(c) == c) {
                triples[i][0] = a;
                triples[i][1] = b;
                triples[i][2] = c;
                    
                i++;

                printf("found: %d, %d, %d\n", a, b, (int) c);
            }
        }
    }
    
    char out[i + 15];

    sprintf(out, "%d\n", limit);
    
    for(int i2 = 0; i2 < i; i2++) {
        for(int j = 0; j < 3; j++) {
            char ln[10];
            sprintf(ln, "%d ", triples[i2][j]);

            strcat(out, ln);
        }

        strcat(out, "\n");
    }

    FILE *f = fopen("cache.txt", "w");
    fprintf(f, "%s", out);
    fclose(f);

    printf("Saved to cache.txt\n");
  return 0;
}

My old algorithm (reproduced, may not be 100% accurate):

#include <stdio.h>
#include <math.h>
#include <string.h>

int main(void) {
    unsigned int limit;

    printf("Max value to bruteforce: ");
    scanf("%d", &limit);

    unsigned int triples[limit][3]; // `0.969 * limit` is exact
    unsigned int i = 0;

    printf("Bruteforcing...\n");

    for(unsigned int b = 1; b < limit; b++) {
        for(unsigned int a = 1; a <= b; a++) {
            for(unsigned int c = 1; c < limit; c++) {
                if(pow(a, 2) + pow(b, 2) == pow(c, 2)) {
                    triples[i][0] = a;
                    triples[i][1] = b;
                    triples[i][2] = c;
                    
                    i++;

                    printf("found: %d, %d, %d\n", a, b, (int) c);
                }
            }
        }
    }
    
    char out[i + 15];

    sprintf(out, "%d\n", limit);
    
    for(int i2 = 0; i2 < i; i2++) {
        for(int j = 0; j < 3; j++) {
            char ln[10];
            sprintf(ln, "%d ", triples[i2][j]);

            strcat(out, ln);
        }

        strcat(out, "\n");
    }

    FILE *f = fopen("cache.txt", "w");
    fprintf(f, "%s", out);
    fclose(f);

    printf("Saved to cache.txt\n");
  return 0;
}

Answer 1

I am unsure what makes 969 such a special number.

Nothing, actually. The posted code has undefined behavior caused by access out of bounds of triplets and the wrong format specifier for limit.

unsigned int limit;
scanf("%d", &limit);
//     ^^  It should be %u

unsigned int triples[limit][3]; // `0.969 * limit` is exact
//                   ^^^^^         ^^^^^^^^^^^^^^^^^^^^^^^^ Nope.

limit is used to bound the value of the bigger cathetus, but it's not enough as count of the number of triplets.

We can write a more simple (IMHO) version of the brute force algorithm, without using any floating-point function from <math.h>:

#include <stdio.h>

unsigned long long sq(unsigned long long x)
{
    return x * x;
}

int main(void)
{
    unsigned int limit;

    printf("Max value to bruteforce: ");
    scanf("%u", &limit);    // input validation is left to the reader.

    puts("\nBruteforcing...");
    unsigned count = 0;

    // The produced triplets will be ordered by the hypotenuse
    // and then the catheti.
    for(unsigned int c = 5; c < limit; ++c)
    {
        for(unsigned int a = 3, b = c - 1; a < b; --b)
        { //                    ^^^^^^^^^  ^^^^^  ^^^
          // The catheti must be greater than 0,
          // different from the hypotenuse and
          // from each other (sqrt(2) is irrational, so...).
            unsigned long long dif = sq(c) - sq(b);
            while ( sq(a) < dif )
            {
                ++a;
            }
            if ( sq(a) == dif )
            {
                ++count;
                printf("%4u %4u %4u\n", a, b, c);
                ++a;
            }
        }
    }
    printf("%u triplets found.\n", count);
    return 0;
}

Live: https://godbolt.org/z/hvTa9zKW4