Време е за НАУКА! v1.1

Три огледала са разположени по стените на равностранен триъгълник, с отразяващите страни навътре. Във всеки връх на триъгълника има безкрайно малък отвор, през който може да мине лазерен лъч.

Да кръстим върховете A, B и C. Има два начина, по които лъч може да влезе в триъгълника през връх C, да се отрази от точно 11 повърхности, и да излезе през същия връх. Единият начин е показан долу, другият е по точно обратния път.

p_201_laserbeam.gif

Има 80840 начина по които лазерен лъч може да влезе през C, да отскочи точно 1000001 пъти, и да излезе пак през C.

По колко начина това може да се получи, ако ограничението е да отскочи точно 12017639147 пъти?

Известно още като Проект Ойлер 202

Математика :)

Математиката е оферта. Освен очевидното "преставаш да си тъп" предимство, го има и моментът, в който решиш нещо.
Да, знам, не съм само аз дето съм го направил това.
Пък и дори само аз да го бях решил, по-голямо усилие е положил този, който е създал пъзела. Той пак е по-умен.
И все пак.

Първи опит. Тъпото решение.

Малко факти - траекторията на един лъч е еднозначно определена от ъгъла, под който влиза в триъгълника. Така де, лъчът си се движи по права линия(тук не се занимаваме с някакви там пост-модернистични квантови физики), а след всяко отражение е точно определено накъде ще отиде. Ъгълът на падане е равен на ъгъла на отражение, физика, девети клас.

Така.
Пускам един лъч. Ако се отрази, му променям посоката. Ако се отрази пак - пак му я променям. И така 12017639147 пъти. И ако излезе пак през върха, в който е влязъл - йей, успели сме! Повтори още безкрайност на брой пъти за всеки един възможен ъгъл под който може да влезе в триъгълника един лъч.

Ммне…

Освен това, с днешната технология това не е особено реално да се направи. С компютър мога да умножа две цели числа и да получа верен резултат.

3 * 3 = 9

Обаче, не по съвсем същия начин стоят нещата с нецелите числа.

0.6 * 0.6 = 0.35999999999999999

Малко е трудно за обяснение, но в общи линии, компютрите не смятат съвсем точно. Така де, грешката е мноого малка - само 0.00000000000000001 нали? Да, ама като умножа три числа, ще стане 0.00000000000000002. Същото е и при събирането, и при другите действия. А за да се изчисли ъгъл на отражение се използват по няколко операции наведнъж. За дванадесетте милиарда отражения, които се искат в задачата, грешката ще се натрупа много сериозно.

Втори дубъл, вече сериозно

С пускане на лъчи тая работа няма да стане.
Тук ще трябва да се помисли.

Всъщност, има само едно нещо, за което човек трябва да се усети. Ама е малко хитро. Цялата работа е да го изкараме така, че въобще да няма отражения. Така де, няма човек на който това да не са му го обяснявали в гимназията. Ама кой ли ти слуша… Разбирам те, аз също не слушах, научих го по-късно.

Вместо да смяташ отражения, просто слагаш един въображаем лъч зад огледалото. И още няколко въображаеми огледала. Ето така:

equivalence.PNG

След това просто слагаме дооооста въображаеми неща.

pyramid.png

Всеки ред е отбелязан с поредния си номер, а в скоби е написано колко пъти трябва да се отрази лъч, за да стигне до връх на този ред. Важно е да се отбележи, че до всяка друга точка трябва да се отрази един път повече. При нулевия и последния връх бройката също е малко по-различна, но те са специални - с тях ще се занимаваме после. Не е много трудно да се докаже, че за да се стигне до връх от ред N, трябва да се мине през 2N - 3 отражения. Използва се математическа индукция. Съответно, за K отражения, решенията ще се намират на ред $\cfrac{K + 3}{2}$.
Това означава, че ако искаме да намерим колко решения има за, например 100001 отражения, трябва да търсим на ред $\cfrac{1000001 + 3}{2} = 500002$

Червената и синята линия изобразяват два пътя на лъч, стигащи до ниво 4. Според условието на задачата, ако стигне до връх на триъгълника, лъчът излиза от него, затова синият лъч не е валиден път. Той ще излезе още на ниво 2, където пресича връх. Тоест, до връх с номер 2 на ниво 4 въобще не може да се стигне. Някои върхове отпадат по този начин - тъжно, но факт.

Червеният път е валиден по тази точка, но има нещо друго. Трябва да се види през кой връх излиза.

В началото имаме три надписани върха - A, B, C. При всяко разгръщане на мрежата, тези върхове отиват някъде в пространството. Трябва да се види кой връх къде отива. Оцветяването на върховете е нещо като решаване на Судоку. Ето как изглежда за първите 4 нива:

nodes.PNG

Тоест, това което трябва да се направи всъщност е:

  1. Дадено е N.
  2. Да се намерят всички лъчи, които започват в най-долния ъгъл на пирамидата и завършват в ред $\cfrac{N + 3}{2}$, като при това не пресичат никой връх(т.е. възел) във вътрешността на пирамидата.
  3. От тези лъчи да се отделят тези, които завършват в черен връх.
  4. Каквото е останало отговаря на решението, освен това няма други лъчи, които да са решения. Йей!

Алгоритъм

Неслучайно избрах редовете да са с номера $1,2 \cdots n$. Ще го използваме това.
По какъв начин можем да разберем дали лъч от началната точка до точка K на ред N минава през някой възел вътре в пирамидата? Малко е дълго за обясняване, но ще използваме най-голям общ делител. Това е най-голямото цяло число, което дели едновременно както K, така и N. Ако това число е 1, казваме че K и N са взаимно прости. На английски функцията се казва greatest common divisor, или още gcd. Няма да се спирам защо точно това е така, но със знания от начален курс по дискретна математика, човек може да се досети, че:

Лъчът до K не минава през вътрешни върхове тогава и само тогава, когато gcd(K,N) = 1

Точно по този критерий отпадат върховете 0 и N за всеки ред N.

Аз си имам един файл utils.py. В него си пиша страничните функцийки. Такива, които може да ползвам за повече от една задача. gcd е такава функция(Виж в края на статията за самия файл).

А как се проверява точка 3) - дали сме в правилния връх? Още едно дребно трикче, нарича се деление по модул. Лесно се доказва - ако $(N + K) \equiv 0 (mod \enspace 3)$, то сме във връх C. Горното означава "N + K е кратно на 3".

При което е готова основната идея на задачата. Ето малко код, който илюстрира точно това, което досега обяснявах:

Това е хубаво. Само че има един недостатък. Извършва N пъти операцията gcd. Която пък всеки път извършва приблизително log(N) (това е логаритъм) операции. Тоест, общо NlogN операции.
За непрограмистите, това означава че програмата ще е бавна. При ограничение N = 12 милиарда, ще трябва да направи към 300 милиарда операции. Можем да разчитаме, че компютърът върши по един милион на секунда (в най-добрия случай и при много оптимизации - 10 милиона/сек). За няколко часа ще изкара верния отговор.

Не. Това не е оптимално, не е красиво. На Дийкстра това не би му харесало.
За какво съм прекарал десетилетие в учене, ако това е най-доброто, на което съм способен?
Не ми се чака толкова. Искам всичко да става за една секунда.

Трети опит. Още веднъж, този път с чувство.

Внимание!
Дотук всичко беше сравнително лесно. Такава програма за известен брой минути ще свърши работата.
Отсега нататък ще ми отнема твърде мого време да обяснявам нещата съвсем подробно.
Затова няма да го правя.

Време е да вкараме тежката артилерия.

Идва време да си зададем въпроса - кога gcd(K,N) е по-голямо от 1?
Отговор: Когато някой от простите делители на N дели K.
M е делител на N когато е изпълнено N = M * L, а (M,N,L) всичките са цели.
Просто число е такова число, което се дели само на 1 и на себе си. 7 е просто. Прост делител на N е делител на N, който е прост :)

Метод на включването и изключването

Та, трябва просто да намерим броя на всички числа, които се делят на някой от простите делители на N.
Как става тази работа?

  1. Първо намираме всички прости делители на N(това става бързо и лесно) :).
  2. След това, по метода на включването и изключването лесно се намира самият брой.

Например N е 30. Делителите му са 2, 3 и 5.
Ясно е, че едно от всеки 2 числа се дели на 2. Едно от всеки 3 се дели на 3 и т.н.
Така че, дотук сумата е $S = [\cfrac{N}{2}] + [\cfrac{N}{3}] + [\cfrac{N}{5}]$.
С квадратни скоби се означава цяла част на число. Тоест, [5.25] = 5, [6.75] = 6.
Само че, по този начин числата, които се делят както на 2, така и на 3 са броени 2 пъти. Затова трябва да ги извадим веднъж.
$S -= [\cfrac{N}{2*3}] + [\cfrac{N}{3*5}] + [\cfrac{N}{5*2}]$
Но, по този начин числата, които се делят и на трите делителя, са извадени 2 пъти. Трябва да ги добавим веднъж…
$S += [\cfrac{N}{2*3*5}]$
И така до края. Докато свършат делителите. Това се нарича метод на включването и изключването. Много е полезен, използва се за броене на разни работи :)

Остава само един проблем - трябва да се гледа остатъка по модул 3 (т.е. "В кой връх съм")

Китайска теорема на остатъците

Да приемем, че делителите на N са a1 a2 … am.
В такъв случай, за всяка една комбинация на тези делители(виж редишната точка), трябва да намерим всички числа X, за които:

  • $X \equiv 0 (mod \enspace a_{c_1})$
  • $X \equiv 0 (mod \enspace a_{c_2})$
  • $X \equiv 0 (mod \enspace a_{c_n})$
  • $X \equiv 3 - N (mod \enspace 3)$

Където с $c_1, c_2, \cdots, c_n$ означаваме кои делители участват в съответната комбинация. Последният ред е просто за пореден път условието за върховете.

Съществува една теорема - Китайската Теорема за Остатъците. Тя гласи, че когато имаме система от сравнения по модул (като тази горе), то тя има решение. При това има единствено решение в интервала $[0, a_{c_1} * a_{c_2} * ... * a_{c_n} * 3)$.
Да кръстим горното произведение P. Та, има единствено решение в $[0,P)$. Няма нужда да го търсим точно кое е, просто знаем, че е там. Нещо повече, КТО ни казва, че във всеки интервал $[k*P, (k+1)*P), k >= 0$ има точно по едно решение на тези сравнения. Тоест, ще разделим N на интервали:

[-----------------------N------------------------]
[---P---][---P---][---P---][---P---][---P---][-Q-]

Като последното Q e по-малко от P. В него ще трябва просто да проверим всички числа и да видим дали някое от тях отговаря на условията. Това в най-лошия случай отнема няколко милиона проверки.
Само дето има още един трик :P
Цялата работа може да стане само с 3 проверки :)
За да разбереш как става, прочети кода на програмата.
За да разбереш наистина как става, нарисувай си пирамидата за 11 отражения и виж, че наистина има точно 2 решения :)

Бързодействие

Като цяло, сложността(бързодействието) на алгоритъма е сума от следните:

  • Намиране на всички прости делители на число: $O(\sqrt{N})$
  • Генериране на всички комбинации на делителите $O(2^{L})$, където L е броят на делителите. За 12 милиарда няма как L да е по-голямо от 20, така че това е по-малко от милион (т.е. по-малко от 1 секунда). А точно за примера от задачата е под 1000.
  • За всяка комбинация се гледа колко пъти участва в крайния резултат и се добавя/вади според гореописания алгоритъм. Това отнема константен брой операции (10-тина най-много), почти за нулево време се изпълнява.

Ето и самият код на програмата.


Работи за 0.6 секунди. От които около половината време е прекарано просто в стартиране на програмата и печатане на резултат. Йей, наука!

q.e.d.

Ето частта от utils.py, която използвах в тази задача:

Unless otherwise stated, the content of this page is licensed under Creative Commons Attribution-ShareAlike 3.0 License