@Markus
Done some tests, and based on these observation, I tend to suggest include the initial coplanar test. See attached file bug55751_inc_coplanar_v2.patch for corresponding patch.
BTW, I forgot to turn off the check: facecolor_is ("none") && edgecolor_is ("none") in previous tests.
The (correct) time cost is:
case 1, case 2, case 3,
scatter3 (no coplanar test) : 0.095 s 0.095 s 0.097 s
__go_patch__ (+coplanar test) : 0.094 s 0.227 s 0.225 s
coplanar split : 0.003 s 0.142 s 0.132 s
coplanar split + init full test : 0.003 s 0.145 s 0.003 s
Case 1: peaks (32).
Case 2: "2D Pseudo Quaterion Rotation Program", the everrunning one.
Case 3: coplanar point set. N = 45001; vecs=rand(N,2)*orth(rand(3,2))';
Test code of _go_patch_ (mimic implementation of scatter3):
tic
hax = newplot ();
edgecolor = "r"; % trigger coplanar test, "none" to suppress.
feval ('__go_patch__', hggroup (),
"xdata", x, "ydata", y, "zdata", z,
"facecolor", "none", "edgecolor", edgecolor,
"faces", 1:numel (x), "vertices", [x, y, z],
"marker", "o",
"markeredgecolor", "none",
"markerfacecolor", "b",
"markersize", 3, "linestyle", "none");
toc
So the impact of initial overall coplanar test is very minor, and the improvement to coplanar input case is huge.
A more detailed time cost of different cases (test code below, time in table below is only for the coplanar test C++ code):
with initial coplanar test
nc 4 5 10 20 50 1e6
random : 1.974 2.297 2.818 3.107 3.307 3.391
coplanar: 1.157 0.919 0.492 0.263 0.139 0.067
no initial coplanar test
nc 4 5 10 20 50 1e6
random : 1.193 1.585 2.396 2.806 3.118 3.227
coplanar: 1.115 1.478 2.185 2.605 2.779 2.830
nc: number of corners in a face, total number of vertex is 1e6.
random: random noncoplanar data.
coplanar: coplanar data.
# random
N = 1e6;
nc = 5;
vertices = rand(N, 3);
faces = 1:N;
faces = reshape (faces, nc, [])';
tic
hp = patch ('Vertices', vertices, 'Faces', faces);
toc
# coplanar
N = 1e6;
nc = 5;
vertices = rand(nc, 2) * orth(rand(3,2))';
faces = 1:nc;
faces = faces .* ones(N / nc, 1);
tic
hp = patch ('Vertices', vertices, 'Faces', faces);
toc
This time table provide information that when should we disable the overall coplanar test. i.e. it should be disabled if nc == 4, otherwise the test should be beneficial on average, assuming more than half cases are coplanar.
The attached file bug55751_inc_coplanar_v2_timing.diff is what I used for the timing measure of coplanar test/face splitting. Apply after the patch.
Further speed up might be possilbe, such as use "rcond" for rank test, instead of "eig". I didn't try it.
Appendix:
Here is some addtional background math in case someone (or me, in the future) want to fully understand the algorithm and fine tune it.
Suppose [x y z] are always in a plane:
a x + b y + c z + d = 0
To find {a,b,c,d}, we can minimize the variance of "epsilon" in
a x + b y + c z + d = epsilon
i.e
min epsilon' * epsilon
= min [a b c d] * [x y z ones()]' * [x y z ones()] * [a b c d]'
Let coor_cov = [x y z ones()]' * [x y z ones()], this is equivalant to
solve the minimum eigenvalue of coor_cov, which equals to min n*var(epsilon),
eigenvector is [a b c d]'. n is the number of points, e.g. length(x).
To avoid the extra "ones()", one need to fit a plane without constant term
a*(xmx) + b*(ymy) + c*(zmz) = 0
where the pivot [mx my mz] = mean([x y z]), then follow the same analysis.
Alternatively, [mx my mz] can be other point (presumably) in the plane,
although different choice introduce diffence bias. Usually "mean([x y z])"
been the most "fair" one.
Other eigenvalues indicate the spreading of the points in other dimensions in the plane.
(file #46498, file #46499)

New idea comes, try the attached patch, it is the method (2) below. The patch is primary for conceptual verification, no thoughtful test or careful cover of every edge case.
With this patch, the time cost for the initial reported sample is reduced to roughly 0.3 second, and the said everrunning example will end in 0.32 second.
Read a bit about graphics.cc and bug #47677.
The current algorithm there, for noncoplanar random data of n points, the time cost is roughly O(n^3): one n for a single coplanar test, one n for search the largest coplanar point set, last n for repeated search of each plane.
I can see at least two possible optimizations.
1. Use binary search instead of the full linear search.
This reduce the time cost to O(n^2 log(n)).
2. Incremental search for coplanar point set.
This is based on the observation that in the 'eig' method: eig(a'*a), the 3x3 matrix a'*a can be calculated by rank one updates corresponding to each point. Thus we can check if the newly added point is inside the plane or not in O(1) (one eig of 3x3 matrix). Time cost is thus reduced to O(n).
To further speed (if necessary), we may need to know the more frequent use case of patch():
a. large list of many non coplanar input. or
b. coplanar input.
If it is (b), probably adding an initial overall coplanar test can help, at the cost of slow down case (a). Current implementation favors (a).
(file #46431)

Not sure if it is proper to post comment when the bug is marked "Fixed"...
@Markus
I tried several coplanar test methods. They all seems slightly faster than your attached checker in comment 14 (file #46341).
N = 1e7;
[Q, ~] = qr(rand(3)); % random orthonormal matrix
vecs = rand(N, 2) * Q(1:2, :); % random coplanar samples
disp(['Normal vector = ', num2str(Q(3, :))]);
[x, y, z] = deal(vecs(:,1), vecs(:,2), vecs(:,3));
tic; is_coplanar (x, y, z); toc
tic
[copl, normal_vector, eig_vals] = is_coplanar_lsm (x, y, z, [], 'svd');
toc
tic
[copl, normal_vector, eig_vals] = is_coplanar_lsm (x, y, z, [], 'eig');
toc
tic
[copl, normal_vector, eig_vals] = is_coplanar_lsm (x, y, z, [], 'lsm');
toc
octavegui:>
Normal vector = 0.72582 0.062478 0.68504
Elapsed time is 0.977858 seconds.
Elapsed time is 0.88521 seconds.
Elapsed time is 0.280604 seconds.
Elapsed time is 0.636894 seconds.
See attachment is_coplanar_lsm.m. There are 3 testers.
 'svd': basically the one you are using in libinterp/corefcn/graphics.cc
 'eig': mathematically equivalent to 'svd' method (svd(a).^2 <=> eig(a'*a)). But instead of doing SVD of a very tall matrix, here only eigenvalue factorization of 3x3 matrix is needed.
 'lsm': find a planar by least square method, i.e. solve [a;b] in ax+by=z.
All of these methods have a clear geometric meaning of tolerance, namely deviation to the best possible planar. In terms of FLOPS, 'lsm' method should be the best, but in practice it is slower than 'eig'.
My feeling is, in order to get better performance, one need some heuristic or greedy algorithm. e.g. first guess the correct planar by a set of good points, then test if all other points are close enough to this planar.
(file #46348)

Attached is a first shot at a coplanarity check using the volume of tetrahedrons.
It is fairly efficient for random data:
>> N=1e7; x = rand (N, 1); y = rand(N, 1); z = rand(N, 1);
>> tic, is_coplanar (x, y, z), toc
ans = 0
Elapsed time is 1.28589 seconds.
And for coplanar (but not colinear) data:
>> N=1e7; x = y = z = 1:N; y = y.*(1).^y; z = z.^3;
>> tic, is_coplanar (x, y, z), toc
ans = 1
Elapsed time is 1.25889 seconds.
But not so efficient if all data is colinear:
>> N=1e5; x = y = z = 1:N;
>> tic, is_coplanar (x, y, z), toc
ans = 1
Elapsed time is 6.73134 seconds.
Is there an efficient way to check whether many 3d vectors are colinear?
Also what criteria could we use to define that two vectors are "sufficiently noncolinear" for the crossproduct?
The results with the matGeom function as a point of comparison:
>> N=1e3; x = rand (N, 1); y = rand(N, 1); z = rand(N, 1);
>> tic, isCoplanar (x, y, z), toc
error: out of memory or dimension too large for Octave's index type
error: called from
nchoosek at line 134 column 9
isCoplanar at line 66 column 5
>> N=1e2; x = rand (N, 1); y = rand(N, 1); z = rand(N, 1);
>> tic, isCoplanar (x, y, z), toc
ans = 0
Elapsed time is 0.337654 seconds.
(file #46341)

Mike Miller said "... you've found a performance regression in Octave 5. Can you file a bug report?"
The example to test the bug he gave is the one below.
[x, y, z] = peaks (32);
[x, y, z] = deal (x(:), y(:), z(:));
tic;
scatter3 (x, y, z);
toc
$ octave4.4.1 t.m
Elapsed time is 0.104946 seconds.
$ octave5.0.91 t.m
Elapsed time is 20.1961 seconds.
My example that had issues with plotting is below (it was never able to finish the plot)
Octave code converted:
%%========================
%%2D Pseudo Quaterion Rotation Program in converted from Scilab 6.01 to Octave 5.0.91
w1=111
w2=333
w3=777
%% No phase differences in this: add in phi1 and phi2 phase differences if desired to the psi and phi angles
n = [0:.001:45];
x1=50.*cos((n./360).*2.*pi.*w1).*(1.144.*n./100); %%w = mu cos theta, the scalar
x2=50.*cos((n./360).*2.*pi.*w2).*(1.144.*n./100); %%cos phi
y1=50.*cos((n./360).*2.*pi.*w3).*(1.144.*n./100); %%cos psi
y2=50.*sin((n./360).*2.*pi.*w1).*(1.144.*n./100); %%sin theta
z1=50.*sin((n./360).*2.*pi.*w2).*(1.144.*n./100); %sin phi
z2=50.*sin((n./360).*2.*pi.*w3).*(1.144.*n./100); %%sin psi
x = y2.*x2;
y=y2.*z1.*y1;
z = y2.*z1.*z2;
scatter3 (x, y, z); %color can be created based on z
xlabel('Xaxis')
ylabel('Yaxis')
zlabel('Zaxis')
grid on
az = 347; %90;
el = 37; %180;
view([az,el]);
