forked from nayuki/Project-Euler-solutions
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathp066.java
More file actions
164 lines (129 loc) · 4.7 KB
/
p066.java
File metadata and controls
164 lines (129 loc) · 4.7 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
/*
* Solution to Project Euler problem 66
* Copyright (c) Project Nayuki. All rights reserved.
*
* https://www.nayuki.io/page/project-euler-solutions
* https://github.com/nayuki/Project-Euler-solutions
*/
import java.math.BigInteger;
import java.util.ArrayList;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
public final class p066 implements EulerSolution {
public static void main(String[] args) {
System.out.println(new p066().run());
}
/*
* Based on this insane theorem: Suppose D > 1 is an integer, non-perfect-square.
*
* Express sqrt(D) as the continued fraction (a0, a1, ..., a_{n-1}, (b0, b1, ..., b_{m-1})),
* where the sequence of b's is the periodic part.
*
* Let p/q (in lowest terms) = (a0, a1, ..., a_{n-1}, b0, b1, ..., b_{m-2}).
* (This is a truncation of the continued fraction with only one period minus the last term.)
*
* Then the minimum solution (x, y) for Pell's equation is given by:
* - (p, q) if m is even
* - (p^2 + D q^2, 2pq) if m is odd
*/
public String run() {
int minN = -1;
BigInteger maxX = BigInteger.ZERO;
for (int n = 2; n <= 1000; n++) {
if (Library.isSquare(n))
continue;
BigInteger x = smallestSolutionX(n);
if (x.compareTo(maxX) > 0) {
minN = n;
maxX = x;
}
}
return Integer.toString(minN);
}
// Returns the smallest x such that x > 0 and there exists some y such that x^2 - n y^2 = 1.
// Requires n to not be a perfect square.
private static BigInteger smallestSolutionX(int n) {
List<BigInteger>[] contFrac = sqrtToContinuedFraction(n);
List<BigInteger> temp = new ArrayList<>();
temp.addAll(contFrac[0]);
temp.addAll(contFrac[1].subList(0, contFrac[1].size() - 1));
Fraction val = new Fraction(temp.get(temp.size() - 1));
for (int i = temp.size() - 2; i >= 0; i--)
val = new Fraction(val.denominator, val.numerator).add(new Fraction(temp.get(i)));
if (contFrac[1].size() % 2 == 0)
return val.numerator;
else
return val.numerator.pow(2).add(val.denominator.pow(2).multiply(BigInteger.valueOf(n)));
}
// Returns the periodic continued fraction of sqrt(n). Requires n to not be a perfect square.
// result[0] is the minimal non-periodic prefix, and result[1] is the minimal periodic tail.
@SuppressWarnings("unchecked")
private static List<BigInteger>[] sqrtToContinuedFraction(int n) {
List<BigInteger> terms = new ArrayList<>();
Map<QuadraticSurd,Integer> seen = new HashMap<>();
QuadraticSurd val = new QuadraticSurd(BigInteger.ZERO, BigInteger.ONE, BigInteger.ONE, BigInteger.valueOf(n));
do {
seen.put(val, seen.size());
BigInteger flr = val.floor();
terms.add(flr);
val = val.subtract(new QuadraticSurd(flr, BigInteger.ZERO, BigInteger.ONE, val.d)).reciprocal();
} while (!seen.containsKey(val));
return new List[]{terms.subList(0, seen.get(val)), terms.subList(seen.get(val), terms.size())};
}
// Represents (a + b * sqrt(d)) / c. d must not be a perfect square.
private static final class QuadraticSurd {
public final BigInteger a, b, c, d;
public QuadraticSurd(BigInteger a, BigInteger b, BigInteger c, BigInteger d) {
if (c.signum() == 0)
throw new IllegalArgumentException();
// Simplify
if (c.signum() == -1) {
a = a.negate();
b = b.negate();
c = c.negate();
}
BigInteger gcd = a.gcd(b).gcd(c);
if (!gcd.equals(BigInteger.ONE)) {
a = a.divide(gcd);
b = b.divide(gcd);
c = c.divide(gcd);
}
this.a = a;
this.b = b;
this.c = c;
this.d = d;
}
public QuadraticSurd subtract(QuadraticSurd other) {
if (!d.equals(other.d))
throw new IllegalArgumentException();
return new QuadraticSurd(a.multiply(other.c).subtract(other.a.multiply(c)), b.multiply(other.c).subtract(other.b.multiply(c)), c.multiply(other.c), d);
}
public QuadraticSurd reciprocal() {
return new QuadraticSurd(a.multiply(c).negate(), b.multiply(c), b.multiply(b).multiply(d).subtract(a.multiply(a)), d);
}
public BigInteger floor() {
BigInteger temp = Library.sqrt(b.multiply(b).multiply(d));
if (b.signum() == -1)
temp = temp.add(BigInteger.ONE).negate();
temp = temp.add(a);
if (temp.signum() == -1)
temp = temp.subtract(c.subtract(BigInteger.ONE));
return temp.divide(c);
}
public boolean equals(Object obj) {
if (!(obj instanceof QuadraticSurd))
return false;
else {
QuadraticSurd other = (QuadraticSurd)obj;
return a.equals(other.a) && b.equals(other.b) && c.equals(other.c) && d.equals(other.d);
}
}
public int hashCode() {
return a.hashCode() + b.hashCode() + c.hashCode() + d.hashCode();
}
public String toString() {
return String.format("(%d + %d*sqrt(%d)) / %d", a, b, d, c);
}
}
}