114 lines
2.7 KiB
C++
114 lines
2.7 KiB
C++
/*
|
||
Consider quadratic Diophantine equations of the form:
|
||
|
||
x2 – Dy2 = 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.
|
||
|
||
-----------------------------------------------------------------------------------------------------
|
||
|
||
The value of x can be found by calculating the convergent of the continued fraction of D,
|
||
right before it repeats.
|
||
|
||
So if y/x is the corresponding approximation of sqrt(D) for the value of a_n of the continued fraction right before it repeats,
|
||
we immediately know x. So in order to solve this we just have to combine the previous two problem solution.
|
||
|
||
https://mathshistory.st-andrews.ac.uk/HistTopics/Pell/
|
||
*/
|
||
|
||
#include <bits/stdc++.h>
|
||
|
||
using namespace std;
|
||
|
||
int gcd(int a, int b){
|
||
if(b == 0){
|
||
return a;
|
||
}
|
||
return gcd(b, a % b);
|
||
}
|
||
|
||
int gcd(int a, int b, int c){
|
||
return gcd(a, gcd(b, c));
|
||
}
|
||
|
||
int p(int n, int p0, int p1){
|
||
if(n == 0){
|
||
return p0;
|
||
} else if(n == 1){
|
||
return p1;
|
||
} else if(n > 1){
|
||
return
|
||
}
|
||
|
||
return 0;
|
||
}
|
||
|
||
// Code yanked from problem 64
|
||
vector<int> findCF(int n){
|
||
float x = sqrt((float)n);
|
||
|
||
vector<int> cf;
|
||
|
||
if(x != sqrt(n)){
|
||
int a, b = 1, c = 1, d = 0;
|
||
int bn, cn, dn, g;
|
||
int b1, c1, d1;
|
||
|
||
for(int i = 0; ; ++i){
|
||
a = floor((floor(b * x) + d) / c);
|
||
cf.push_back(a);
|
||
|
||
bn = b*c;
|
||
cn = b*b*n - d*d - a*a*c*c + 2*a*c*d;
|
||
dn = a*c*c - c*d;
|
||
|
||
g = gcd(bn, cn, dn);
|
||
|
||
b = bn / g;
|
||
c = cn / g;
|
||
d = dn / g;
|
||
|
||
if(i == 0){
|
||
b1 = b;
|
||
c1 = c;
|
||
d1 = d;
|
||
} else if(b1 == b && c1 == c && d1 == d){
|
||
break;
|
||
}
|
||
}
|
||
}
|
||
|
||
return cf;
|
||
}
|
||
|
||
int main(){
|
||
|
||
cout << "Hello this is Patrick" << endl;
|
||
auto start = chrono::high_resolution_clock::now();
|
||
|
||
int maxX = 0;
|
||
for(int d = 1; d <= 1000; ++d){
|
||
auto cf = findCF(d);
|
||
|
||
// Insert problem 065 code once I have that fixed for big integers (cpp needs a dedicated internal library for that)
|
||
}
|
||
|
||
auto duration = chrono::duration_cast<chrono::milliseconds>(chrono::high_resolution_clock::now() - start);
|
||
cout << (float)duration.count()/1000 << endl;
|
||
return 0;
|
||
}
|