diff -r d9ff240f1339 scripts/set/unique.m
--- a/scripts/set/unique.m	Tue Feb 20 14:02:11 2024 -0800
+++ b/scripts/set/unique.m	Wed Feb 21 11:40:46 2024 +1100
@@ -84,9 +84,6 @@
 ## outputs @var{i}, @var{j} will follow the shape of the input @var{x} rather
 ## than always being column vectors.
 ##
-## The third output, @var{j}, has not been implemented yet when the sort
-## order is @qcode{"stable"}.
-##
 ## @seealso{uniquetol, union, intersect, setdiff, setxor, ismember}
 ## @end deftypefn
 
@@ -232,11 +229,55 @@
 
   ## Calculate j output (3rd output)
   if (nargout > 2)
-    j = i;  # cheap way to copy dimensions
-    j(i) = cumsum ([1; ! match(:)]);
-    if (! optsorted)
-      warning ("unique: third output J is not yet implemented");
-      j = [];
+    if (optsorted)
+      j = i;  # cheap way to copy dimensions
+      j(i) = cumsum ([1; ! match(:)]);
+    else
+      ## Find the stable output order reverse mapping, x=y(j), and avoid
+      ## any element comparisons. The following mappings apply equally for
+      ## a vector, x(i), and for matrix rows, x(i,:).
+      ## An example:
+      ##  a = [3 2 1 2; 2 1 2 1];
+      ##  [o1, o2, o3] = unique (a, "stable");
+      ##  x = vec(a) = [3;2;2;1;1;2;2;1]
+      ##  [y,i] = sort(x)
+      ##  y = [1;1;1;2;2;2;2;3] -> [3;2;1]
+      ##  i = [4;5;8;2;3;6;7;1]
+      ##  match = [1;1;0;1;1;1;0]
+      ##  u =    [1;4;8]
+      ##  i(u) = [4;2;1]
+      ##  s =    [3;2;1]
+      ##  l = [1;1;1;4;4;4;4;8]
+      ## For k=1, i(8)=1 -> l(8)=8 -> y(1) -> j(1)=1
+      ## For k=2, i(4)=2 -> l(4)=4 -> y(2) -> j(2)=2
+      ## For k=3, i(5)=3 -> l(5)=4 -> y(2) -> j(3)=2
+      ## For k=4, i(1)=4 -> l(1)=1 -> y(3) -> j(4)=3
+      ## For k=5, i(2)=5 -> l(2)=1 -> y(3) -> j(5)=3
+      ## For k=6, i(6)=6 -> l(6)=4 -> y(2) -> j(6)=2
+      ## For k=7, i(7)=7 -> l(7)=4 -> y(2) -> j(7)=2
+      ## For k=8, i(3)=8 -> l(3)=1 -> y(3) -> j(8)=3
+
+      ## Find the indexes of the unique elements of sort(x)
+      u = find ([0; match(:)] == 0);
+      nu = numel (u);
+
+      ## The sorted output is y=x(i(u)) and the stable output is y=x(i(u(s)))
+      [~, s] = sort (i(u));
+
+      ## Make a mapping of the unique elements of sort(x)
+      l = u(cumsum ([1; ! match(:)]));
+
+      ## Now find the reverse mapping x=y(j)
+      ## For the above example sort(i) returns:
+      ##  j = [ 8; 4; 5; 1; 2; 6; 7; 3 ]
+      ## and (u(s) == l(j).') returns:
+      ##  [ 1  0  0  0  0  0  0  0
+      ##    0  1  1  0  0  1  1  0
+      ##    0  0  0  1  1  0  0  1 ]
+      [~, j] = sort (i);
+      j = rem (find (u(s) == l(j).'), nu);
+      j(j==0) = nu;
+
     endif
 
     if (optlegacy && isrowvec)
@@ -302,11 +343,10 @@
 %! assert (j, [1;1;2;3;3;3;4]);
 
 %!test
-%! [y,i,~] = unique ([4,4,2,2,2,3,1], "stable");
+%! [y,i,j] = unique ([4,4,2,2,2,3,1], "stable");
 %! assert (y, [4,2,3,1]);
 %! assert (i, [1;3;6;7]);
-%! ## FIXME: 'j' input not calculated with stable
-%! ##assert (j, []);
+%! assert (j, [1;1;2;2;2;3;4]);
 
 %!test
 %! [y,i,j] = unique ([1,1,2,3,3,3,4]', "last");
@@ -335,11 +375,10 @@
 
 %!test
 %! A = [4,5,6; 1,2,3; 4,5,6];
-%! [y,i,~] = unique (A, "rows", "stable");
+%! [y,i,j] = unique (A, "rows", "stable");
 %! assert (y, [4,5,6; 1,2,3]);
 %! assert (A(i,:), y);
-%! ## FIXME: 'j' output not calculated correctly with "stable"
-%! ##assert (y(j,:), A);
+%! assert (y(j,:), A);
 
 ## Test "legacy" option
 %!test
@@ -355,6 +394,52 @@
 %! assert (i, [2; 5; 4; 3]);
 %! assert (j, [4; 1; 4; 3; 2]);
 
+%!test <*65176>
+%! a = [3 2 1 2; 2 1 2 1];
+%! [o1, o2, o3] = unique (a);
+%! assert ({o1, o2, o3}, {[1;2;3], [4;2;1], [3;2;2;1;1;2;2;1]});
+%! [o1, o2, o3] = unique (a, "stable");
+%! assert ({o1, o2, o3}, {[3;2;1], [1;2;4], [1;2;2;3;3;2;2;3]})
+
+%!test <*65176>
+%! a = [3 2 1 2; 2 1 2 1];
+%! [o1, o2, o3] = unique (a(1,:), "rows");
+%! assert ({o1, o2, o3}, {a(1,:), 1, 1});
+%! [o1, o2, o3] = unique (a(1,:), "rows", "stable");
+%! assert ({o1, o2, o3}, {a(1,:), 1, 1});
+%! [o1, o2, o3] = unique (a, "rows");
+%! assert ({o1, o2, o3}, {[a(2,:); a(1,:)], [2;1], [2;1]});
+%! [o1, o2, o3] = unique (a, "rows", "stable");
+%! assert ({o1, o2, o3}, {a, [1;2], [1;2]});
+%! [o1, o2, o3] = unique ([a;a], "rows");
+%! assert ({o1, o2, o3}, {[a(2,:); a(1,:)], [2;1], [2;1;2;1]});
+%! [o1, o2, o3] = unique ([a;a], "rows", "stable");
+%! assert ({o1, o2, o3}, {a, [1;2], [1;2;1;2]});
+
+%!test <*65176>
+%! a = gallery ("integerdata", [-100, 100], 6, 6);
+%! a = [a(2,:); a(1:5,:); a(2:6,:)];
+%! [o1, o2, o3] = unique (a);
+%! assert ({o1, o1(o3), o2, o3}, {a(:)(o2), a(:), ...
+%! [26;22;34;45;57; 6;11;17;33;28;35;15;56; 2;59; 4;66; ...
+%!  16;50;49;27;24;37;44;48;39;38;13;23; 5;12;46;55; 1], ...
+%! [34;14;34;16;30; 6;34;16;30; 6; 7;31;28;31;12;18; 8;31;12;18; 8; 2;29; ...
+%!  22;29; 1;21;10;29; 1;21;10; 9; 3;11; 3;23;27;26; 3;23;27;26;24; 4;32; ...
+%!  4; 25;20;19; 4;25;20;19;33;13; 5;13;15; 2;24;13;15; 2;24;17]});
+%! [o1, o2, o3] = unique (a, "stable");
+%! assert ({o1, o1(o3), o2, o3}, {a(:)(o2), a(:), ...
+%! [ 1; 2; 4; 5; 6;11;12;13;15;16;17;22;23;24;26;27;28; ...
+%!  33;34;35;37;38;39;44;45;46;48;49;50;55;56;57;59;66], ...
+%! [ 1; 2; 1; 3; 4; 5; 1; 3; 4; 5; 6; 7; 8; 7; 9;10;11; 7; 9;10;11;12;13; ...
+%!  14;13;15;16;17;13;15;16;17;18;19;20;19;21;22;23;19;21;22;23;24;25;26;...
+%!  25;27;28;29;25;27;28;29;30;31;32;31;33;12;24;31;33;12;24;34]});
+%! [o1, o2, o3] = unique (a, "rows");
+%! assert ({o1, o1(o3,:), o2, o3}, {a(o2,:), a, ...
+%! [6;11;2;4;5;1], [6;3;6;4;5;1;6;4;5;1;2]});
+%! [o1, o2, o3] = unique (a, "rows", "stable");
+%! assert ({o1, o1(o3,:), o2, o3}, {a(o2,:), a, ...
+%! [1;2;4;5;6;11], [1;2;1;3;4;5;1;3;4;5;6]});
+
 ## Test input validation
 %!error <Invalid call> unique ()
 %!error <X must be an array or cell array of strings> unique ({1})
@@ -376,6 +461,4 @@
 %!error <invalid option> unique ({"a", "b", "c"}, "rows", "UnknownOption2")
 %!error <invalid option> unique ({"a", "b", "c"}, "UnknownOption1", "last")
 %!warning <"rows" is ignored for cell arrays> unique ({"1"}, "rows");
-%!warning <third output J is not yet implemented>
-%! [y,i,j] = unique ([2,1], "stable");
-%! assert (j, []);
+
