Solved problem 12.

This commit is contained in:
2018-02-04 23:11:24 +01:00
parent 230956e89c
commit b3a5bfac52
2 changed files with 392 additions and 1 deletions

View File

@@ -197,6 +197,211 @@ div#notebook {
<p>We can see that 28 is the first triangle number to have over five divisors.</p>
<p>What is the value of the first triangle number to have over five hundred divisors?</p>
</div>
</div>
</div>
<div class="cell border-box-sizing text_cell rendered">
<div class="prompt input_prompt">
</div>
<div class="inner_cell">
<div class="text_cell_render border-box-sizing rendered_html">
<p>Let's try brute forcing.</p>
</div>
</div>
</div>
<div class="cell border-box-sizing code_cell rendered">
<div class="input">
<div class="prompt input_prompt">In&nbsp;[5]:</div>
<div class="inner_cell">
<div class="input_area">
<div class=" highlight hl-ipython3"><pre><span></span><span class="k">def</span> <span class="nf">get_divisors</span><span class="p">(</span><span class="n">n</span><span class="p">):</span>
<span class="k">return</span> <span class="p">[</span><span class="n">i</span> <span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="nb">range</span><span class="p">(</span><span class="mi">1</span><span class="p">,</span> <span class="nb">int</span><span class="p">(</span><span class="n">n</span><span class="o">/</span><span class="mi">2</span><span class="p">)</span> <span class="o">+</span> <span class="mi">1</span><span class="p">)</span> <span class="k">if</span> <span class="n">n</span> <span class="o">%</span> <span class="n">i</span> <span class="o">==</span> <span class="mi">0</span><span class="p">]</span> <span class="o">+</span> <span class="p">[</span><span class="n">n</span><span class="p">]</span>
<span class="k">assert</span><span class="p">(</span><span class="n">get_divisors</span><span class="p">(</span><span class="mi">28</span><span class="p">)</span> <span class="o">==</span> <span class="p">[</span><span class="mi">1</span><span class="p">,</span><span class="mi">2</span><span class="p">,</span><span class="mi">4</span><span class="p">,</span><span class="mi">7</span><span class="p">,</span><span class="mi">14</span><span class="p">,</span><span class="mi">28</span><span class="p">])</span>
<span class="k">def</span> <span class="nf">triangle_number_generator_function</span><span class="p">():</span>
<span class="n">c</span> <span class="o">=</span> <span class="mi">0</span>
<span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="nb">range</span><span class="p">(</span><span class="mi">1</span><span class="p">,</span> <span class="mi">100000000000</span><span class="p">):</span>
<span class="n">c</span> <span class="o">+=</span> <span class="n">i</span>
<span class="k">yield</span> <span class="n">c</span>
<span class="c1">#ts = triangle_number_generator_function()</span>
<span class="c1">#for t in ts:</span>
<span class="c1"># if len(get_divisors(t)) &gt; 500:</span>
<span class="c1"># print(t)</span>
<span class="c1"># break</span>
</pre></div>
</div>
</div>
</div>
</div>
<div class="cell border-box-sizing text_cell rendered">
<div class="prompt input_prompt">
</div>
<div class="inner_cell">
<div class="text_cell_render border-box-sizing rendered_html">
<p>That failed miserably. We have to come up with something more efficient. I looked at my old Haskell solution which looks like this:</p>
<pre><code>divisor_count' x = product . map (succ . length) . group $ prim_factors x</code></pre>
<p>Add first I did not understand what it does at all. But the algorithm is as follows:</p>
<pre><code># get prime factors, for example for 28
2 * 2 * 7
# group the primes
[[2, 2], [7]]
# get the length of each group and increment it by one
[3, 2]
# get the product of all values
6</code></pre>
<p>Honestly, I have no idea why this gives as the number of divisors. I will implement the solution and then try to come up with an explanation.</p>
</div>
</div>
</div>
<div class="cell border-box-sizing code_cell rendered">
<div class="input">
<div class="prompt input_prompt">In&nbsp;[6]:</div>
<div class="inner_cell">
<div class="input_area">
<div class=" highlight hl-ipython3"><pre><span></span><span class="k">def</span> <span class="nf">is_prime</span><span class="p">(</span><span class="n">n</span><span class="p">,</span> <span class="n">smaller_primes</span><span class="p">):</span>
<span class="k">for</span> <span class="n">s</span> <span class="ow">in</span> <span class="n">smaller_primes</span><span class="p">:</span>
<span class="k">if</span> <span class="n">n</span> <span class="o">%</span> <span class="n">s</span> <span class="o">==</span> <span class="mi">0</span><span class="p">:</span>
<span class="k">return</span> <span class="kc">False</span>
<span class="k">if</span> <span class="n">s</span> <span class="o">*</span> <span class="n">s</span> <span class="o">&gt;</span> <span class="n">n</span><span class="p">:</span>
<span class="k">return</span> <span class="kc">True</span>
<span class="k">return</span> <span class="kc">True</span>
<span class="k">def</span> <span class="nf">prime_generator_function</span><span class="p">():</span>
<span class="n">primes</span> <span class="o">=</span> <span class="p">[</span><span class="mi">2</span><span class="p">,</span> <span class="mi">3</span><span class="p">,</span> <span class="mi">5</span><span class="p">,</span> <span class="mi">7</span><span class="p">]</span>
<span class="k">for</span> <span class="n">p</span> <span class="ow">in</span> <span class="n">primes</span><span class="p">:</span>
<span class="k">yield</span> <span class="n">p</span>
<span class="k">while</span> <span class="kc">True</span><span class="p">:</span>
<span class="n">p</span> <span class="o">+=</span> <span class="mi">2</span>
<span class="k">if</span> <span class="n">is_prime</span><span class="p">(</span><span class="n">p</span><span class="p">,</span> <span class="n">primes</span><span class="p">):</span>
<span class="n">primes</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="n">p</span><span class="p">)</span>
<span class="k">yield</span> <span class="n">p</span>
<span class="k">def</span> <span class="nf">get_prime_factors</span><span class="p">(</span><span class="n">number</span><span class="p">):</span>
<span class="n">prime_generator</span> <span class="o">=</span> <span class="n">prime_generator_function</span><span class="p">()</span>
<span class="n">remainder</span> <span class="o">=</span> <span class="n">number</span>
<span class="n">factors</span> <span class="o">=</span> <span class="p">[]</span>
<span class="k">for</span> <span class="n">p</span> <span class="ow">in</span> <span class="n">prime_generator</span><span class="p">:</span>
<span class="k">while</span> <span class="n">remainder</span> <span class="o">%</span> <span class="n">p</span> <span class="o">==</span> <span class="mi">0</span><span class="p">:</span>
<span class="n">remainder</span> <span class="o">/=</span> <span class="n">p</span>
<span class="n">factors</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="n">p</span><span class="p">)</span>
<span class="k">if</span> <span class="n">remainder</span> <span class="o">==</span> <span class="mi">1</span> <span class="ow">or</span> <span class="n">p</span> <span class="o">*</span> <span class="n">p</span> <span class="o">&gt;</span> <span class="n">number</span><span class="p">:</span>
<span class="k">break</span>
<span class="k">if</span> <span class="n">remainder</span> <span class="o">!=</span> <span class="mi">1</span><span class="p">:</span>
<span class="n">factors</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="n">remainder</span><span class="p">)</span>
<span class="k">return</span> <span class="n">factors</span>
</pre></div>
</div>
</div>
</div>
</div>
<div class="cell border-box-sizing text_cell rendered">
<div class="prompt input_prompt">
</div>
<div class="inner_cell">
<div class="text_cell_render border-box-sizing rendered_html">
<p>This are the prime numbers related functions we already now. Now we implement a group function, the product function we already know, and based on that the algorithm to get the number of divisors mentioned above.</p>
</div>
</div>
</div>
<div class="cell border-box-sizing code_cell rendered">
<div class="input">
<div class="prompt input_prompt">In&nbsp;[7]:</div>
<div class="inner_cell">
<div class="input_area">
<div class=" highlight hl-ipython3"><pre><span></span><span class="k">def</span> <span class="nf">group</span><span class="p">(</span><span class="n">xs</span><span class="p">):</span>
<span class="kn">from</span> <span class="nn">functools</span> <span class="k">import</span> <span class="n">reduce</span>
<span class="k">def</span> <span class="nf">f</span><span class="p">(</span><span class="n">xss</span><span class="p">,</span> <span class="n">x</span><span class="p">):</span>
<span class="k">if</span> <span class="n">xss</span> <span class="ow">and</span> <span class="n">x</span> <span class="ow">in</span> <span class="n">xss</span><span class="p">[</span><span class="o">-</span><span class="mi">1</span><span class="p">]:</span>
<span class="n">xss</span><span class="p">[</span><span class="o">-</span><span class="mi">1</span><span class="p">]</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="n">x</span><span class="p">)</span>
<span class="k">else</span><span class="p">:</span>
<span class="n">xss</span><span class="o">.</span><span class="n">append</span><span class="p">([</span><span class="n">x</span><span class="p">])</span>
<span class="k">return</span> <span class="n">xss</span>
<span class="k">return</span> <span class="n">reduce</span><span class="p">(</span><span class="n">f</span><span class="p">,</span> <span class="n">xs</span><span class="p">,</span> <span class="p">[])</span>
<span class="k">def</span> <span class="nf">product</span><span class="p">(</span><span class="n">xs</span><span class="p">):</span>
<span class="kn">from</span> <span class="nn">functools</span> <span class="k">import</span> <span class="n">reduce</span>
<span class="kn">from</span> <span class="nn">operator</span> <span class="k">import</span> <span class="n">mul</span>
<span class="k">return</span> <span class="n">reduce</span><span class="p">(</span><span class="n">mul</span><span class="p">,</span> <span class="n">xs</span><span class="p">,</span> <span class="mi">1</span><span class="p">)</span>
<span class="k">assert</span><span class="p">(</span><span class="n">group</span><span class="p">([</span><span class="mi">1</span><span class="p">,</span> <span class="mi">1</span><span class="p">,</span> <span class="mi">3</span><span class="p">,</span> <span class="mi">2</span><span class="p">,</span> <span class="mi">2</span><span class="p">])</span> <span class="o">==</span> <span class="p">[[</span><span class="mi">1</span><span class="p">,</span> <span class="mi">1</span><span class="p">],</span> <span class="p">[</span><span class="mi">3</span><span class="p">],</span> <span class="p">[</span><span class="mi">2</span><span class="p">,</span> <span class="mi">2</span><span class="p">]])</span>
<span class="k">def</span> <span class="nf">get_number_of_divisors</span><span class="p">(</span><span class="n">n</span><span class="p">):</span>
<span class="c1"># get prime factors, for example for 28</span>
<span class="n">ps</span> <span class="o">=</span> <span class="n">get_prime_factors</span><span class="p">(</span><span class="n">n</span><span class="p">)</span>
<span class="c1"># group the primes</span>
<span class="n">ps</span> <span class="o">=</span> <span class="n">group</span><span class="p">(</span><span class="n">ps</span><span class="p">)</span>
<span class="c1"># get the length of each group and increment it by one</span>
<span class="n">ps</span> <span class="o">=</span> <span class="nb">map</span><span class="p">(</span><span class="k">lambda</span> <span class="n">x</span><span class="p">:</span> <span class="nb">len</span><span class="p">(</span><span class="n">x</span><span class="p">)</span> <span class="o">+</span> <span class="mi">1</span><span class="p">,</span> <span class="n">ps</span><span class="p">)</span>
<span class="c1"># get the product of all values</span>
<span class="k">return</span> <span class="n">product</span><span class="p">(</span><span class="n">ps</span><span class="p">)</span>
<span class="k">assert</span><span class="p">(</span><span class="n">get_number_of_divisors</span><span class="p">(</span><span class="mi">28</span><span class="p">)</span> <span class="o">==</span> <span class="mi">6</span><span class="p">)</span>
</pre></div>
</div>
</div>
</div>
</div>
<div class="cell border-box-sizing text_cell rendered">
<div class="prompt input_prompt">
</div>
<div class="inner_cell">
<div class="text_cell_render border-box-sizing rendered_html">
<p>Now we are ready to do another brute force attempt.</p>
</div>
</div>
</div>
<div class="cell border-box-sizing code_cell rendered">
<div class="input">
<div class="prompt input_prompt">In&nbsp;[8]:</div>
<div class="inner_cell">
<div class="input_area">
<div class=" highlight hl-ipython3"><pre><span></span><span class="n">ts</span> <span class="o">=</span> <span class="n">triangle_number_generator_function</span><span class="p">()</span>
<span class="k">for</span> <span class="n">t</span> <span class="ow">in</span> <span class="n">ts</span><span class="p">:</span>
<span class="k">if</span> <span class="n">get_number_of_divisors</span><span class="p">(</span><span class="n">t</span><span class="p">)</span> <span class="o">&gt;</span> <span class="mi">500</span><span class="p">:</span>
<span class="nb">print</span><span class="p">(</span><span class="n">t</span><span class="p">)</span>
<span class="k">break</span>
</pre></div>
</div>
</div>
</div>
<div class="output_wrapper">
<div class="output">
<div class="output_area"><div class="prompt"></div>
<div class="output_subarea output_stream output_stdout output_text">
<pre>76576500
</pre>
</div>
</div>
</div>
</div>
</div>
<div class="cell border-box-sizing text_cell rendered">
<div class="prompt input_prompt">
</div>
<div class="inner_cell">
<div class="text_cell_render border-box-sizing rendered_html">
<p>Now the only question is why this crazy algorithm works. Okay, I got it with the help of <a href="https://www.math.upenn.edu/~deturck/m170/wk2/numdivisors.html">this</a> page. The problem is actually an instance of the multiplication principle for counting things. Each prime can be used to calculate (n + 1) other divisors. The incrementation by one is required because we can also choose to not use a certain prime. Of course, to get the potential combinations we have to get the product of the potential divisors for all primes. The web page explains it better. The only question is whether I came up with this algorithm myself back then.</p>
</div>
</div>
</div>

View File

@@ -27,6 +27,188 @@
"What is the value of the first triangle number to have over five hundred divisors?"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's try brute forcing."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"def get_divisors(n):\n",
" return [i for i in range(1, int(n/2) + 1) if n % i == 0] + [n]\n",
"\n",
"assert(get_divisors(28) == [1,2,4,7,14,28])\n",
"\n",
"def triangle_number_generator_function():\n",
" c = 0\n",
" for i in range(1, 100000000000):\n",
" c += i\n",
" yield c\n",
"\n",
"#ts = triangle_number_generator_function()\n",
"#for t in ts:\n",
"# if len(get_divisors(t)) > 500:\n",
"# print(t)\n",
"# break\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"That failed miserably. We have to come up with something more efficient. I looked at my old Haskell solution which looks like this:\n",
"\n",
"~~~\n",
"divisor_count' x = product . map (succ . length) . group $ prim_factors x\n",
"~~~\n",
"\n",
"Add first I did not understand what it does at all. But the algorithm is as follows:\n",
"\n",
"~~~\n",
"# get prime factors, for example for 28\n",
"2 * 2 * 7\n",
"# group the primes\n",
"[[2, 2], [7]]\n",
"# get the length of each group and increment it by one\n",
"[3, 2]\n",
"# get the product of all values\n",
"6\n",
"~~~\n",
"\n",
"Honestly, I have no idea why this gives as the number of divisors. I will implement the solution and then try to come up with an explanation."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"def is_prime(n, smaller_primes):\n",
" for s in smaller_primes:\n",
" if n % s == 0:\n",
" return False\n",
" if s * s > n:\n",
" return True\n",
" return True\n",
"\n",
"def prime_generator_function():\n",
" primes = [2, 3, 5, 7]\n",
" for p in primes:\n",
" yield p\n",
" while True:\n",
" p += 2\n",
" if is_prime(p, primes):\n",
" primes.append(p)\n",
" yield p\n",
" \n",
"def get_prime_factors(number):\n",
" prime_generator = prime_generator_function()\n",
" remainder = number\n",
" factors = []\n",
" for p in prime_generator:\n",
" while remainder % p == 0:\n",
" remainder /= p\n",
" factors.append(p)\n",
" if remainder == 1 or p * p > number:\n",
" break\n",
" if remainder != 1:\n",
" factors.append(remainder)\n",
" return factors"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This are the prime numbers related functions we already now. Now we implement a group function, the product function we already know, and based on that the algorithm to get the number of divisors mentioned above."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"def group(xs):\n",
" from functools import reduce\n",
" def f(xss, x):\n",
" if xss and x in xss[-1]:\n",
" xss[-1].append(x)\n",
" else:\n",
" xss.append([x])\n",
" return xss\n",
" return reduce(f, xs, [])\n",
"\n",
"def product(xs):\n",
" from functools import reduce\n",
" from operator import mul\n",
" return reduce(mul, xs, 1)\n",
"\n",
"assert(group([1, 1, 3, 2, 2]) == [[1, 1], [3], [2, 2]])\n",
"\n",
"def get_number_of_divisors(n):\n",
" # get prime factors, for example for 28\n",
" ps = get_prime_factors(n)\n",
" # group the primes\n",
" ps = group(ps)\n",
" # get the length of each group and increment it by one\n",
" ps = map(lambda x: len(x) + 1, ps)\n",
" # get the product of all values\n",
" return product(ps)\n",
" \n",
"assert(get_number_of_divisors(28) == 6)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now we are ready to do another brute force attempt."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"76576500\n"
]
}
],
"source": [
"ts = triangle_number_generator_function()\n",
"for t in ts:\n",
" if get_number_of_divisors(t) > 500:\n",
" print(t)\n",
" break"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now the only question is why this crazy algorithm works. Okay, I got it with the help of [this](https://www.math.upenn.edu/~deturck/m170/wk2/numdivisors.html) page. The problem is actually an instance of the multiplication principle for counting things. Each prime can be used to calculate (n + 1) other divisors. The incrementation by one is required because we can also choose to not use a certain prime. Of course, to get the potential combinations we have to get the product of the potential divisors for all primes. The web page explains it better. The only question is whether I came up with this algorithm myself back then."
]
},
{
"cell_type": "code",
"execution_count": null,
@@ -58,7 +240,11 @@
},
"tags": [
"triangular",
"number"
"number",
"divisors",
"factorization",
"prime",
"group"
]
},
"nbformat": 4,