c – Babbage Problem – squares ending in digits 269696

//the strategy of take the rest of division by 1e06 is
//to take the a number how 6 last digits are 269696
    while (((square=current*current) % 1000000 != 269696) && (square<INT_MAX)) {
        current++;
    }

This is part of the “C” example on rosettacode.org. The comments, and not only them, more belong to towerofbabel.disorg!

But rosettacode has a very nice presentation:

What is the smallest positive integer whose square ends in the digits 269,696?

…Babbage asked in a letter, to give an example of what his yet-to-build engine could be working on.

He thought the answer might be 99,736, whose square is 9,947,269,696;
but he couldn’t be certain.

I hope this is all true, because the story is almost too good: I really wonder how he got at this solution; it is correct, but there is a much smaller one.

It is suggested you write it as if for Mr. Babbage himself, who seems to have a clever pen-and-paper method, and who knows the basics (only the basics, but very well).

I found some rare programs that only check roots ending on 4 or 6, because the end of the 269696 ending is a 6. I extended this keenly to endings 64 and 36, which both only give squares ending on 96.

Here is my code. I turned it around and first calculate an approximate root for every number with that ending. At least I treat it as approximate.

My main Q is about the rounding and the way I assign and use babb and diff. Before I used round() it was not really working; now I removed all casts and it seems to work.

/* "Babbage Problem"
    = (Smallest) number whose square ends in ...269696 ? */
/* The root must end on 4 or 6 to give ending 6,
   but also on 36 or 64 to give ending 96 ?!? */

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

int main() {
    const int ENDING = 269696;
    const int EXPMAX = 38;
    const double mu = 0.001;

    long n,
         babb;       /* nearest integer to root */

    double root,     /* approximation, unless sqrt() is used */
           diff,     /* between root and babb */
           modroot;  /* root wuth last int digits only */

    for (n = ENDING; n < 1L<<EXPMAX; n += 1000*1000) {

        /* sqrt() is faster than exp(log()/2) and has no fraction/diff at all if really integer */
        root = exp(log(n)/2);
        //root = sqrt((double)n); 

        babb = round(root);
        diff = root - babb;
        /* mod 100 with 36 and 64, or mod 10 with 4 and 6 */
        modroot = babb % 100 + diff;
        if (fabs(36 - modroot) < mu ||
            fabs(64 - modroot) < mu   ) {

            /* Check with integer division */
            if (n % babb == 0)
                putchar('*');
            else
                putchar(' ');

            printf("%16ld %20.12f %12ld %20.12fn", n, root, babb, diff);
        }
    }
    return 0;
}

Output:

*       638269696   25264.000000000004        25264       0.000000000004
*      9947269696   99735.999999999927        99736      -0.000000000073
*     22579269696  150263.999999999971       150264      -0.000000000029
*     50506269696  224735.999999999884       224736      -0.000000000116
      55020269696  234563.999147354130       234564      -0.000852645870
      70456269696  265435.999246522610       265436      -0.000753477390
*     75770269696  275263.999999999942       275264      -0.000000000058
*    122315269696  349736.000000000175       349736       0.000000000175
     129286269696  359563.999443770794       359564      -0.000556229206
     152440269696  390435.999487752211       390436      -0.000512247789
*    160211269696  400264.000000000058       400264       0.000000000058
*    225374269696  474735.999999999651       474736      -0.000000000349
     234802269696  484563.999587257451       484564      -0.000412742549
     265674269696  515435.999611978768       515436      -0.000388021232

So Babbage’s pen-and-paper 99736 is the second one, being at 9947 million. But the first one is after 638 iterations/millions. He came from somewhere else. This is with exp(log()).

The near misses (no stars) are also interesting. With a smaller mu, n up to 2^46 and sqrt() the last lines are:

    *  67646282269696 8224736.000000000000      8224736       0.000000000000
       67808044269696 8234563.999975712039      8234564      -0.000024287961
       68317432269696 8265435.999975803308      8265436      -0.000024196692
    *  68479994269696 8275264.000000000000      8275264       0.000000000000
    *  69718091269696 8349736.000000000000      8349736       0.000000000000
       69882310269696 8359563.999976075254      8359564      -0.000023924746

Does sqrt() always return .00000 ?

Is there something like exp(log()/2) but even simpler?

I only need 3 or 4 significant digits around the decimal point. An imprecise but specialized square root (or log) function. I don’t really understand that “binary estimation” for the seed value. How can I think about 2^n without log2()?