Project Euler 153

/**
 * Copyright uwi ( http://wonderfl.net/user/uwi )
 * MIT License ( http://www.opensource.org/licenses/mit-license.php )
 * Downloaded from: http://wonderfl.net/c/3Zi8
 */

package {
    import flash.display.Sprite;
    import flash.text.TextField;
    import flash.utils.getTimer;
    // @see http://projecteuler.net/index.php?section=problems&id=153
    public class Euler153 extends Sprite {
        private var _tf : TextField;
  
        public function Euler153() {
            _tf = new TextField();
            _tf.width = 465;
            _tf.height = 465;
            addChild(_tf);
            
            var s : int = getTimer();
            t(solve(100000));
            var g : int = getTimer();
            t((g - s) + " ms");
        }

        // N/(a+bi)=N(a-bi)/(a^2+b^2)
        // N,a,bが互いに素とすると、a^2+b^2はNの約数。
        // (a,b)=mとすると、(a^2+b^2)/mはNの約数。
        // KがNの約数で、1<=N<=10^8のとき、該当するNの個数は[10^8/K]個
        // よって求める値は、上限をSとして
        // Σ_{a>0,b,a^2+b^2<=S,(a,b)=1}(aΣ_m[S/m(a^2+b^2)]m) [[複素数]]
        // + Σ_c [S/c]c  [[整数]]
        // S=10^8は30秒以上かかる・・
        private function solve(S : int) : MPI
        {
            var cache : Array = [];
            var primes : Array = doEratosthenes(Math.sqrt(S));
            
            var sum : MPI = new MPI(0);
            for(var a : int = 1;a * a <= S;a++){
                var bsup : int = Math.sqrt(S - a * a);
                var bs : Array = new Array(bsup + 1);
                // 互いに素マップ
                if(a > 1){
                    var sa : int = Math.sqrt(a);
                    for each(var p : int in primes){
                        if(p > a)break;
                        if(a % p == 0){
                            var j : int;
                            for(j = p;j <= bsup;j+=p){
                                bs[j] = 1;
                            }
                        }
                    }
                    /*
                    for(var i : int = 1;i <= sa;i++){
                        if(a % i == 0){
                            if(i > 1){
                                var j : int;
                                for(j = i;j <= bsup;j+=i){
                                    bs[j] = 1;
                                }
                            }
                            var q : int = a / i;
                            for(j = q;j <= bsup;j+=q){
                                bs[j] = 1;
                            }
                        }
                    }
                    */
                }
                
                var llsum : Number = 0;
                for(var b : int = 1;b <= bsup;b++){
                    if(bs[b])continue;
                    var lsum : Number;
                    var msup : int = S / (a * a + b * b);
                    if(cache[msup]){
                        lsum = cache[msup];
                    }else{
                        lsum = 0;
                        for(var m : int = 1;m <= msup;m++){
                            lsum += int(msup / m) * m;
                        }
                        cache[msup] = lsum;
                    }
                    llsum += lsum;
                }
                llsum *= a;
                sum.add(new MPI(llsum));
            }
            sum.mul(new MPI(2));
            
            var lllsum : Number = 0;
            for(var c : int = 1;c <= S;c++){
                lllsum += int(S / c) * c;
            }
            sum.add(new MPI(lllsum));
            
            return sum;
        }
        
        private function GCD(a : int, b : int) : int
        {
            return b == 0 ? a : GCD(b, a % b);
        }

	private function t(o : *) : void
	{
            _tf.appendText("" + o + "\n");
	}
    
        private static function doEratosthenes(n : int) : Array
        {
            var ar : Vector.<uint> = new Vector.<uint>(n / 2 - 1);
            var i : int;
            for(i = 0;i < ar.length;i++)ar[i] = 1;
            
            var sq : int = (Math.sqrt(n) - 3) >> 1;
            for(var p : int = 0;p <= sq;p++){
                if(ar[p] == 1){
                    var m : int = (p << 1) + 3;
                    var m2 : int = m << 1;
                    for(var mm : int = m * m;mm <= n;mm += m2){
                        ar[(mm - 3) >> 1] = 0;
                    }
                }
            }
            var ret : Array = [2];
            for(i = 0;i < ar.length;i++){
                if(ar[i] == 1)ret.push((i << 1) + 3);
            }
            return ret;
        }
    }
}

class MPI
{
    public var d : Vector.<uint>;
    
    public function MPI(v : Number = 0)
    {
        if(v == 0){
            d = Vector.<uint>([0]);
        }else{
            d = Vector.<uint>([]);
            for(var i : Number = v;i != 0;i = Math.floor(i / 10))d.push(i % 10);
        }
    }
    
    public function carry() : void
    {
        for(var i : int = 0;i < d.length;i++){
            var c : uint = d[i] / 10;
            if(c > 0){
                if(i == d.length - 1){
                    d.push(0);
                }
                d[i + 1] += c;
            }
            d[i] %= 10;
        }
    }
    
    public function add(n : MPI) : MPI
    {
        var i : int;
        for(i = 0;i < Math.min(d.length, n.d.length);i++){
            d[i] += n.d[i];
        }
        for(i = d.length;i < n.d.length;i++){
            d.push(n.d[i]);
        }
        carry();
        return this;
    }
    
    public function mul(n : MPI) : MPI
    {
        var ret : Vector.<uint> = new Vector.<uint>(d.length + n.d.length);
        var i : int, j : int;
        for(i = 0;i < ret.length;i++)ret[i] = 0;
        for(i = 0;i < d.length;i++){
            for(j = 0;j < n.d.length;j++){
                ret[i + j] += d[i] * n.d[j];
            }
        }
        d = ret;
        carry();
        shave();
        return this;
    }
    
    public function shave() : void
    {
        for(var i : int = d.length - 1;i >= 1;i--){
            if(d[i] > 0)break;
            d.pop();
        }
    }
    
    public function toString() : String
    {
        var ret : String = "";
        for(var i : int = d.length - 1;i >= 0;i--){
            ret += d[i].toString();
        }
        return ret;
    } 
}