project euler 66

Problem 66

---

Diophantine equation

Consider quadratic Diophantine equations of the form:

x 2 – Dy 2 = 1

For example, when D=13, the minimal solution in x is 6492 – 13×1802 = 1.

It can be assumed that there are no solutions in positive integers when D is square.

By finding minimal solutions in x for D = {2, 3, 5, 6, 7}, we obtain the following:

32 – 2×22 = 1
22 – 3×12 = 1
92 – 5×42 = 1
52 – 6×22 = 1
82 – 7×32 = 1

Hence, by considering minimal solutions in x for D ≤ 7, the largest x is obtained when D=5.

Find the value of D ≤ 1000 in minimal solutions of x for which the largest value of x is obtained.

---

丢番图方程
考虑如下形式的二次丢番图方程：

x 2 – Dy 2 = 1

举例而言，当D=13时，x的最小值出现在6492 – 13×1802 = 1。

可以断定，当D是平方数时，这个方程不存在正整数解。

对于D= {2, 3, 5, 6, 7}分别求出x取最小值的解，我们得到：

32 – 2×22 = 1
22 – 3×12 = 1
92 – 5×42 = 1
52 – 6×22 = 1
82 – 7×32 = 1

因此，对于所有D ≤ 7，当D=5时x的最小值最大。

对于D ≤ 1000，求使得x的最小值最大的D值。

package cuppics;

import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;

import junit.framework.TestCase;

public class Prj66 extends TestCase {

	public static final int NUM_OF_CHOOSE = 1000;
	public static final int D_LIMIT = 1000;

	
	/**
	 * http://mathworld.wolfram.com/PellEquation.html
	 */
	public void testDiophantineEquation() {

		int maxId = 0;
		String maxStr = "";

		for (int i = 2; i <= D_LIMIT; i++) {
			if (isSquare(i)) {
				continue;
			}
			String[] pq = solvePell(getSqure(i), getCircle(i));
			if (compStr(pq[0], maxStr)) {
				maxStr = pq[0];
				maxId = i;
			}

			String fstr = "(x,y)=(%s,%s)";
			fstr = String.format(fstr, pq[0], pq[1]);
			System.out.println(fstr);

		}
		
		String fstr = "D=%d, x =%s";
		fstr = String.format(fstr, maxId, maxStr);
		System.out.println(fstr);

	}

	boolean compStr(String lf, String rt) {
		if (lf.length() > rt.length()) {
			return true;
		}

		if (lf.length() == rt.length()) {
			if (lf.compareTo(rt) >= 0) {
				return true;
			}
			return false;
		}

		return false;

	}

	public String[] solvePell(int a0, List<Integer> cycle) {
		int sz = cycle.size();

		int r = sz - 1;
		if (r % 2 == 0) {
			r = 2 * r + 1;
		}

		r = r + 1;

		return iterArr(r, a0, cycle);
	}

	int getSqure(int val) {
		return (int) Math.sqrt(val);
	}

	boolean isSquare(int val) {
		int _val = (int) Math.sqrt(val);
		return _val * _val == val;
	}

	public int getBits() {
		return NUM_OF_CHOOSE * 3 + 1;
	}

	public int[] power(int[] arr, int multiplier, int b) {

		int[] ret = Arrays.copyOf(arr, arr.length);
		int count = b;
		while (count > 1) {

			ret = multiply(ret, multiplier);
			count--;
		}

		return ret;
	}

	public int[] multiply(int[] arr, int b) {

		int[] ret = Arrays.copyOf(arr, arr.length);

		for (int i = 0; i < arr.length; i++) {
			ret[i] *= b;
		}

		for (int i = arr.length - 1; i > 0; i--) {

			ret[i - 1] = ret[i - 1] + ret[i] / 10;
			ret[i] = ret[i] % 10;
		}

		return ret;
	}

	public int[] add(int[] a, int[] b) {
		assert (a.length == b.length);

		int[] ret = new int[a.length];

		for (int i = 0; i < a.length; i++) {
			ret[i] = a[i] + b[i];
		}

		for (int i = a.length - 1; i > 0; i--) {
			ret[i - 1] = ret[i - 1] + ret[i] / 10;
			ret[i] %= 10;
		}

		return ret;
	}

	public int[] int2Arr(int val) {

		String str = Integer.toString(val, 10);
		int[] ret = new int[str.length()];
		for (int i = 0; i < ret.length; i++) {
			ret[i] = Integer.parseInt(String.valueOf(str.charAt(i)));
		}
		return ret;
	}

	public int[] enLarge(int[] val, int bits) {

		assert (val.length <= bits);

		int[] ret = new int[bits];

		for (int i = bits - val.length, j = 0; i < bits; i++, j++) {
			ret[i] = val[j];
		}
		return ret;

	}

	public class CycleCounter {

		private int count = 0;

		private List<Integer> data;

		public CycleCounter(List<Integer> bits) {
			this.data = bits;
			this.count = 0;
		}

		public int getNext() {

			++count;
			if (count >= data.size()) {
				count = 0;
			}

			return data.get(count);
		}

	}

	/**
	 * [a0;a1,a2....] a0/1, (a1a0 + 1)/a1,
	 * 
	 * @param n
	 * @return
	 */
	public String[] iterArr(int n, int init_a0, List<Integer> cycle) {

		int bits = getBits();

		// start
		int a0 = init_a0;
		int a1 = cycle.get(0);
		// int an = 2;
		int p0 = a0;
		int p1 = a1 * a0 + 1;
		int q0 = 1;
		int q1 = a1;

		int[] pn_1 = enLarge(int2Arr(p1), bits);
		int[] pn_2 = enLarge(int2Arr(p0), bits);
		int[] qn_1 = enLarge(int2Arr(q1), bits);
		int[] qn_2 = enLarge(int2Arr(q0), bits);

		if (n == 1) {
			return new String[] { Integer.toString(a0), Integer.toString(1) };
		} else if (n == 2) {
			return new String[] { Integer.toString(p1), Integer.toString(q1) };
		}

		assert (n >= 2);
		int[] pn = enLarge(int2Arr(0), bits);
		int[] qn = enLarge(int2Arr(0), bits);

		CycleCounter cc = new CycleCounter(cycle);

		for (int i = 3; i <= n; i++) {
			int ccc = cc.getNext();
			pn = add(multiply(pn_1, ccc), pn_2);
			qn = add(multiply(qn_1, ccc), qn_2);

			// String pnStr = intArr2Str(pn);
			// String qnStr = intArr2Str(qn);
			// System.out.println(pnStr + "," + qnStr);
			pn_2 = Arrays.copyOf(pn_1, bits);
			pn_1 = Arrays.copyOf(pn, bits);
			qn_2 = Arrays.copyOf(qn_1, bits);
			qn_1 = Arrays.copyOf(qn, bits);
		}

		return new String[] { intArr2Str(pn), intArr2Str(qn) };
	}

	private String intArr2Str(int[] pn) {
		StringBuilder sb = new StringBuilder();

		int i = 0;
		while (pn[i] == 0) {
			i++;
		}
		for (; i < pn.length; i++) {
			sb.append(pn[i]);
		}
		return sb.toString();

	}

	/**
	 * https://en.wikipedia.org/wiki/Methods_of_computing_square_roots#
	 * Continued_fraction_expansion
	 */
	public List<Integer> getCircle(int S) {

		int val = (int) Math.sqrt(S);
		// square
		if (val * val == S) {
			return new ArrayList<Integer>();
		}

		List<Integer> ret = new ArrayList<Integer>();

		int m0 = 0;
		int d0 = 1;
		int a0 = (int) Math.floor(Math.sqrt(S));
		int _a0 = a0;

		int mn = 0;
		int dn = 0;
		int an = 0;
		while (true) {
			mn = d0 * a0 - m0;
			dn = (S - mn * mn) / d0;
			an = (int) Math.floor((_a0 + mn) / dn);

			m0 = mn;
			d0 = dn;
			a0 = an;
			ret.add(a0);
			if (an == 2 * _a0) {
				break;
			}
		}
		return ret;
	}

}