RosettaCodeData/Task/Miller-Rabin-primality-test/Phix/miller-rabin-primality-test...

105 lines
26 KiB
Plaintext

(phixonline)-->
<span style="color: #000080;font-style:italic;">-- this is transpiled (then manually copied) to mpz_prime() in mpfr.js:</span>
<span style="color: #004080;">mpz</span> <span style="color: #000000;">modp47</span> <span style="color: #0000FF;">=</span> <span style="color: #004600;">NULL</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">w</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">witness_ranges</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">mpz_prime_mr</span><span style="color: #0000FF;">(</span><span style="color: #004080;">mpz</span> <span style="color: #000000;">p</span><span style="color: #0000FF;">,</span> <span style="color: #004080;">integer</span> <span style="color: #000000;">k</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">10</span><span style="color: #0000FF;">)</span>
<span style="color: #000080;font-style:italic;">-- deterministic to 3,317,044,064,679,887,385,961,981</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">primes</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">5</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">7</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">11</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">13</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">17</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">19</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">23</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">29</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">31</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">37</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">41</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">43</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">47</span><span style="color: #0000FF;">}</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">,</span><span style="color: #000000;">primes</span><span style="color: #0000FF;">[$])<=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">return</span> <span style="color: #7060A8;">find</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">mpz_get_integer</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">),</span><span style="color: #000000;">primes</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">modp47</span><span style="color: #0000FF;">=</span><span style="color: #004600;">NULL</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">modp47</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"614_889_782_588_491_410"</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- === product(primes), largest &lt; 2^64</span>
<span style="color: #000000;">w</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">()</span>
<span style="color: #000080;font-style:italic;">-- Best known deterministic witnesses for given range and set of bases
-- https://miller-rabin.appspot.com/
-- https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test</span>
<span style="color: #000000;">witness_ranges</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{{</span><span style="color: #008000;">"341_531"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"9345883071009581737"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"1_050_535_501"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"336781006125"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"9639812373923155"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"350_269_456_337"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"4230279247111683200"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"14694767155120705706"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"16641139526367750375"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"55_245_642_489_451"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"2"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"141889084524735"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"1199124725622454117"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"11096072698276303650"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"7_999_252_175_582_851"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"2"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"4130806001517"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"149795463772692060"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"186635894390467037"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"3967304179347715805"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"585_226_005_592_931_977"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"2"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"123635709730000"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"9233062284813009"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"43835965440333360"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"761179012939631437"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"1263739024124850375"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"18_446_744_073_709_551_615"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"2"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"325"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"9375"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"28178"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"450775"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"9780504"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"1795265022"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"318_665_857_834_031_151_167_461"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"2"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"3"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"5"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"7"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"11"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"13"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"17"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"19"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"23"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"29"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"31"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"37"</span><span style="color: #0000FF;">}},</span>
<span style="color: #0000FF;">{</span><span style="color: #008000;">"3_317_044_064_679_887_385_961_981"</span><span style="color: #0000FF;">,{</span><span style="color: #008000;">"2"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"3"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"5"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"7"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"11"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"13"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"17"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"19"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"23"</span><span style="color: #0000FF;">,</span>
<span style="color: #008000;">"29"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"31"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"37"</span><span style="color: #0000FF;">,</span> <span style="color: #008000;">"41"</span><span style="color: #0000FF;">}}}</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">(</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">1</span><span style="color: #0000FF;">])</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">2</span><span style="color: #0000FF;">])</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">2</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">(</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">2</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">])</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #7060A8;">mpz_gcd</span><span style="color: #0000FF;">(</span><span style="color: #000000;">w</span><span style="color: #0000FF;">,</span><span style="color: #000000;">p</span><span style="color: #0000FF;">,</span><span style="color: #000000;">modp47</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">w</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)!=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">return</span> <span style="color: #004600;">false</span> <span style="color: #000080;font-style:italic;">-- eliminates 86.2% of all integers</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #000080;font-style:italic;">--
-- Choose input witness bases:
--</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">witnesses</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">,</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[$][</span><span style="color: #000000;">1</span><span style="color: #0000FF;">])>=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">witnesses</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">k</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">k</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">mpz</span> <span style="color: #000000;">a</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">()</span>
<span style="color: #7060A8;">mpz_sub_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">p</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpz_rand</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- a := 0..a-1 (cf rand(n) yields 1..n)</span>
<span style="color: #7060A8;">mpz_add_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">witnesses</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">a</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">else</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">,</span><span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">1</span><span style="color: #0000FF;">])<</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">witnesses</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">witness_ranges</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]</span>
<span style="color: #008080;">exit</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #004080;">mpz</span> <span style="color: #000000;">d</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">()</span>
<span style="color: #7060A8;">mpz_sub_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">p</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">mpz</span> <span style="color: #000000;">nm1</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">)</span>
<span style="color: #000080;font-style:italic;">-- d &gt;&gt;= 4 while (d & 0xf) == 0 # suck out factors of 2
-- (d &gt;&gt;= (d & 3)^2; d &gt;&gt;= (d & 1)^1) if d.even? # 4 bits at a time</span>
<span style="color: #008080;">while</span> <span style="color: #7060A8;">mpz_even</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #7060A8;">mpz_fdiv_q_2exp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">witnesses</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">mpz</span> <span style="color: #000000;">b</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">witnesses</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span>
<span style="color: #008080;">if</span> <span style="color: #008080;">not</span> <span style="color: #7060A8;">mpz_divisible_p</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">,</span><span style="color: #000000;">p</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span> <span style="color: #000080;font-style:italic;">-- skip multiples of input</span>
<span style="color: #004080;">mpz</span> <span style="color: #000000;">s</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">),</span>
<span style="color: #000000;">y</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">()</span>
<span style="color: #7060A8;">mpz_powm</span><span style="color: #0000FF;">(</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">p</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- y := b^d % p</span>
<span style="color: #008080;">while</span> <span style="color: #7060A8;">mpz_cmp_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)!=</span><span style="color: #000000;">0</span>
<span style="color: #008080;">and</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">nm1</span><span style="color: #0000FF;">)!=</span><span style="color: #000000;">0</span>
<span style="color: #008080;">and</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">,</span><span style="color: #000000;">nm1</span><span style="color: #0000FF;">)!=</span><span style="color: #000000;">0</span> <span style="color: #008080;">do</span>
<span style="color: #7060A8;">mpz_powm_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">y</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">p</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- y := y^2 mod p</span>
<span style="color: #7060A8;">mpz_mul_2exp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">s</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- s &lt;&lt; 1</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">nm1</span><span style="color: #0000FF;">)!=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_even</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span> <span style="color: #008080;">return</span> <span style="color: #004600;">false</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">return</span> <span style="color: #004600;">true</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<!--