112 lines
2.8 KiB
C
112 lines
2.8 KiB
C
/* All square roots are periodic when written as continued fractions and can be written in the form:
|
||
*
|
||
* √N=a0+1/(a1+1/(a2+1/(a3+…
|
||
*
|
||
* For example, let us consider √23:
|
||
*
|
||
* √23=4+√23−4=4+1/(1/(√23−4))=4+1/(1+(√23−3)/7)
|
||
*
|
||
* If we continue we would get the following expansion:
|
||
*
|
||
* √23=4+1/(1+1/(3+1/(1+1/(8+…
|
||
*
|
||
* The process can be summarised as follows:
|
||
* a0=4,1/(√23−4)=(√23+4)/7=1+(√23−3)/7
|
||
* a1=1,7/(√23−3)=7(√23+3)/14=3+(√23−3)/2
|
||
* a2=3,2/(√23−3)=2(√23+3)/14=1+(√23−4)/7
|
||
* a3=1,7/(√23−4)=7(√23+4)/7=8+√23−4
|
||
* a4=8,1/(√23−4)=(√23+4)/7=1+(√23−3)/7
|
||
* a5=1,7/(√23−3)=7(√23+3)/14=3+(√23−3)/2
|
||
* a6=3,2/(√23−3)=2(√23+3)/14=1+(√23−4)/7
|
||
* a7=1,7/(√23−4)=7(√23+4)/7=8+√23−4
|
||
*
|
||
* It can be seen that the sequence is repeating. For conciseness, we use the notation √23=[4;(1,3,1,8)], to indicate that the block (1,3,1,8)
|
||
* repeats indefinitely.
|
||
*
|
||
* The first ten continued fraction representations of (irrational) square roots are:
|
||
*
|
||
* √2=[1;(2)], period=1
|
||
* √3=[1;(1,2)], period=2
|
||
* √5=[2;(4)], period=1
|
||
* √6=[2;(2,4)], period=2
|
||
* √7=[2;(1,1,1,4)], period=4
|
||
* √8=[2;(1,4)], period=2
|
||
* √10=[3;(6)], period=1
|
||
* √11=[3;(3,6)], period=2
|
||
* √12=[3;(2,6)], period=2
|
||
* √13=[3;(1,1,1,1,6)], period=5
|
||
*
|
||
* Exactly four continued fractions, for N≤13, have an odd period.
|
||
*
|
||
* How many continued fractions for N≤10000 have an odd period?*/
|
||
|
||
#define _POSIX_C_SOURCE 199309L
|
||
|
||
#include <stdio.h>
|
||
#include <stdlib.h>
|
||
#include <math.h>
|
||
#include <time.h>
|
||
#include "projecteuler.h"
|
||
|
||
int is_square(int n);
|
||
|
||
int main(int argc, char **argv)
|
||
{
|
||
int i, count = 0, period;
|
||
int *fraction;
|
||
double elapsed;
|
||
struct timespec start, end;
|
||
|
||
clock_gettime(CLOCK_MONOTONIC, &start);
|
||
|
||
for(i = 2; i <= 10000; i++)
|
||
{
|
||
/* Perfect squares are obviously not represented as continued fractions.
|
||
* For all other numbers, calculate their period and check if it's odd.*/
|
||
if(!is_square(i))
|
||
{
|
||
if((fraction = build_sqrt_cont_fraction(i, &period, 300)) == NULL)
|
||
{
|
||
fprintf(stderr, "Error! Build_cont_fraction function returned NULL\n");
|
||
return 1;
|
||
}
|
||
|
||
if(period % 2 != 0)
|
||
{
|
||
count++;
|
||
}
|
||
}
|
||
}
|
||
|
||
free(fraction);
|
||
|
||
clock_gettime(CLOCK_MONOTONIC, &end);
|
||
|
||
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
|
||
|
||
printf("Project Euler, Problem 64\n");
|
||
printf("Answer: %d\n", count);
|
||
|
||
printf("Elapsed time: %.9lf seconds\n", elapsed);
|
||
|
||
return 0;
|
||
}
|
||
|
||
int is_square(int n)
|
||
{
|
||
int m;
|
||
double p;
|
||
|
||
p = sqrt(n);
|
||
m = p;
|
||
|
||
if(p == m)
|
||
{
|
||
return 1;
|
||
}
|
||
else
|
||
{
|
||
return 0;
|
||
}
|
||
}
|