bigsqrt - Python vs Pike

[Skaperen]

this is a commandline program i wrote many years ago in both Python (my first Python program) and Pike. both are given below.  the Pike code is substantially faster.  this program does its work using big ints, so that ends up being what is compared.  give it one argument, the number to find the square root of.  if no arguments are given, it will compute the golden ratio.  it just keeps outputting digits to get more and more precise (never ending).

My code

#!/usr/bin/python
# -*- coding: ISO-8859-1 -*-
#-----------------------------------------------------------------------------
# Copyright � 2009 - Philip Howard - All rights reserved
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 2
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
#-----------------------------------------------------------------------------
# program    sqrt_big
#
# purpose    Calculate lots of digits for the square root of a given number
#        using a scaled integer adaptation of Newton's method.
#
# usage        sqrt_big [number [digits]]
#
# note        The number defaults to zero.  Since the square root of zero
#        is not defined, this program calculates the Golden Ratio as
#        a special case.
#
# note        The number of digits defaults to zero.  If it is zero, this
#        program keeps running until something stops it, like user
#        intervention or out of memory.
#-----------------------------------------------------------------------------
import sys
import math

stderr = sys.stderr
stdout = sys.stdout

def main():
    write_limit = 10000

    # Get 0 or 1 or 2 arguments (defaults: 0 and 0)
    this_arg = '0'
    digits_max = 0
    if len(sys.argv) > 1:
        this_arg = sys.argv[1]
        if len(sys.argv) > 2:
            digits_max = int(sys.argv[2])

    num_parts = this_arg.split('.')
    if len(num_parts) == 1:
        num_parts = num_parts + ['0']
    if len(num_parts) != 2:
        stderr.write( "Number needs 0 or 1 '.'\n" )
    return 1
    if len(num_parts[0]) < 1:
        num_parts[0] = '0'
    if len(num_parts[1]) < 8:
        num_parts[1] = (num_parts[1] + '00000000')[0:8]
    num_whole = int( num_parts[0] )
    num_fract = int( num_parts[1] )
    num_scale = 10 ** len( num_parts[1] )

    if num_whole < 0 or num_fract < 0:
    stderr.write( "Number is negative\n" )
    return 1

    if num_whole == 0 and num_fract == 0:
        goldratio =         1
        num       = 125000000
        num_scale = 100000000
    else:
        goldratio =         0
        num       = num_whole * num_scale + num_fract

    # Use floating point sqrt() function to make an estimate.
    num_float = float(num)
    num_float /= float(num_scale)
    num_float = math.sqrt(num_float)
    num_float *= float(num_scale)
    root = int(num_float)
    root_scale = num_scale
    root_str = str(root)[0:1]

    # Do 12 cycles with short precision growth.
    for n in range( 0, 12 ):
        root = ( ( num * root_scale * root_scale + num_scale * root * root ) * root_scale ) / ( root + root )
        root /= root_scale
        root_scale *= num_scale
        root_str = str(root)

    digits = 1 + len(str(root)) - len(str(root_scale))
    if goldratio:
        # Fake the addition of 0.5 for the Golden Ratio.
        stdout.write( '1.6' )
        digits_done = 2
    elif digits > 0:
        # Output just the whole digits to get the fraction point in there.
        stdout.write( root_str[0:digits] )
        stdout.write( '.' )
        digits_done = digits
    elif digits <= 0:
    # Output enough 0 digits to align the fraction.
        stdout.write( '0.' )
    stdout.write( '0'*(-digits) )
        digits_done = 0

    # Do remaining cycles with long precision growth.
    for n in range( 0, 60 ):
        root = ( ( num * root_scale * root_scale + num_scale * root * root ) * root_scale ) / ( root + root )
        root_scale = root_scale * root_scale * num_scale
        root_str = str(root)
        digits = len(root_str) * 7 / 8;

        # Output the digits, with a limit per write, up to the maximum.
        while digits_done < digits:
            count = digits - digits_done
            if count > write_limit:
                count = write_limit
            if digits_max > 0:
                remain = digits_max - digits_done
                if count > remain:
                    count = remain
            stdout.write( root_str[ digits_done : digits_done + count ] )
            digits_done += count
            if digits_max > 0 and digits_done >= digits_max:
                break
        if digits_max > 0 and digits_done >= digits_max:
            break
    stdout.write( '\n' )
    return 0

main()

#!/usr/bin/env pike
//-----------------------------------------------------------------------------
// Copyright � 2009 - Philip Howard - All rights reserved
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License
// as published by the Free Software Foundation; either version 2
// of the License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
//-----------------------------------------------------------------------------
// program    sqrt_big.pike
//
// purpose    Calculate lots of digits for the square root of a given number
//        using a scaled integer adaptation of Newton's method.
//
// usage    [pike] sqrt_big.pike [number [digits]]
//
// note        The number defaults to zero.  If it is zero, this program
//        calculates the Golden Ratio as a special case.
//
// note        The number of digits defaults to zero.  If it is zero, this
//        program keeps running until it is out of memory.
//-----------------------------------------------------------------------------
int main( int argc, array(string) arga )
{
    array    num_parts;
    string    this_arg,root_str;
    float    num_float;
    int        num,num_whole,num_fract,num_scale;
    int        root,root_scale;
    int        digits,digits_done,digits_max;
    int        goldratio;
    int        write_limit;


    write_limit = 10000;

    // Get 0 or 1 or 2 arguments (defaults: 0 and 0)
    this_arg = argc <= 1 ? "0" : arga[1];
    digits_max = argc <= 2 ? 0 : (int) arga[2];

    // Break apart the argument number.
    num_parts = this_arg / ".";
    if ( sizeof( num_parts ) == 1 ) num_parts = num_parts + ({ "0" });
    if ( sizeof( num_parts ) != 2 ) return 1;
    if ( strlen( num_parts[0] ) < 1 ) num_parts[0] = "0";
    if ( strlen( num_parts[1] ) < 8 ) num_parts[1] = ( num_parts[1] + "00000000" )[0..7];
    num_whole = (int) num_parts[0];
    num_fract = (int) num_parts[1];
    num_scale = 10->pow( strlen( num_parts[1] ) );

    // Do not calculate square root of negative number.
    if ( num_whole < 0 ) return 1;
    if ( num_fract < 0 ) return 1;

    // If zero is specified, calculate the golden ratio instead.
    if ( num_whole == 0 && num_fract == 0 ) {
    // To find golden ratio, find square root of 1.25 and add 0.5 to result.
    goldratio =         1;
    num       = 125000000;
    num_scale = 100000000;
    } else {
    goldratio =         0;
    num       = num_whole * num_scale + num_fract;
    }

    // Use floating point sqrt() function to make an estimate.
    num_float = (float) num;
    num_float /= (float) num_scale;
    num_float = sqrt( num_float );
    num_float *= (float) num_scale;
    root = (int) num_float;
    root_scale = num_scale;
    root_str = ( (string) root )[0..0];

    // Do 13 cycles with short precision growth.
    for ( int n = 0; n < 13; ++ n ) {
    root = ( ( num * root_scale * root_scale + num_scale * root * root ) * root_scale ) / ( root + root );
    root /= root_scale;
    root_scale *= num_scale;
    root_str = (string) root;
    }

    // Display the whole part first.
    digits = 1 + strlen( (string) root ) - strlen( (string) root_scale );
    if ( goldratio ) {
    // Fake the addition of 0.5 for the Golden Ratio.
    write( "1.6" );
    digits_done = 2;
    } else if ( digits > 0 ) {
    // Output just the whole digits to get the fraction point in there.
    write( "%s.", root_str[ 0 .. digits - 1 ] );
    digits_done = digits;
    } else if ( digits <= 0 ) {
    // Output enough 0 digits to align the fraction.
    write( "0.%s", "0" * (-digits) );
    digits_done = 0;
    }

    // Do remaining cycles with long precision growth.
    for ( int n = 0; n < 60; ++ n ) {
    string root_str;
    root = ( ( num * root_scale * root_scale + num_scale * root * root ) * root_scale ) / ( root + root );
    root_scale = root_scale * root_scale * num_scale;

    root_str = (string) root;
    digits = ( strlen( root_str ) * 14 ) / 15;

    // Output the digits, with a limit per write, up to the maximum.
    while ( digits_done < digits ) {
        int count;

        count = digits - digits_done;
        if ( count > write_limit ) count = write_limit;
        if ( digits_max > 0 ) {
        int remain;
        remain = digits_max - digits_done;
        if ( count > remain ) count = remain;
        }
        write( "%s", root_str[ digits_done .. digits_done + count - 1 ] );
        digits_done += count;
        if ( digits_max > 0 && digits_done >= digits_max ) break;
    }
    if ( digits_max > 0 && digits_done >= digits_max ) break;
    }

    write( "\n" );
    return 0;
}