<< problem 133 - Repunit nonfactors | Same differences - problem 135 >> |
Problem 134: Prime pair connection
(see projecteuler.net/problem=134)
Consider the consecutive primes p_1 = 19 and p_2 = 23.
It can be verified that 1219 is the smallest number such that the last digits are formed by p_1 whilst also being divisible by p_2.
In fact, with the exception of p_1 = 3 and p_2 = 5, for every pair of consecutive primes, p_2 > p_1, there exist values of n for which the last digits
are formed by p_1 and n is divisible by p_2. Let S be the smallest of these values of n.
Find sum{S} for every pair of consecutive primes with 5 <= p_1 <= 1000000.
My Algorithm
I solved this problem twice:
1. a bruteForce
approach (finishes in about 140 seconds)
2. a smarter chineseRemainderTheorem
solution (finishes in about 0.1 second)
Both call tens(x)
which returns the smallest number 10^k that is bigger than x, e.g. tens(456) = 1000
.
My brute-force algorithm consists of a simple loop starting at tens(smallPrime) + smallPrime
(e.g. 100+19=119 when checking the pair 19,23)
and increments until the division by largePrime
has a zero remainder.
Short and simple code but very, very slow ...
I must admit that I heard the name en.wikipedia.org/wiki/Chinese_remainder_theorem years ago but had no idea what it is about.
Wikipedia's explanations weren't exactly clear to me and it took me 3 hours to come up with correct code.
Nevertheless, now I have a small class to compute the extended Euclidean algorithm (see ExtendedGcd
) and learnt
a few things about Monsieur Bezout (en.wikipedia.org/wiki/Étienne_Bézout), too ...
In order to find S I have to solve
n == 0 mod p_2 and
n == p_1 mod t where t_1 = tens(p_1)
The extended Euclidean algorithm of (p_2,t) returns two values x and y.
One solution n is 0 * y * t + p_1 * x * p_2 which can be reduced to p_1 x p_2.
The smallest positive solutions S is n mod (p_2 t).
If this is negative, then I have to add the modulo p_2 t once.
Note
Calling tens()
every time is overkill: when the current prime exceeds the previous value, then just multiply it by 10.
But I doubt that you observe any noticeable performance gains when optimizing that aspect.
Interactive test
This feature is not available for the current problem.
My code
… was written in C++11 and can be compiled with G++, Clang++, Visual C++. You can download it, too.
#include <iostream>
#include <vector>
// return the smallest 10^k bigger than x
// e.g. tens(456) = 1000 => "a 1 followed by as many 0s as x has digits"
unsigned long long tens(unsigned long long x)
{
unsigned long long result = 1;
while (result <= x)
result *= 10;
return result;
}
// when you are too lazy to type and your CPU is too fast ...
// no, seriously: quite useful for verifying small solutions
unsigned long long bruteForce(unsigned long long smallPrime, unsigned long long largePrime)
{
// solve x == 0 mod largePrime
// and x == smallPrime mod tens(smallPrime)
// find 10^k with the minimum amount of zeros
auto shift = tens(smallPrime);
// start with smallest possible value
auto result = shift + smallPrime;
// is it a multiple of b ?
while (result % largePrime != 0)
result += shift; // no, keep going ...
return result;
}
// extended Euclidean algorithm
struct ExtendedGcd
{
// this typedef allows me to switch easily between int and long long
typedef long long Number;
// find solutions x and y (so-called Bezout coefficients)
ExtendedGcd(Number a, Number b)
{
// iterative algorithm from https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm
Number s = 0, lastS = 1;
Number t = 1, lastT = 0;
Number r = b, lastR = a; // remainder
while (r != 0)
{
Number quotient = lastR / r;
Number tmp;
tmp = lastR; lastR = r; r = tmp - quotient * r;
tmp = lastS; lastS = s; s = tmp - quotient * s;
tmp = lastT; lastT = t; t = tmp - quotient * t;
}
// fill members
gcd = lastR;
x = lastS;
y = lastT;
}
// Bezout coefficients
Number x;
Number y;
// just in case we need it, too ...
Number gcd;
};
// use Chinese Remainder Theorem
unsigned long long chineseRemainderTheorem(unsigned int smallPrime, unsigned int largePrime)
{
// solve x == 0 mod largePrime
// and x == smallPrime mod tens(smallPrime)
// e.g. Wolfram Alpha "x = 0 mod 23, x = 19 mod 100" => 1219
auto modulo1 = largePrime;
auto modulo2 = tens(smallPrime);
// apply extended Euclidean algorithm
ExtendedGcd gcd(modulo1, modulo2);
//long long result = 0 * gcd.y * modulo2 + (long long)smallPrime * gcd.x * modulo1;
// => multiplying by zero cancels half of the equation
auto result = smallPrime * gcd.x * modulo1;
// reduce to smallest solution
auto product = modulo1 * modulo2;
result %= (long long)product; // make sure it's a signed modulo
// "too small" ?
if (result < 0)
result += product;
return result;
}
int main()
{
// result
unsigned long long sum = 0;
// sieve based on trial division
unsigned int lastPrime = 2;
std::vector<unsigned int> primes = { lastPrime };
for (unsigned int i = 3; ; i += 2)
{
bool isPrime = true;
// test against all prime numbers we have so far (in ascending order)
for (auto p : primes)
{
// next prime is too large to be a divisor ?
if (p*p > i)
break;
// divisible ? => not prime
if (i % p == 0)
{
isPrime = false;
break;
}
}
// no prime ?
if (!isPrime)
continue;
auto lastPrime = primes.back();
primes.push_back(i);
// find solution
if (lastPrime >= 5)
{
//sum += bruteForce(lastPrime, i);
sum += chineseRemainderTheorem(lastPrime, i);
}
// done ?
if (i > 1000000)
break;
}
std::cout << sum << std::endl;
return 0;
}
This solution contains 23 empty lines, 34 comments and 2 preprocessor commands.
Benchmark
The correct solution to the original Project Euler problem was found in 0.12 seconds on an Intel® Core™ i7-2600K CPU @ 3.40GHz.
(compiled for x86_64 / Linux, GCC flags: -O3 -march=native -fno-exceptions -fno-rtti -std=gnu++11 -DORIGINAL
)
See here for a comparison of all solutions.
Note: interactive tests run on a weaker (=slower) computer. Some interactive tests are compiled without -DORIGINAL
.
Changelog
May 24, 2017 submitted solution
May 24, 2017 added comments
Difficulty
Project Euler ranks this problem at 45% (out of 100%).
Links
projecteuler.net/thread=134 - the best forum on the subject (note: you have to submit the correct solution first)
Code in various languages:
C# www.mathblog.dk/project-euler-134-prime-pair-connection/ (written by Kristian Edlund)
Java github.com/nayuki/Project-Euler-solutions/blob/master/java/p134.java (written by Nayuki)
Those links are just an unordered selection of source code I found with a semi-automatic search script on Google/Bing/GitHub/whatever.
You will probably stumble upon better solutions when searching on your own.
Maybe not all linked resources produce the correct result and/or exceed time/memory limits.
Heatmap
Please click on a problem's number to open my solution to that problem:
green | solutions solve the original Project Euler problem and have a perfect score of 100% at Hackerrank, too | |
yellow | solutions score less than 100% at Hackerrank (but still solve the original problem easily) | |
gray | problems are already solved but I haven't published my solution yet | |
blue | solutions are relevant for Project Euler only: there wasn't a Hackerrank version of it (at the time I solved it) or it differed too much | |
orange | problems are solved but exceed the time limit of one minute or the memory limit of 256 MByte | |
red | problems are not solved yet but I wrote a simulation to approximate the result or verified at least the given example - usually I sketched a few ideas, too | |
black | problems are solved but access to the solution is blocked for a few days until the next problem is published | |
[new] | the flashing problem is the one I solved most recently |
I stopped working on Project Euler problems around the time they released 617.
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 | 165 | 166 | 167 | 168 | 169 | 170 | 171 | 172 | 173 | 174 | 175 |
176 | 177 | 178 | 179 | 180 | 181 | 182 | 183 | 184 | 185 | 186 | 187 | 188 | 189 | 190 | 191 | 192 | 193 | 194 | 195 | 196 | 197 | 198 | 199 | 200 |
201 | 202 | 203 | 204 | 205 | 206 | 207 | 208 | 209 | 210 | 211 | 212 | 213 | 214 | 215 | 216 | 217 | 218 | 219 | 220 | 221 | 222 | 223 | 224 | 225 |
226 | 227 | 228 | 229 | 230 | 231 | 232 | 233 | 234 | 235 | 236 | 237 | 238 | 239 | 240 | 241 | 242 | 243 | 244 | 245 | 246 | 247 | 248 | 249 | 250 |
251 | 252 | 253 | 254 | 255 | 256 | 257 | 258 | 259 | 260 | 261 | 262 | 263 | 264 | 265 | 266 | 267 | 268 | 269 | 270 | 271 | 272 | 273 | 274 | 275 |
276 | 277 | 278 | 279 | 280 | 281 | 282 | 283 | 284 | 285 | 286 | 287 | 288 | 289 | 290 | 291 | 292 | 293 | 294 | 295 | 296 | 297 | 298 | 299 | 300 |
301 | 302 | 303 | 304 | 305 | 306 | 307 | 308 | 309 | 310 | 311 | 312 | 313 | 314 | 315 | 316 | 317 | 318 | 319 | 320 | 321 | 322 | 323 | 324 | 325 |
326 | 327 | 328 | 329 | 330 | 331 | 332 | 333 | 334 | 335 | 336 | 337 | 338 | 339 | 340 | 341 | 342 | 343 | 344 | 345 | 346 | 347 | 348 | 349 | 350 |
351 | 352 | 353 | 354 | 355 | 356 | 357 | 358 | 359 | 360 | 361 | 362 | 363 | 364 | 365 | 366 | 367 | 368 | 369 | 370 | 371 | 372 | 373 | 374 | 375 |
376 | 377 | 378 | 379 | 380 | 381 | 382 | 383 | 384 | 385 | 386 | 387 | 388 | 389 | 390 | 391 | 392 | 393 | 394 | 395 | 396 | 397 | 398 | 399 | 400 |
401 | 402 | 403 | 404 | 405 | 406 | 407 | 408 | 409 | 410 | 411 | 412 | 413 | 414 | 415 | 416 | 417 | 418 | 419 | 420 | 421 | 422 | 423 | 424 | 425 |
426 | 427 | 428 | 429 | 430 | 431 | 432 | 433 | 434 | 435 | 436 | 437 | 438 | 439 | 440 | 441 | 442 | 443 | 444 | 445 | 446 | 447 | 448 | 449 | 450 |
451 | 452 | 453 | 454 | 455 | 456 | 457 | 458 | 459 | 460 | 461 | 462 | 463 | 464 | 465 | 466 | 467 | 468 | 469 | 470 | 471 | 472 | 473 | 474 | 475 |
476 | 477 | 478 | 479 | 480 | 481 | 482 | 483 | 484 | 485 | 486 | 487 | 488 | 489 | 490 | 491 | 492 | 493 | 494 | 495 | 496 | 497 | 498 | 499 | 500 |
501 | 502 | 503 | 504 | 505 | 506 | 507 | 508 | 509 | 510 | 511 | 512 | 513 | 514 | 515 | 516 | 517 | 518 | 519 | 520 | 521 | 522 | 523 | 524 | 525 |
526 | 527 | 528 | 529 | 530 | 531 | 532 | 533 | 534 | 535 | 536 | 537 | 538 | 539 | 540 | 541 | 542 | 543 | 544 | 545 | 546 | 547 | 548 | 549 | 550 |
551 | 552 | 553 | 554 | 555 | 556 | 557 | 558 | 559 | 560 | 561 | 562 | 563 | 564 | 565 | 566 | 567 | 568 | 569 | 570 | 571 | 572 | 573 | 574 | 575 |
576 | 577 | 578 | 579 | 580 | 581 | 582 | 583 | 584 | 585 | 586 | 587 | 588 | 589 | 590 | 591 | 592 | 593 | 594 | 595 | 596 | 597 | 598 | 599 | 600 |
601 | 602 | 603 | 604 | 605 | 606 | 607 | 608 | 609 | 610 | 611 | 612 | 613 | 614 | 615 | 616 | 617 | 618 | 619 | 620 | 621 | 622 | 623 | 624 | 625 |
626 | 627 | 628 | 629 | 630 | 631 | 632 | 633 | 634 | 635 | 636 | 637 | 638 | 639 | 640 | 641 | 642 | 643 | 644 | 645 | 646 | 647 | 648 | 649 | 650 |
651 | 652 | 653 | 654 | 655 | 656 | 657 | 658 | 659 | 660 | 661 | 662 | 663 | 664 | 665 | 666 | 667 | 668 | 669 | 670 | 671 | 672 | 673 | 674 | 675 |
676 | 677 | 678 | 679 | 680 |
I scored 13526 points (out of 15700 possible points, top rank was 17 out of ≈60000 in August 2017) at Hackerrank's Project Euler+.
My username at Project Euler is stephanbrumme while it's stbrumme at Hackerrank.
Look at my progress and performance pages to get more details.
Copyright
I hope you enjoy my code and learn something - or give me feedback how I can improve my solutions.
All of my solutions can be used for any purpose and I am in no way liable for any damages caused.
You can even remove my name and claim it's yours. But then you shall burn in hell.
The problems and most of the problems' images were created by Project Euler.
Thanks for all their endless effort !!!
<< problem 133 - Repunit nonfactors | Same differences - problem 135 >> |