comparison scripts/specfun/perms.m @ 31117:46e15523ca06

perms.m: Small cleanups for Octave coding conventions (bug #60364) * perms.m: Wrap long lines in documentation to < 80 characters. Change output in documentation example to match what Octave actually produces. Use true/false for boolean variable "unique_v" rather than 0/1. Cuddle parentheses when doing indexing and use a space when calling a function. Add FIXME notes requesting an explanation of the apparently complicated algorithm being used for permutations and unque permutations. Remove period at end of error() message text per Octave conventions. Change BIST input validation to more precisely check error() message.
author Rik <rik@octave.org>
date Tue, 05 Jul 2022 08:57:15 -0700
parents 05729b0e39e4
children fd29c7a50a78
comparison
equal deleted inserted replaced
31116:7d3bda173b63 31117:46e15523ca06
30 ## 30 ##
31 ## Results are returned in inverse lexicographic order. The result has size 31 ## Results are returned in inverse lexicographic order. The result has size
32 ## @code{factorial (@var{n}) * @var{n}}, where @var{n} is the length of 32 ## @code{factorial (@var{n}) * @var{n}}, where @var{n} is the length of
33 ## @var{v}. Any repeated elements are included in the output. 33 ## @var{v}. Any repeated elements are included in the output.
34 ## 34 ##
35 ## If the optional argument @qcode{"unique"} is given, then only unique 35 ## If the optional argument @qcode{"unique"} is given then only unique
36 ## permutations are returned, using less memory and generally taking less time 36 ## permutations are returned, using less memory and generally taking less time
37 ## than calling @code{unique (perms (@var{v}), "rows")}. 37 ## than calling @code{unique (perms (@var{v}), "rows")}.
38 ## 38 ##
39 ## Example 1 39 ## Example 1
40 ## 40 ##
55 ## 55 ##
56 ## @example 56 ## @example
57 ## @group 57 ## @group
58 ## perms ([1, 1, 2, 2], "unique") 58 ## perms ([1, 1, 2, 2], "unique")
59 ## @result{} 59 ## @result{}
60 ## 1 1 2 2
61 ## 1 2 1 2
62 ## 1 2 2 1
63 ## 2 1 1 2
64 ## 2 1 2 1
60 ## 2 2 1 1 65 ## 2 2 1 1
61 ## 2 1 2 1
62 ## 2 1 1 2
63 ## 1 2 2 1
64 ## 1 2 1 2
65 ## 1 1 2 2
66 ## @end group 66 ## @end group
67 ## @end example 67 ## @end example
68 ## 68 ##
69 ## Programming Note: If the @qcode{"unique"} option is not used, the length of 69 ## Programming Note: If the @qcode{"unique"} option is not used, the length of
70 ## @var{v} should be no more than 10-12 to limit memory consumption. Even with 70 ## @var{v} should be no more than 10-12 to limit memory consumption. Even with
71 ## @qcode{"unique"}, there should be no more than 10-12 unique elements in @var{v}. 71 ## @qcode{"unique"}, there should be no more than 10-12 unique elements in
72 ## @var{v}.
72 ## @seealso{permute, randperm, nchoosek} 73 ## @seealso{permute, randperm, nchoosek}
73 ## @end deftypefn 74 ## @end deftypefn
74 75
75 ## FIXME: In principle it should be more efficient to do indexing using uint8 76 ## FIXME: In principle it should be more efficient to do indexing using uint8
76 ## type. However, benchmarking shows doubles are faster. If this changes in 77 ## type. However, benchmarking shows doubles are faster. If this changes in
80 81
81 if (nargin < 1) 82 if (nargin < 1)
82 print_usage (); 83 print_usage ();
83 endif 84 endif
84 85
85 unique_v = 0; 86 unique_v = false;
86 if (nargin == 2) 87 if (nargin == 2)
87 if (! strcmpi (opt, "unique")) 88 if (! strcmpi (opt, "unique"))
88 error ('perms: option must be the string "unique".'); 89 error ('perms: option must be the string "unique"');
89 else
90 unique_v = 1;
91 endif 90 endif
91 unique_v = true;
92 endif 92 endif
93 93
94 v = v(:).'; # convert to row vector 94 v = v(:).'; # convert to row vector
95 if (isnumeric (v) || ischar (v)) 95 if (isnumeric (v) || ischar (v))
96 ## Order of output is only dependent on the actual values for 96 ## Order of output is only dependent on the actual values for
111 P = v([3 2 1; 3 1 2; 2 3 1; 2 1 3; 1 3 2; 1 2 3]); 111 P = v([3 2 1; 3 1 2; 2 3 1; 2 1 3; 1 3 2; 1 2 3]);
112 endswitch 112 endswitch
113 if (unique_v) 113 if (unique_v)
114 P = unique (P, "rows"); 114 P = unique (P, "rows");
115 endif 115 endif
116 else 116
117 if (! unique_v) 117 elseif (! unique_v)
118 v = v(end:-1:1); 118 ## FIXME: Brief explanation of the algorithm being used would be useful.
119 n-= 1; 119 v = v(end:-1:1);
120 120 n-= 1;
121 idx = zeros (factorial (n), n); 121
122 idx(1:6, n-2:n) = [1, 2, 3;1, 3, 2;2, 1, 3;2, 3, 1;3, 1, 2;3, 2, 1]+(n-3); 122 idx = zeros (factorial (n), n);
123 f = 2; # jump-start for efficiency with medium n 123 idx(1:6, n-2:n) = [1, 2, 3;1, 3, 2;2, 1, 3;2, 3, 1;3, 1, 2;3, 2, 1]+(n-3);
124 for j = 3:n-1 124 f = 2; # jump-start for efficiency with medium n
125 b = 1:n; 125 for j = 3:n-1
126 f *= j; 126 b = 1:n;
127 perm = idx(1:f, n-(j-1):n); 127 f *= j;
128 idx(1:(j+1)*f, n-j) = (n-j:n)(ones (f, 1),:)(:); 128 perm = idx(1:f, n-(j-1):n);
129 for i=0:j 129 idx(1:(j+1)*f, n-j) = (n-j:n)(ones (f, 1),:)(:);
130 b(i+n-j) -= 1; 130 for i = 0:j
131 idx((1:f)+i*f, n-(j-1):n) = b(perm); 131 b(i+n-j) -= 1;
132 endfor 132 idx((1:f)+i*f, n-(j-1):n) = b(perm);
133 endfor 133 endfor
134 134 endfor
135 n += 1; 135
136 f *= n-1; 136 n += 1;
137 P = v(1)(ones (factorial (n), n)); 137 f *= n-1;
138 P(:,1) = v(ones (f, 1),:)(:); 138 P = v(1)(ones (factorial (n), n));
139 139 P(:,1) = v(ones (f, 1),:)(:);
140 for i = 1:n 140
141 b = v([1:i-1 i+1:n]); 141 for i = 1:n
142 P((1:f)+(i-1)*f, 2:end) = b(idx); 142 b = v([1:i-1 i+1:n]);
143 endfor 143 P((1:f)+(i-1)*f, 2:end) = b(idx);
144 else # user wants unique permutations 144 endfor
145 [v, i, j] = unique(v); 145
146 h = accumarray (j, 1)'; 146 else # unique permutations
147 idx = m_perms (h); 147 [v, ~, j] = unique (v);
148 P = v (sortrows (idx')); 148 h = accumarray (j, 1)';
149 endif 149 idx = m_perms (h);
150 P = v(sortrows (idx'));
150 endif 151 endif
151 152
152 endfunction 153 endfunction
153 154
154 function out_perms = m_perms (multis) 155 function out_perms = m_perms (multis)
156 ## FIXME: Brief explanation of the algorithm being used would be useful.
155 157
156 l = numel (multis); 158 l = numel (multis);
157 if (l == 1) 159 if (l == 1)
158 out_perms = uint8 (ones (multis, 1)); 160 out_perms = uint8 (ones (multis, 1));
159 else 161 else
160 p1 = m_perms (multis (1:floor (l/2))); 162 p1 = m_perms (multis (1:floor (l/2)));
161 p2 = m_perms (multis (floor(l/2)+1:l)) + max (p1 (:, 1)); 163 p2 = m_perms (multis (floor (l/2+1):l)) + max (p1(:, 1));
162 l1 = rows (p1); 164 l1 = rows (p1);
163 l2 = rows (p2); 165 l2 = rows (p2);
164 cp1 = columns (p1); 166 cp1 = columns (p1);
165 cp2 = columns (p2); 167 cp2 = columns (p2);
166 168
167 p = nchoosek (1:l1+l2, l1); 169 p = nchoosek (1:l1+l2, l1);
168 rp = rows(p); 170 rp = rows (p);
169 171
170 ii = false (l1+l2, rp); 172 ii = false (l1+l2, rp);
171 ii(p + (0:rp - 1)' * (l1 + l2)) = true; 173 ii(p + (0:rp - 1)' * (l1 + l2)) = true;
172 out_perms = zeros (l1 + l2, cp1 * cp2 * rp, "uint8"); 174 out_perms = zeros (l1 + l2, cp1 * cp2 * rp, "uint8");
173 out_perms (repmat ( ii, cp1 * cp2, 1)(:)) = repmat (p1, cp2, rp)(:); 175 out_perms(repmat ( ii, cp1 * cp2, 1)(:)) = repmat (p1, cp2, rp)(:);
174 out_perms (repmat (!ii, cp1 * cp2, 1)(:)) = repmat (p2, 1, cp1 * rp)(:); 176 out_perms(repmat (!ii, cp1 * cp2, 1)(:)) = repmat (p2, 1, cp1 * rp)(:);
175 endif 177 endif
176 178
177 endfunction 179 endfunction
180
178 181
179 %!assert (rows (perms (1:6)), factorial (6)) 182 %!assert (rows (perms (1:6)), factorial (6))
180 %!assert (perms (pi), pi) 183 %!assert (perms (pi), pi)
181 %!assert (perms ([pi, e]), [pi, e; e, pi]) 184 %!assert (perms ([pi, e]), [pi, e; e, pi])
182 %!assert (perms ([1,2,3]), [3,2,1;3,1,2;2,3,1;2,1,3;1,3,2;1,2,3]) 185 %!assert (perms ([1,2,3]), [3,2,1;3,1,2;2,3,1;2,1,3;1,3,2;1,2,3])
188 191
189 %!assert (sortrows (perms ("abb", "unique")), ["abb"; "bab"; "bba"]); 192 %!assert (sortrows (perms ("abb", "unique")), ["abb"; "bab"; "bba"]);
190 %!assert (size (perms ([1 1 1 1 2 2 2 3 3], "unique")), [1260 9]) 193 %!assert (size (perms ([1 1 1 1 2 2 2 3 3], "unique")), [1260 9])
191 %!assert (size (perms (int8([1 1 1 1 1 2 2 2 2 3 3 3]), "unique")), [27720 12]) 194 %!assert (size (perms (int8([1 1 1 1 1 2 2 2 2 3 3 3]), "unique")), [27720 12])
192 195
193 ## Should work for any array type, such as cells and structs, and not 196 ## Should work for any array type, such as cells and structs,
194 ## only for numeric data. 197 ## and not only for numeric data.
195 198
196 %!assert <*52431> (perms ({1}), {1}) 199 %!assert <*52431> (perms ({1}), {1})
197 %!assert <*52431> (perms ({0.1, "foo"}), {"foo", 0.1; 0.1, "foo"}) 200 %!assert <*52431> (perms ({0.1, "foo"}), {"foo", 0.1; 0.1, "foo"})
198 %!assert <*52431> (perms ({"foo", 0.1}), {0.1, "foo"; "foo", 0.1}) 201 %!assert <*52431> (perms ({"foo", 0.1}), {0.1, "foo"; "foo", 0.1})
199 %!assert <*52431> (perms ({"foo"; 0.1}), {0.1, "foo"; "foo", 0.1}) 202 %!assert <*52431> (perms ({"foo"; 0.1}), {0.1, "foo"; "foo", 0.1})
228 %! assert (perms (reshape (s, 0, 0)), reshape (s, 1, 0)); 231 %! assert (perms (reshape (s, 0, 0)), reshape (s, 1, 0));
229 %! assert (perms (reshape (s, 0, 1)), reshape (s, 1, 0)); 232 %! assert (perms (reshape (s, 0, 1)), reshape (s, 1, 0));
230 233
231 ## Test input validation 234 ## Test input validation
232 %!error <Invalid call> perms () 235 %!error <Invalid call> perms ()
233 %!error <option must be the string> perms (1:5, "foobar") 236 %!error <option must be the string "unique"> perms (1:5, "foobar")
234 %!error <option must be the string> perms (1:5, {"foo"}) 237 %!error <option must be the string "unique"> perms (1:5, {"foo"})
235 238