Add more solutions and minor corrections

Added solutions for problem 71, 72, 73, 74 and 75 in C and python
and made a few corrections in the code for other problems.
This commit is contained in:
daniele 2019-09-29 19:38:39 +02:00
parent 2cdbca922a
commit 7c055b749e
Signed by: fuxino
GPG Key ID: 6FE25B4A3EE16FDA
21 changed files with 962 additions and 6 deletions

1
.gitignore vendored
View File

@ -1 +1,2 @@
Python/__pycache__
ProblemsOverviews/*

View File

@ -27,7 +27,7 @@ int main(int argc, char **argv)
{
for(n = 1; n < m && !found; n++)
{
if(gcd(m, n) == 1 && (m % 2 == 0 && n % 2 != 0 || m % 2 != 0 && n % 2 == 0))
if(gcd(m, n) == 1 && ((m % 2 == 0 && n % 2 != 0) || (m % 2 != 0 && n % 2 == 0)))
{
a = m * m - n * n;
b = 2 * m * n;

View File

@ -77,7 +77,7 @@ int main(int argc, char **argv)
/* Use the function implemented in projecteuler.c to find the maximum path.*/
max = find_max_path(triang, 15);
for(i = 0; i < 100; i++)
for(i = 0; i < 15; i++)
{
free(triang[i]);
}

View File

@ -49,7 +49,7 @@ int count(int value, int n, int i)
if(value == 200)
{
return n+1;
return n + 1;
}
else if(value > 200)
{

View File

@ -85,7 +85,7 @@ int count_digits(long int a)
return count;
}
int is_permutation(int a, int b)
int is_permutation(long int a, long int b)
{
int i;
int digits1[10] = {0}, digits2[10] = {0};

View File

@ -27,7 +27,11 @@ int main(int argc, char **argv)
clock_gettime(CLOCK_MONOTONIC, &start);
primes = sieve(N);
if((primes = sieve(N)) == NULL)
{
fprintf(stderr, "Error! Sieve function returned NULL\n");
return 1;
}
for(i = 2; i < N; i++)
{

64
C/p071.c Normal file
View File

@ -0,0 +1,64 @@
/* Consider the fraction, n/d, where n and d are positive integers. If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
*
* If we list the set of reduced proper fractions for d 8 in ascending order of size, we get:
*
* 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
*
* It can be seen that 2/5 is the fraction immediately to the left of 3/7.
*
* By listing the set of reduced proper fractions for d 1,000,000 in ascending order of size, find the numerator
* of the fraction immediately to the left of 3/7.*/
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include "projecteuler.h"
#define N 1000000
int main(int argc, char **argv)
{
int i, j, n, d, max_n;
double elapsed, max = 0.0;
struct timespec start, end;
clock_gettime(CLOCK_MONOTONIC, &start);
/* For each denominator q, we need to find the biggest numerator p for which
* p/q<a/b, where a/b is 3/7 for this problem. So:
* pb<aq
* pb<=aq-1
* p<=(aq-1)/b
* So we can set p=(3*q-1)/7 (using integer division).*/
for(i = 2; i <= N; i++)
{
j = (3 * i - 1) / 7;
if((double)j / i > max)
{
n = j;
d = i;
max = (double)n / d;
/* Reduce the fraction if it's not already reduced.*/
if(gcd(i, j) > 1)
{
n /= gcd(i, j);
d /= gcd(i, j);
}
max_n = n;
}
}
clock_gettime(CLOCK_MONOTONIC, &end);
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
printf("Project Euler, Problem 71\n");
printf("Answer: %d\n", max_n);
printf("Elapsed time: %.9lf seconds\n", elapsed);
return 0;
}

54
C/p072.c Normal file
View File

@ -0,0 +1,54 @@
/* Consider the fraction, n/d, where n and d are positive integers. If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
*
* If we list the set of reduced proper fractions for d 8 in ascending order of size, we get:
*
* 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
*
* It can be seen that there are 21 elements in this set.
*
* How many elements would be contained in the set of reduced proper fractions for d 1,000,000?*/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include "projecteuler.h"
#define N 1000001
int main(int argc, char **argv)
{
int i;
int *primes;
long int count = 0;
double elapsed;
struct timespec start, end;
clock_gettime(CLOCK_MONOTONIC, &start);
if((primes = sieve(N)) == NULL)
{
fprintf(stderr, "Error! Sieve function returned NULL\n");
return 1;
}
/* For any denominator d, the number of reduced proper fractions is
* the number of fractions n/d where gcd(n, d)=1, which is the definition
* of Euler's Totient Function phi. It's sufficient to calculate phi for each
* denominator and sum the value.*/
for(i = 2; i < N; i++)
{
count += phi(i, primes);
}
clock_gettime(CLOCK_MONOTONIC, &end);
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
printf("Project Euler, Problem 72\n");
printf("Answer: %ld\n", count);
printf("Elapsed time: %.9lf seconds\n", elapsed);
return 0;
}

51
C/p073.c Normal file
View File

@ -0,0 +1,51 @@
/* Consider the fraction, n/d, where n and d are positive integers. If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
*
* If we list the set of reduced proper fractions for d 8 in ascending order of size, we get:
*
* 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
*
* It can be seen that there are 3 fractions between 1/3 and 1/2.
*
* How many fractions lie between 1/3 and 1/2 in the sorted set of reduced proper fractions for d 12,000?*/
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include "projecteuler.h"
int main(int argc, char **argv)
{
int i, j, limit, count = 0;
double elapsed;
struct timespec start, end;
clock_gettime(CLOCK_MONOTONIC, &start);
/* For each denominator q, we need to find the fractions p/q for which
* 1/3<p/q<1/2. For the lower limit, if p=q/3, then p/q=1/3, so we take
* p=q/3+1. For the upper limit, if p=q/2 p/q=1/2, so we take p=(q-1)/2.*/
for(i = 2; i <= 12000; i++)
{
limit = (i - 1) / 2;
for(j = i / 3 + 1; j <= limit; j++)
{
/* Increment the counter if the current fraction is reduced.*/
if(gcd(j, i) == 1)
{
count++;
}
}
}
clock_gettime(CLOCK_MONOTONIC, &end);
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
printf("Project Euler, Problem 73\n");
printf("Answer: %d\n", count);
printf("Elapsed time: %.9lf seconds\n", elapsed);
return 0;
}

132
C/p074.c Normal file
View File

@ -0,0 +1,132 @@
/* The number 145 is well known for the property that the sum of the factorial of its digits is equal to 145:
*
* 1! + 4! + 5! = 1 + 24 + 120 = 145
*
* Perhaps less well known is 169, in that it produces the longest chain of numbers that link back to 169; it turns out that there are
* only three such loops that exist:
*
* 169 363601 1454 169
* 871 45361 871
* 872 45362 872
*
* It is not difficult to prove that EVERY starting number will eventually get stuck in a loop. For example,
*
* 69 363600 1454 169 363601 ( 1454)
* 78 45360 871 45361 ( 871)
* 540 145 ( 145)
*
* Starting with 69 produces a chain of five non-repeating terms, but the longest non-repeating chain with a starting number
* below one million is sixty terms.
*
* How many chains, with a starting number below one million, contain exactly sixty non-repeating terms?*/
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#define N 1000000
int factorial(int n);
int len_chain(int n);
int factorials[10] = {0};
int chains[N] = {0};
int main(int argc, char **argv)
{
int i, count = 0;
double elapsed;
struct timespec start, end;
clock_gettime(CLOCK_MONOTONIC, &start);
/* Simple brute force approach, for every number calculate
* the length of the chain.*/
for(i = 3; i < N; i++)
{
if(len_chain(i) == 60)
{
count++;
}
}
clock_gettime(CLOCK_MONOTONIC, &end);
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
printf("Project Euler, Problem 74\n");
printf("Answer: %d\n", count);
printf("Elapsed time: %.9lf seconds\n", elapsed);
return 0;
}
/* Recursively calculate the factorial, of a number, saving the result
* so that it can be reused later.*/
int factorial(int n)
{
if(n == 0 || n == 1)
{
return 1;
}
else if(factorials[n] != 0)
{
return factorials[n];
}
else
{
factorials[n]=n*factorial(n-1);
return factorials[n];
}
}
int len_chain(int n)
{
int i, count = 0, finished = 0, value, tmp;
int chain[60];
value = n;
chain[count] = value;
while(!finished)
{
tmp = 0;
count++;
/* Generate the next number of the chain by taking
* the digits of the current value, calculating the
* factorials and adding them.*/
while(value != 0)
{
tmp += factorial(value % 10);
value /= 10;
}
/* If the chain length for the new value has already been
* calculated before, use the saved value (only chains for
* values smaller than N are saved).*/
if(tmp < N && chains[tmp] != 0)
{
return count + chains[tmp];
}
value = tmp;
/* If the current value is already present in the chain,
* the chain is finished.*/
for(i = 0; i < count; i++)
{
if(chain[i] == value)
{
finished = 1;
}
}
chain[count] = value;
}
chains[n] = count;
return count;
}

102
C/p075.c Normal file
View File

@ -0,0 +1,102 @@
/* It turns out that 12 cm is the smallest length of wire that can be bent to form an integer sided right angle triangle in exactly one way,
* but there are many more examples.
*
* 12 cm: (3,4,5)
* 24 cm: (6,8,10)
* 30 cm: (5,12,13)
* 36 cm: (9,12,15)
* 40 cm: (8,15,17)
* 48 cm: (12,16,20)
*
* In contrast, some lengths of wire, like 20 cm, cannot be bent to form an integer sided right angle triangle, and other lengths allow
* more than one solution to be found; for example, using 120 cm it is possible to form exactly three different integer sided right angle triangles.
*
* 120 cm: (30,40,50), (20,48,52), (24,45,51)
*
* Given that L is the length of the wire, for how many values of L 1,500,000 can exactly one integer sided right angle triangle be formed?*/
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include "projecteuler.h"
#define N 1500000
int main(int argc, char **argv)
{
int i, a, b, c, tmpa, tmpb, tmpc, m, n, count = 0;
int *l;
double elapsed;
struct timespec start, end;
clock_gettime(CLOCK_MONOTONIC, &start);
if((l = (int *)calloc(N + 1, sizeof(int))) == NULL)
{
fprintf(stderr, "Error while allocating memory\n");
return 1;
}
/* Generate all Pythagorean triplets using Euclid's algorithm:
* For m>=2 and n<m:
* a=m*m-n*n
* b=2*m*n
* c=m*m+n*n
* This gives a primitive triple if gcd(m, n)=1 and exactly one
* of m and n is odd. To generate all the triples, generate all
* the primitive one and multiply them by i=2,3, ..., n until the
* perimeter is larger than the limit. The limit for m is 865, because
* when m=866 even with the smaller n (i.e. 1) the perimeter is greater
* than the given limit.*/
for(m = 2; m < 866; m++)
{
for(n = 1; n < m; n++)
{
if(gcd(m, n) == 1 && ((m % 2 == 0 && n % 2 != 0) || (m % 2 != 0 && n % 2 == 0)))
{
a = m * m - n * n;
b = 2 * m * n;
c = m * m + n * n;
if(a + b + c <= N)
{
l[a+b+c]++;
}
i = 2;
tmpa = i * a;
tmpb = i * b;
tmpc = i * c;
while(tmpa + tmpb + tmpc <= N)
{
l[tmpa+tmpb+tmpc]++;
i++;
tmpa = i * a;
tmpb = i * b;
tmpc = i * c;
}
}
}
}
for(i = 0; i <= N; i++)
{
if(l[i] == 1)
count++;
}
free(l);
clock_gettime(CLOCK_MONOTONIC, &end);
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
printf("Project Euler, Problem 75\n");
printf("Answer: %d\n", count);
printf("Elapsed time: %.9lf seconds\n", elapsed);
return 0;
}

69
C/p076.c Normal file
View File

@ -0,0 +1,69 @@
/* It is possible to write five as a sum in exactly six different ways:
*
* 4 + 1
* 3 + 2
* 3 + 1 + 1
* 2 + 2 + 1
* 2 + 1 + 1 + 1
* 1 + 1 + 1 + 1 + 1
How many different ways can one hundred be written as a sum of at least two positive integers?*/
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
int count(int value, int n, int i);
int integers[99];
int main(int argc, char **argv)
{
int i, n;
double elapsed;
struct timespec start, end;
clock_gettime(CLOCK_MONOTONIC, &start);
for(i = 0; i < 99; i++)
integers[i] = i + 1;
n = count(0, 0, 0);
clock_gettime(CLOCK_MONOTONIC, &end);
elapsed = (end.tv_sec - start.tv_sec) + (double)(end.tv_nsec - start.tv_nsec) / 1000000000;
printf("Project Euler, Problem 76\n");
printf("Answer: %d\n", n);
printf("Elapsed time: %.9lf seconds\n", elapsed);
return 0;
}
int count(int value, int n, int i)
{
int j;
for(j = i; j < 99; j++)
{
value += integers[j];
if(value == 100)
{
return n + 1;
}
else if(value > 100)
{
return n;
}
else
{
n = count(value, n, j);
value -= integers[j];
}
}
return n;
}

View File

@ -690,3 +690,96 @@ int phi_semiprime(int n, int p, int q)
return n - (p + q) + 1;
}
}
/* Function to calculate phi(n) for any n. If n is prime, phi(n)=n-1,
* is n is semiprime, use the above function. In any other case,
* phi(n)=n*prod(1-1/p) for every distinct prime p that divides n.*/
int phi(int n, int *primes)
{
int i, p, q, limit;
double ph = (double)n;
/* If n is prime, phi(n)=n-1*/
if(primes[n])
{
return n - 1;
}
/* If n is semiprime, use above function.*/
if(is_semiprime(n, &p, &q, primes))
{
return phi_semiprime(n, p, q);
}
/* If 2 is a factor of n, multiply the current ph (which now is n)
* by 1-1/2, then divide all factors 2.*/
if(n % 2 == 0)
{
ph *= (1.0 - 1.0 / 2.0);
do
{
n /= 2;
}while(n % 2 == 0);
}
/* If 3 is a factor of n, multiply the current ph by 1-1/3,
* then divide all factors 3.*/
if(n % 3 == 0)
{
ph *= (1.0 - 1.0 / 3.0);
do
{
n /= 3;
}while(n % 3 == 0);
}
/*Any number can have only one prime factor greater than its
* square root, so we can stop checking at this point and deal
* with the only factor larger than sqrt(n), if present, at the end.*/
limit = floor(sqrt(n));
/* Every prime other than 2 and 3 is in the form 6k+1 or 6k-1.
* If I check all those value no prime factors of the number
* will be missed. For each of these possible primes, check if
* they are prime, then check if the number divides n, in which
* case update the current ph.*/
for(i = 5; i <= limit; i += 6)
{
if(primes[i])
{
if(n % i == 0)
{
ph *= (1.0 - 1.0 / i);
do
{
n /= i;
}while(n % i == 0);
}
}
if(primes[i+2])
{
if(n % (i + 2) == 0)
{
ph *= (1.0 - 1.0 /(i + 2));
do
{
n /= (i + 2);
}while(n % (i + 2) == 0);
}
}
}
/* After dividing all prime factors smaller than sqrt(n), n is either 1
* or is equal to the only prime factor greater than sqrt(n). In this
* second case, we need to update ph with the last prime factor.*/
if(n > 1)
{
ph *= (1.0 - 1.0 / n);
}
return (int)ph;
}

View File

@ -23,5 +23,6 @@ int *build_sqrt_cont_fraction(int i, int *period, int l);
int pell_eq(int i, mpz_t x);
int is_semiprime(int n, int *p, int *q, int *primes);
int phi_semiprime(int n, int p, int q);
int phi(int n, int *primes);
#endif

View File

@ -10,8 +10,9 @@
#
# If the product of these four fractions is given in its lowest common terms, find the value of the denominator.
from math import gcd
from timeit import default_timer
from projecteuler import gcd
def main():
start = default_timer()

54
Python/p071.py Normal file
View File

@ -0,0 +1,54 @@
#!/usr/bin/python3
# Consider the fraction, n/d, where n and d are positive integers. If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
#
# If we list the set of reduced proper fractions for d ≤ 8 in ascending order of size, we get:
#
# 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
#
# It can be seen that 2/5 is the fraction immediately to the left of 3/7.
#
# By listing the set of reduced proper fractions for d ≤ 1,000,000 in ascending order of size, find the numerator
# of the fraction immediately to the left of 3/7.
from math import gcd
from timeit import default_timer
def main():
start = default_timer()
N = 1000000
max_ = 0
# For each denominator q, we need to find the biggest numerator p for which
# p/q<a/b, where a/b is 3/7 for this problem. So:
# pb<aq
# pb<=aq-1
# p<=(aq-1)/b
# So we can set p=(3*q-1)/7 (using integer division).
for i in range(2, N+1):
j = (3 * i - 1) // 7
if j / i > max_:
n = j
d = i
max_ = n / d
# Reduce the fractions if it's not already reduced.
if gcd(i, j) > 1:
n /= gcd(i, j)
d /= gcd(i, j)
max_n = n
end = default_timer()
print('Project Euler, Problem 71')
print('Answer: {}'.format(max_n))
print('Elapsed time: {:.9f} seconds'.format(end - start))
if __name__ == '__main__':
main()

39
Python/p072.py Normal file
View File

@ -0,0 +1,39 @@
#!/usr/bin/python3
# Consider the fraction, n/d, where n and d are positive integers. If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
#
# If we list the set of reduced proper fractions for d ≤ 8 in ascending order of size, we get:
#
# 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
#
# It can be seen that there are 21 elements in this set.
#
# How many elements would be contained in the set of reduced proper fractions for d ≤ 1,000,000?
from timeit import default_timer
from projecteuler import sieve, phi
def main():
start = default_timer()
N = 1000001
count = 0
primes = sieve(N)
# For any denominator d, the number of reduced proper fractions is
# the number of fractions n/d where gcd(n, d)=1, which is the definition
# of Euler's Totient Function phi. It's sufficient to calculate phi for each
# denominator and sum the value.
for i in range(2, N):
count = count + phi(i, primes)
end = default_timer()
print('Project Euler, Problem 72')
print('Answer: {}'.format(count))
print('Elapsed time: {:.9f} seconds'.format(end - start))
if __name__ == '__main__':
main()

41
Python/p073.py Normal file
View File

@ -0,0 +1,41 @@
#!/usr/bin/python3
# Consider the fraction, n/d, where n and d are positive integers. If n<d and HCF(n,d)=1, it is called a reduced proper fraction.
#
# If we list the set of reduced proper fractions for d ≤ 8 in ascending order of size, we get:
#
# 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
#
# It can be seen that there are 3 fractions between 1/3 and 1/2.
#
# How many fractions lie between 1/3 and 1/2 in the sorted set of reduced proper fractions for d ≤ 12,000?
from math import gcd
from timeit import default_timer
def main():
start = default_timer()
count = 0
# For each denominator q, we need to find the fractions p/q for which
# 1/3<p/q<1/2. For the lower limit, if p=q/3, then p/q=1/3, so we take
# p=q/3+1. For the upper limit, if p=q/2 p/q=1/2, so we take p=(q-1)/2.
for i in range(2, 12001):
limit = (i - 1) // 2 + 1
for j in range(i//3+1, limit):
# Increment the counter if the current fraction is reduced.*/
if gcd(j, i) == 1:
count = count + 1
end = default_timer()
print('Project Euler, Problem 73')
print('Answer: {}'.format(count))
print('Elapsed time: {:.9f} seconds'.format(end - start))
if __name__ == '__main__':
main()

97
Python/p074.py Normal file
View File

@ -0,0 +1,97 @@
#!/usr/bin/python3
# The number 145 is well known for the property that the sum of the factorial of its digits is equal to 145:
#
# 1! + 4! + 5! = 1 + 24 + 120 = 145
#
# Perhaps less well known is 169, in that it produces the longest chain of numbers that link back to 169; it turns out that there are
# only three such loops that exist:
#
# 169 → 363601 → 1454 → 169
# 871 → 45361 → 871
# 872 → 45362 → 872
#
# It is not difficult to prove that EVERY starting number will eventually get stuck in a loop. For example,
#
# 69 → 363600 → 1454 → 169 → 363601 (→ 1454)
# 78 → 45360 → 871 → 45361 (→ 871)
# 540 → 145 (→ 145)
#
# Starting with 69 produces a chain of five non-repeating terms, but the longest non-repeating chain with a starting number
# below one million is sixty terms.
#
# How many chains, with a starting number below one million, contain exactly sixty non-repeating terms?
from math import factorial
from timeit import default_timer
def len_chain(n):
global N
global chains
chain = [0] * 60
count = 0
finished = 0
value = n
chain[count] = value
while not finished:
tmp = 0
count = count + 1
# Generate the next number of the chain by taking
# the digits of the current value, calculating the
# factorials and adding them.*/
while value != 0:
tmp = tmp + factorial(value % 10)
value = value // 10
# If the chain length for the new value has already been
# calculated before, use the saved value (only chains for
# values smaller than N are saved).*/
if tmp < N and chains[tmp] != 0:
return count + chains[tmp]
value = tmp
# If the current value is already present in the chain,
# the chain is finished.*/
for i in range(count):
if chain[i] == value:
finished = 1
chain[count] = value
chains[n] = count
return count
def main():
start = default_timer()
global N
global chains
N = 1000000
chains = [0] * N
count = 0
# Simple brute force approach, for every number calculate
# the length of the chain.
for i in range(3, N):
if len_chain(i) == 60:
count = count + 1
end = default_timer()
print('Project Euler, Problem 74')
print('Answer: {}'.format(count))
print('Elapsed time: {:.9f} seconds'.format(end - start))
if __name__ == '__main__':
main()

80
Python/p075.py Normal file
View File

@ -0,0 +1,80 @@
#!/usr/bin/python3
# It turns out that 12 cm is the smallest length of wire that can be bent to form an integer sided right angle triangle in exactly one way,
# but there are many more examples.
#
# 12 cm: (3,4,5)
# 24 cm: (6,8,10)
# 30 cm: (5,12,13)
# 36 cm: (9,12,15)
# 40 cm: (8,15,17)
# 48 cm: (12,16,20)
#
# In contrast, some lengths of wire, like 20 cm, cannot be bent to form an integer sided right angle triangle, and other lengths allow
# more than one solution to be found; for example, using 120 cm it is possible to form exactly three different integer sided right angle triangles.
#
# 120 cm: (30,40,50), (20,48,52), (24,45,51)
#
# Given that L is the length of the wire, for how many values of L ≤ 1,500,000 can exactly one integer sided right angle triangle be formed?
from math import gcd
from timeit import default_timer
def main():
start = default_timer()
N = 1500000
l = [0] * (N+1)
# Generate all Pythagorean triplets using Euclid's algorithm:
# For m>=2 and n<m:
# a=m*m-n*n
# b=2*m*n
# c=m*m+n*n
# This gives a primitive triple if gcd(m, n)=1 and exactly one
# of m and n is odd. To generate all the triples, generate all
# the primitive one and multiply them by i=2,3, ..., n until the
# perimeter is larger than the limit. The limit for m is 865, because
# when m=866 even with the smaller n (i.e. 1) the perimeter is greater
# than the given limit.
for m in range(2, 866):
for n in range(1, m):
if gcd(m, n) == 1 and ((m % 2 == 0 and n % 2 != 0) or (m % 2 != 0 and n % 2 == 0)):
a = m * m - n * n
b = 2 * m * n
c = m * m + n * n
if(a + b + c <= N):
l[a+b+c] = l[a+b+c] + 1
i = 2
tmpa = i * a
tmpb = i * b
tmpc = i * c
while(tmpa + tmpb + tmpc <= N):
l[tmpa+tmpb+tmpc] = l[tmpa+tmpb+tmpc] + 1
i = i + 1
tmpa = i * a
tmpb = i * b
tmpc = i * c
count = 0
for i in range(N+1):
if l[i] == 1:
count = count + 1
end = default_timer()
print('Project Euler, Problem 75')
print('Answer: {}'.format(count))
print('Elapsed time: {:.9f} seconds'.format(end - start))
if __name__ == '__main__':
main()

View File

@ -324,3 +324,76 @@ def phi_semiprime(n, p, q):
return n - p
else:
return n - (p + q) + 1
def phi(n, primes):
# If n is primes, phi(n)=n-1.
if primes[n] == 1:
return n - 1
# If n is semiprime, use above function.
semi_p, p, q = is_semiprime(n, primes)
if semi_p:
return phi_semiprime(n, p, q)
ph = n
# If 2 is a factor of n, multiply the current ph (which now is n)
# by 1-1/2, then divide all factors 2.
if n % 2 == 0:
ph = ph * (1 - 1 / 2)
while True:
n = n // 2
if n % 2 != 0:
break
# If 3 is a factor of n, multiply the current ph by 1-1/3,
# then divide all factors 3.
if n % 3 == 0:
ph = ph * (1 - 1 / 3)
while True:
n = n // 3
if n % 3 != 0:
break
# Any number can have only one prime factor greater than its
# square root, so we can stop checking at this point and deal
# with the only factor larger than sqrt(n), if present, at the end
limit = floor(sqrt(n)) + 1
# Every prime other than 2 and 3 is in the form 6k+1 or 6k-1.
# If I check all those value no prime factors of the number
# will be missed. For each of these possible primes, check if
# they are prime, then check if the number divides n, in which
# case update the current ph.
for i in range(5, limit, 6):
if primes[i]:
if n % i == 0:
ph = ph * (1 - 1 / i)
while True:
n = n // i
if n % i != 0:
break
if primes[i+2]:
if n % (i + 2) == 0:
ph = ph * (1 - 1 / (i + 2))
while True:
n = n // (i + 2)
if n % (i + 2) != 0:
break
# After dividing all prime factors smaller than sqrt(n), n is either 1
# or is equal to the only prime factor greater than sqrt(n). In this
# second case, we need to update ph with the last prime factor.
if n > 1:
ph = ph * (1 - 1 / n)
return ph