Mercurial > octave
annotate scripts/specfun/isprime.m @ 31051:5f015f2829b7
maint: Merge stable to default
author | Arun Giridhar <arungiridhar@gmail.com> |
---|---|
date | Wed, 01 Jun 2022 17:20:22 -0400 |
parents | 83f9f8bda883 3a93a77dd4cf |
children | 597f3ee61a48 |
rev | line source |
---|---|
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
1 ######################################################################## |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
2 ## |
30564
796f54d4ddbf
update Octave Project Developers copyright for the new year
John W. Eaton <jwe@octave.org>
parents:
30400
diff
changeset
|
3 ## Copyright (C) 2000-2022 The Octave Project Developers |
27918
b442ec6dda5c
use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents:
27803
diff
changeset
|
4 ## |
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
5 ## See the file COPYRIGHT.md in the top-level directory of this |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
6 ## distribution or <https://octave.org/copyright/>. |
5827 | 7 ## |
8 ## This file is part of Octave. | |
9 ## | |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
10 ## Octave is free software: you can redistribute it and/or modify it |
5827 | 11 ## under the terms of the GNU General Public License as published by |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
12 ## the Free Software Foundation, either version 3 of the License, or |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
13 ## (at your option) any later version. |
5827 | 14 ## |
15 ## Octave is distributed in the hope that it will be useful, but | |
16 ## WITHOUT ANY WARRANTY; without even the implied warranty of | |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
17 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22323
diff
changeset
|
18 ## GNU General Public License for more details. |
5827 | 19 ## |
20 ## You should have received a copy of the GNU General Public License | |
7016 | 21 ## along with Octave; see the file COPYING. If not, see |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23220
diff
changeset
|
22 ## <https://www.gnu.org/licenses/>. |
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
23 ## |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
24 ######################################################################## |
5827 | 25 |
26 ## -*- texinfo -*- | |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
27 ## @deftypefn {} {@var{tf} =} isprime (@var{x}) |
20162
2645f9ef8c88
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
28 ## Return a logical array which is true where the elements of @var{x} are prime |
2645f9ef8c88
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19833
diff
changeset
|
29 ## numbers and false where they are not. |
11431
0d9640d755b1
Improve docstrings for all isXXX functions.
Rik <octave@nomad.inbox5.com>
parents:
10664
diff
changeset
|
30 ## |
19052
bb0c5e182c12
isprime.m: Update docstring to note that 1 is not prime.
Rik <rik@octave.org>
parents:
19050
diff
changeset
|
31 ## A prime number is conventionally defined as a positive integer greater than |
bb0c5e182c12
isprime.m: Update docstring to note that 1 is not prime.
Rik <rik@octave.org>
parents:
19050
diff
changeset
|
32 ## 1 (e.g., 2, 3, @dots{}) which is divisible only by itself and 1. Octave |
bb0c5e182c12
isprime.m: Update docstring to note that 1 is not prime.
Rik <rik@octave.org>
parents:
19050
diff
changeset
|
33 ## extends this definition to include both negative integers and complex |
bb0c5e182c12
isprime.m: Update docstring to note that 1 is not prime.
Rik <rik@octave.org>
parents:
19050
diff
changeset
|
34 ## values. A negative integer is prime if its positive counterpart is prime. |
bb0c5e182c12
isprime.m: Update docstring to note that 1 is not prime.
Rik <rik@octave.org>
parents:
19050
diff
changeset
|
35 ## This is equivalent to @code{isprime (abs (x))}. |
19597
db92e7e28e1f
strip trailing whitespace from most source files
John W. Eaton <jwe@octave.org>
parents:
19052
diff
changeset
|
36 ## |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
37 ## If @code{class (@var{x})} is complex, then primality is tested in the domain |
25143
13fd0610480f
doc: Use https whenever possible in @url entries.
Rik <rik@octave.org>
parents:
25054
diff
changeset
|
38 ## of Gaussian integers (@url{https://en.wikipedia.org/wiki/Gaussian_integer}). |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
39 ## Some non-complex integers are prime in the ordinary sense, but not in the |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
40 ## domain of Gaussian integers. For example, @math{5 = (1+2i)*(1-2i)} shows |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
41 ## that 5 is not prime because it has a factor other than itself and 1. |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
42 ## Exercise caution when testing complex and real values together in the same |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
43 ## matrix. |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
44 ## |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
45 ## Examples: |
5827 | 46 ## |
11431
0d9640d755b1
Improve docstrings for all isXXX functions.
Rik <octave@nomad.inbox5.com>
parents:
10664
diff
changeset
|
47 ## @example |
0d9640d755b1
Improve docstrings for all isXXX functions.
Rik <octave@nomad.inbox5.com>
parents:
10664
diff
changeset
|
48 ## @group |
0d9640d755b1
Improve docstrings for all isXXX functions.
Rik <octave@nomad.inbox5.com>
parents:
10664
diff
changeset
|
49 ## isprime (1:6) |
26600
f6730533820e
doc: clean up doc example blocks for accuracy and consistent formatting
Mike Miller <mtmiller@octave.org>
parents:
26376
diff
changeset
|
50 ## @result{} 0 1 1 0 1 0 |
11431
0d9640d755b1
Improve docstrings for all isXXX functions.
Rik <octave@nomad.inbox5.com>
parents:
10664
diff
changeset
|
51 ## @end group |
0d9640d755b1
Improve docstrings for all isXXX functions.
Rik <octave@nomad.inbox5.com>
parents:
10664
diff
changeset
|
52 ## @end example |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
53 ## |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
54 ## @example |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
55 ## @group |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
56 ## isprime ([i, 2, 3, 5]) |
26600
f6730533820e
doc: clean up doc example blocks for accuracy and consistent formatting
Mike Miller <mtmiller@octave.org>
parents:
26376
diff
changeset
|
57 ## @result{} 0 0 1 0 |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
58 ## @end group |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
59 ## @end example |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
60 ## |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
61 ## Programming Note: @code{isprime} is suitable for all @var{x} |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
62 ## in the range abs(@var{x}) |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
63 ## @tex |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
64 ## $ < 2^{64}$. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
65 ## @end tex |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
66 ## @ifnottex |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
67 ## < 2^64. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
68 ## @end ifnottex |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
69 ## |
27803
e8823425a60b
doc: Change markup for "Matlab".
Markus Mützel <markus.muetzel@gmx.de>
parents:
26600
diff
changeset
|
70 ## Compatibility Note: @sc{matlab} does not extend the definition of prime |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
71 ## numbers and will produce an error if given negative or complex inputs. |
5827 | 72 ## @seealso{primes, factor, gcd, lcm} |
73 ## @end deftypefn | |
74 | |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
75 function tf = isprime (x) |
7125 | 76 |
28891
de5f2f9a64ff
maint: Use same coding style when checking for a minimum of 1 input.
Rik <rik@octave.org>
parents:
27923
diff
changeset
|
77 if (nargin < 1) |
19041
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
78 print_usage (); |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
79 elseif (any (fix (x) != x)) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
80 error ("isprime: X contains non-integer entries"); |
19041
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
81 endif |
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
82 |
19045
920a400929ca
* isprime.m: Return an empty array for empty input
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19041
diff
changeset
|
83 if (isempty (x)) |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
84 tf = x; |
19045
920a400929ca
* isprime.m: Return an empty array for empty input
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19041
diff
changeset
|
85 return; |
920a400929ca
* isprime.m: Return an empty array for empty input
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19041
diff
changeset
|
86 endif |
920a400929ca
* isprime.m: Return an empty array for empty input
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19041
diff
changeset
|
87 |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
88 if (iscomplex (x)) |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
89 tf = isgaussianprime (x); |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
90 return; |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
91 endif |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
92 |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
93 ## Code strategy is to quickly compare entries in x with small primes |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
94 ## using lookup(), then do direct division on larger numbers up to |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
95 ## a threshold, then call Miller-Rabin for numbers over the threshold. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
96 |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
97 x = abs (x); # handle negative entries |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
98 |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
99 ## Generate prime table of suitable length up to maxp. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
100 ## The value of maxp needs to be at least 37, |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
101 ## because of the method used by __isprimelarge__ below. |
30400
168da23530b4
maint: strip trailing spaces from code base.
Rik <rik@octave.org>
parents:
30261
diff
changeset
|
102 maxp = 37; |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
103 pr = [2 3 5 7 11 13 17 19 23 29 31 37]; |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
104 tf = lookup (pr, x, "b"); # quick search for table matches. |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
105 |
31050
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
106 THRESHOLD = 195e8; |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
107 ## FIXME: THRESHOLD is the input value at which Miller-Rabin |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
108 ## becomes more efficient than direct division. For smaller numbers, |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
109 ## use direct division. For larger numbers, use Miller-Rabin. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
110 ## |
31050
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
111 ## From numerical experiments in Jun 2022, this was observed: |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
112 ## THRESHOLD Division Miller-Rabin |
31050
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
113 ## 29e9 29.8196s 26.2484s (previous value) |
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
114 ## 20e9 26.7445s 26.0161s |
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
115 ## 10e9 20.9330s 25.3247s |
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
116 ## 19e9 26.5397s 26.8987s |
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
117 ## 195e8 26.5735s 26.4749s |
3a93a77dd4cf
Minor performance tweaks to isprime.m and __isprimelarge__.cc
Arun Giridhar <arungiridhar@gmail.com>
parents:
30564
diff
changeset
|
118 ## which is close enough, so new threshold = 195e8. |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
119 ## |
30400
168da23530b4
maint: strip trailing spaces from code base.
Rik <rik@octave.org>
parents:
30261
diff
changeset
|
120 ## The test code was this: |
168da23530b4
maint: strip trailing spaces from code base.
Rik <rik@octave.org>
parents:
30261
diff
changeset
|
121 ## n = THRESHOLD - (1:1e7); tic; isprime(n); toc |
168da23530b4
maint: strip trailing spaces from code base.
Rik <rik@octave.org>
parents:
30261
diff
changeset
|
122 ## n = THRESHOLD + (1:1e7); tic; isprime(n); toc |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
123 ## |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
124 ## Two notes for future programmers: |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
125 ## |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
126 ## 1. Test and tune THRESHOLD periodically. Miller-Rabin is only CPU-limited, |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
127 ## while factorization by division is very memory-intensive. This is |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
128 ## plainly noticeable from the loudness of the computer fans when running |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
129 ## the division technique! CPU speed and RAM speed scale differently over |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
130 ## time, so test and tune THRESHOLD periodically. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
131 ## |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
132 ## 2. If you make improvements elsewhere in the code that favor one over |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
133 ## the other (not symmetric), you should also retune THRESHOLD afterwards. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
134 ## If the Miller-Rabin part is sped up, the optimum THRESHOLD will |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
135 ## decrease, and if factorization is sped up, it will increase. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
136 ## |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
137 |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
138 ## Process large entries that are still suitable for direct division |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
139 m = x (x > maxp & x <= THRESHOLD); |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
140 if ( ! isempty (m)) |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
141 ## Start by dividing through by the small primes until the remaining list |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
142 ## of entries is small (and most likely prime themselves). |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
143 pr2 = primes (sqrt (max (m))); |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
144 tf |= lookup (pr2, x, "b"); |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
145 for p = pr2 |
19041
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
146 m = m(rem (m, p) != 0); |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
147 if (numel (m) < numel (pr) / 10) |
19041
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
148 break; |
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
149 endif |
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
150 endfor |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
151 |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
152 ## Check the remaining list of possible primes against the remaining |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
153 ## prime factors which were not tested in the for loop. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
154 ## This is just an optimization to use arrayfun over for loop. |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
155 pr2 = pr2 (pr2 > p); |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
156 mm = arrayfun (@(x) all (rem (x, pr2)), m); |
19041
0976f9fccbbd
isprime.m: Verify that input is a real number (bug #43041)
Rik <rik@octave.org>
parents:
17744
diff
changeset
|
157 m = m(mm); |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
158 |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
159 ## Add any remaining entries, which are truly prime, to the results. |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
160 if ( ! isempty (m)) |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
161 tf |= lookup (sort (m), x, "b"); |
10664
faff5367cc05
second isprime rewrite
Jaroslav Hajek <highegg@gmail.com>
parents:
10657
diff
changeset
|
162 endif |
7125 | 163 endif |
164 | |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
165 ## Process remaining entries (everything above THRESHOLD) with Miller-Rabin |
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
166 ii = (x(:)' > THRESHOLD); |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
167 tf(ii) = __isprimelarge__ (x(ii)); |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
168 |
5827 | 169 endfunction |
9984
d1cc2e0ddf55
isprime: produce logical result
John W. Eaton <jwe@octave.org>
parents:
7125
diff
changeset
|
170 |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
171 function tf = isgaussianprime (z) |
21758
ffad2baa90f7
maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
172 |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
173 ## Assume prime unless proven otherwise |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
174 tf = true (size (z)); |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
175 |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
176 x = real (z); |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
177 y = imag (z); |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
178 |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
179 ## If purely real or purely imaginary, ordinary prime test for |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
180 ## that complex part if that part is 3 mod 4. |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
181 xidx = y==0 & mod (x, 4) == 3; |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
182 yidx = x==0 & mod (y, 4) == 3; |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
183 |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
184 tf(xidx) &= isprime (x(xidx)); |
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
185 tf(yidx) &= isprime (y(yidx)); |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
186 |
19050
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
187 ## Otherwise, prime if x^2 + y^2 is prime |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
188 zidx = ! (xidx | yidx); # Skip entries that were already evaluated |
2f117c4b5cb0
isprime.m: Document negative and complex inputs.
Rik <rik@octave.org>
parents:
19046
diff
changeset
|
189 zabs = x(zidx).^2 + y(zidx).^2; |
30558
83aeaba707d8
doc: Use TF for output variable in documentation for isXXX functions in scripts/ directory.
Rik <rik@octave.org>
parents:
30400
diff
changeset
|
190 tf(zidx) &= isprime (zabs); |
21758
ffad2baa90f7
maint: Use newlines to make code more readable.
Rik <rik@octave.org>
parents:
20852
diff
changeset
|
191 |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
192 endfunction |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
193 |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
194 |
14363
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
195 %!assert (isprime (3), true) |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
196 %!assert (isprime (4), false) |
30261
a49c635b179d
isprime.m: Speed up for large integers using Miller-Rabin test (bug #61312).
Arun Giridhar <arungiridhar@gmail.com>
parents:
29359
diff
changeset
|
197 %!assert (isprime (uint64 (18446744073709551557)), true) |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
198 %!assert (isprime (5i), false) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
199 %!assert (isprime (7i), true) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
200 %!assert (isprime ([1+2i, (2+3i)*(-1+2i)]), [true, false]) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
201 %!assert (isprime (-2), true) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
202 %!assert (isprime (complex (-2)), false) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
203 %!assert (isprime (2i), false) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
204 %!assert (isprime ([i, 2, 3, 5]), [false, false, true, false]) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
205 %!assert (isprime (0), false) |
14363
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
206 %!assert (isprime (magic (3)), logical ([0, 0, 0; 1, 1, 1; 0, 0, 1])) |
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
207 |
19833
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19697
diff
changeset
|
208 ## Test input validation |
28896
90fea9cc9caa
test: Add expected error message <Invalid call> to BIST tests for nargin.
Rik <rik@octave.org>
parents:
28891
diff
changeset
|
209 %!error <Invalid call> isprime () |
19046
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
210 %!error <X contains non-integer entries> isprime (0.5i) |
89e275a4f6f6
Allow isprime to handle complex numbers (bug #43041)
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
19045
diff
changeset
|
211 %!error <X contains non-integer entries> isprime (0.5) |