Project Euler 66

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

package {
    import flash.display.Sprite;
    import flash.text.TextField;
    import flash.utils.getTimer;
    // @see http://projecteuler.net/index.php?section=problems&id=66
    public class Euler66 extends Sprite {
        private var _tf : TextField;
  
        public function Euler66() {
            _tf = new TextField();
            _tf.width = 465;
            _tf.height = 465;
            addChild(_tf);
            
            var s : int = getTimer();
            tr(solve(1000));
            var g : int = getTimer();
            tr((g - s) + " ms");
        }
        
        // @see http://www004.upp.so-net.ne.jp/s_honma/pell/pell.htm
        // まず解を与えるインデックスを求めて、そこまでのkの配列を求める。 (bsec)
        // kから、xの値を求め、x^2-Dy^2=-1の場合は、x^2+Dy^2が最小解、1の場合はそのままxが最小解となる。
        private function solve(M : int) : int
        {
            var max : Number = 0;
            var maxd : int = 0;
            
            var maxrep : MPI = new MPI(0);
            var maxD : int = -1;
            
            for(var D : int = 1, sq : int = 1;D <= M;D++){
                if(sq * sq == D){
                    sq++;
                    continue;
                }
               var ks : Array = bseq(D);
                var rep : Array = repA(ks);
                // 簡易的に1の位だけみる
                var r : int = (rep[0][0] % 10) * (rep[0][0] % 10) - D * (rep[1][0] % 10) * (rep[1][0] % 10);
                var val : MPI;
                if(r % 10 == -1){
                		// x^2-Dy^2=-1の場合。
                		// 解の整数項は、x^2+Dy^2
                		var xx : MPI = new MPI(rep[0][1] + "000000000000000");
                		xx.addint_(rep[0][0]);
                		var yy : MPI = new MPI(rep[1][1] + "000000000000000");
                		yy.addint_(rep[1][0]);
                		val = MPI.add(MPI.mul(xx, xx), MPI.mul(MPI.mul(yy, yy), new MPI(D)));
               }else{
                		// x^2-Dy^2=1の場合。
                		val = new MPI(rep[0][1] + "000000000000000");
                		val.addint_(rep[0][0]);
               }
               if(MPI.comp(val, maxrep) > 0){
               	maxrep = val;
               	maxD = D;
               }
            }
            tr(maxrep);
            return maxD;
        }
        
        private function repA(ks : Array) : Array
        {
        		var pr : Array = [0, 0];
        		var cr : Array = [1, 0];
        		var ps : Array = [1, 0];
        		var cs : Array = [0, 0];
        		for each(var k : uint in ks){
        			var nr : Array = pr.concat();
        			var ns : Array = ps.concat();
        			for(var i : uint = 0;i < k;i++){
        				nr = N2.add(nr, cr);
        				ns = N2.add(ns, cs);
        			}
        			pr = cr; cr = nr;
        			ps = cs; cs = ns;
        		}
        		return [cr, cs];
        }
        
        private function bseq(n : int) : Array
        {
        		var r : int = 0;
        		var s : int = 1;
        		var t : int = 1;
        		var ret : Array = [];
        		var sq : Number = Math.sqrt(n);
        		for(var i : int = 1;i < 100;i++){
        			var k : int = Math.floor((r + s * sq) / t);
        			ret.push(k);
        			var nr : int = t * (r - k * t);
        			var ns : int = -t * s;
        			var nt : int = (r - k * t) * (r - k * t) - s * s * n;
        			if(nt < 0){
        				nr = -nr;
        				ns = -ns;
        				nt = -nt;
        			}
        			var g : int = GCD(GCD(nr, ns), nt);
        			nr /= g;
        			ns /= g;
        			nt /= g;
//        			tr([nr, ns, nt]);
        			if(nt == 1)break;
        			r = nr;
        			s = ns;
        			t = nt;
        		}
        		return ret;
        }
        
        private function GCD(a : int, b : int) : int
        {
        		while(b > 0){
        			var c : int = a;
        			a = b;
        			b = c % b;
        		}
        		return a;
        }
        
        private function tr(...o : Array) : void
        {
            _tf.appendText(o + "\n");
            _tf.scrollV = _tf.maxScrollV;
        }
    }
}

class N2{
	public static const N : Number = 1E15;	
	public static function add(a : Array, b : Array) : Array
	{
		var r0 : Number = a[0] + b[0];
		var rp : int = Math.floor(r0 / N);
		r0 -= rp * N;
		var r1 : Number = a[1] + b[1] + rp;
		return [r0, r1];
	}
	
	public static function sub(a : Array, b : Array) : Array
	{
		var r0 : Number = a[0] - b[0];
		var rp : int = Math.floor(r0 / N);
		r0 -= rp * N;
		var r1 : Number = a[1] - b[1] + rp;
		return [r0, r1];
	}
	
	public static function inc(a : Array, b : Number) : Array
	{
		var r0 : Number = a[0] + b;
		var rp : int = Math.floor(r0 / N);
		r0 -= rp * N;
		var r1 : Number = a[1] + rp;
		
		a[0] = r0;
		a[1] = r1;
		return a;
	}
	
	public static function comp(a : Array, b : Array) : Number
	{
		if(a[1] != b[1])return a[1] - b[1];
		return a[0] - b[0];
	}
}

class MPI
{
    public var d : Vector.<uint>;
    
    public function MPI(v : * = 0)
    {
        if(v is Number){
            d = Vector.<uint>([]);
            for(var n : Number = v;n != 0;n = Math.floor(n / 10))d.push(n % 10);
            if(v == 0)d.push(0);
        }else if(v is String){
            d = Vector.<uint>([]);
            for(var i : int = v.length - 1;i >= 0;i--){
                d.push(uint(v.charAt(i)));
            }
        }else if(v is Array){
            d = Vector.<uint>(v);
        }else if(v is Vector.<uint>){
            d = v.concat();
        }else if(v is MPI){
            d = v.d.concat();
        }
        if(d)shave();
    }
    
    public static function add(a : MPI, b : MPI) : MPI
    {
        if(b == null)return a;
        
        var ret : MPI = new MPI(a);
        var i : uint;
        var minl : uint = Math.min(ret.d.length, b.d.length);
        for(i = 0;i < minl;i++){
            ret.d[i] += b.d[i];
        }
        for(i = ret.d.length;i < b.d.length;i++){
            ret.d.push(b.d[i]);
        }
        ret.carry();
        return ret;
    }
    
    public function add_(b : MPI) : MPI
    {
        if(b == null)return this;
        
        var i : uint;
        var minl : uint = Math.min(this.d.length, b.d.length);
        for(i = 0;i < minl;i++){
            this.d[i] += b.d[i];
        }
        for(i = this.d.length;i < b.d.length;i++){
            this.d.push(b.d[i]);
        }
        this.carry();
        return this;
    }
    
    public static function addint(a : MPI, b : Number) : MPI
    {
        var ret : MPI = new MPI(a);
        ret.d[0] += b;
        ret.carry();
        return ret;
    }
    
    public function addint_(b : Number) : MPI
    {
        this.d[0] += b;
        this.carry();
        return this;
    }
    
    public function inc() : void
    {
        d[0]++;
        for(var i : uint = 0;i < d.length;i++){
            if(d[i] >= 10){
                var c : uint = d[i] / 10;
                if(i == d.length - 1){
                    d.push(0);
                }
                d[i + 1] += c;
                d[i] %= 10;
            }else{
                break;
            }
        }
    }
    
    public static function sub(a: MPI, b : MPI) : MPI
    {
        if(b == null)return a;
        if(comp(a, b) < 0){
            var m : MPI = a; a = b; b = m;
        }
        
        // a - b
        var ret : MPI = new MPI(a);
        var i : uint;
        for(i = 0;i < b.d.length;i++){
            ret.d[i] -= b.d[i];
        }
        
        ret.minuscarry();
        return ret;
    }
    
    public function sub_(b : MPI) : MPI
    {
        if(b == null)return this;
        
        // a - b
        var i : uint;
        for(i = 0;i < b.d.length;i++){
            this.d[i] -= b.d[i];
        }
        
        this.minuscarry();
        return this;
    }
    
    public static function subint(a : MPI, b : Number) : MPI
    {
        var ret : MPI = new MPI(a);
        ret.d[0] -= b;
        ret.minuscarry();
        return ret;
    }
    
    public function subint_(b : Number) : MPI
    {
        this.d[0] -= b;
        this.minuscarry();
        return this;
    }
    
    public static function comp(a : MPI, b : MPI) : int
    {
        if(a.d.length != b.d.length)return a.d.length - b.d.length;
        for(var i : int = a.d.length - 1;i >= 0;i--){
            if(a.d[i] != b.d[i]){
                return a.d[i] - b.d[i];
            }
        }
        return 0;
    }
    
    public static function mul(a : MPI, b : MPI) : MPI
    {
        var retv : Vector.<uint> = new Vector.<uint>(a.d.length + b.d.length);
        var i : uint, j : uint;
        for(i = 0;i < retv.length;i++)retv[i] = 0;
        for(i = 0;i < a.d.length;i++){
            for(j = 0;j < b.d.length;j++){
                retv[i + j] += a.d[i] * b.d[j];
            }
        }
        var ret : MPI = new MPI(retv);
        ret.carry();
        ret.shave();
        return ret;
    }
    
    public static function mulint(a : MPI, b : int) : MPI
    {
        var ret : MPI = new MPI(a);
        for(var i : uint = 0;i < ret.d.length;i++){
            ret.d[i] *= b;
        }
        ret.carry();
        return ret;
    }
    
    public function mulint_(b : int) : MPI
    {
        for(var i : uint = 0;i < d.length;i++){
            this.d[i] *= b;
        }
        this.carry();
        return this;
    }
    
    public static function div(a : MPI, b : MPI) : Array
    {
        var bm : Array = new Array(10);
        var i : int, j : uint;
        for(i = 1;i <= 9;i++)bm[i] = mulint(b, i);
        
        var retr : Vector.<uint> = Vector.<uint>([]);
        var m : MPI = new MPI(0);
        for(i = a.d.length - 1;i >= 0;i--){
            m.mul10(1);
            m = addint(m, a.d[i]);
            for(j = 1;j <= 9;j++){
                if(comp(m, bm[j]) < 0)break;
            }
            retr.push(j - 1);
            if(j >= 2)m.sub_(bm[j - 1]);
        }
        
        var ret : MPI = new MPI(retr.reverse());
        ret.shave();
        m.shave();
        return [ret, m];
    }
    
    public static function divint(a : MPI, b : int) : Array
    {
        var i : int;
        
        var retr : Vector.<uint> = Vector.<uint>([]);
        var m : int = 0;
        for(i = a.d.length - 1;i >= 0;i--){
            m = m * 10 + a.d[i];
            var dv : int = m / b;
            retr.push(dv);
            m %= b;
        }
        
        var ret : MPI = new MPI(retr.reverse());
        ret.shave();
        return [ret, m];
    }
    
    // 末尾の0の個数
    public function nLast0() : int
    {
        for(var i : uint = 0;i < d.length && d[i] == 0;i++);
        return i;
    }
    
    // 位上げ
    public function carry() : void
    {
        for(var i : uint = 0;i < d.length;i++){
            if(d[i] >= 10){
                var c : uint = d[i] / 10;
                if(i == d.length - 1){
                    d.push(0);
                }
                d[i + 1] += c;
                d[i] %= 10;
            }
        }
    }
    
    // 位下げ
    public function minuscarry() : void
    {
        for(var i : uint = 0;i < d.length;i++){
            if(d[i] <= 2147483648)continue; // やっつけ
            var v : int = d[i] - 4294967296;
            var c : int = (-v + 9) / 10;
            d[i + 1] -= c;
            d[i] += 10 * c;
        }
        shave();
    }
    
    // 上位の桁から0を削る
    public function shave() : void
    {
        for(var i : uint = d.length - 1;i >= 1;i--){
            if(d[i] > 0)break;
            d.pop();
        }
    }
    
    public function mul10(x : int) : void
    {
        var r : Vector.<uint> = d.reverse();
        for(var i : uint = 0;i < x;i++){
            r.push(0);
        }
        d = r.reverse();
    }
    
    public function isZero() : Boolean
    {
        return d.length == 1 && d[0] == 0;
    }

    public function toNumber() : Number
    {
        return Number(toString());
    }

    public function toString() : String
    {
        return d.concat().reverse().join('');
    } 
}