HP Forums - minimum spanning tree lua HPPL

Good morning Albert Chan, here is the program, which seems to work, but it gives wrong results. For instance:
Points: = {{3,4.5},{1,5},{5,6.5},{2.5,5},{7,3},{1.5,1},{4,2.5},{2,2},{8,4},{3,9},{4.5,0.5},{5,2}};
Edges: = {{3,4.5,1,5,3,4.5,5,6.5,3,4.5,2.5,5,3,4.5,7,3,3,4.5,1.5,1,3,4.5,4,2.5,3,4.5,2,2,3,4.5,8,4,3,4.5,3,9,3,4.5,4.5,0.5,3,4.5,5,2},{1,5,5,6.5,1,5,2.5,5,1,5,7,3,1,5,1.5,1,1,5,4,2.5,1,5,2,2,1,5,8,4,1,5,3,9,1,5,4.5,0.5,1,5,5,2},{5,6.5,2.5,5,5,6.5,7,3,5,6.5,1.5,1,5,6.5,4,2.5,5,6.5,2,2,5,6.5,8,4,5,6.5,3,9,5,6.5,4.5,0.5,5,6.5,5,2},{2.5,5,7,3,2.5,5,1.5,1,2.5,5,4,2.5,2.5,5,2,2,2.5,5,8,4,2.5,5,3,9,2.5,5,4.5,0.5,2.5,5,5,2},{7,3,1.5,1,7,3,4,2.5,7,3,2,2,7,3,8,4,7,3,3,9,7,3,4.5,0.5,7,3,5,2},{1.5,1,4,2.5,1.5,1,2,2,1.5,1,8,4,1.5,1,3,9,1.5,1,4.5,0.5,1.5,1,5,2},{4,2.5,2,2,4,2.5,8,4,4,2.5,3,9,4,2.5,4.5,0.5,4,2.5,5,2},{2,2,8,4,2,2,3,9,2,2,4.5,0.5,2,2,5,2},{8,4,3,9,8,4,4.5,0.5,8,4,5,2},{3,9,4.5,0.5,3,9,5,2},{4.5,0.5,5,2}};

To compare the wrong result to the right result, I would like to send you an image, but I cannot find the "attach image" option. So I'm sending you a message on the forum.
Thanks for your availability, greetings, Roberto.

Code:

#cas

dist(x1,y1,x2,y2):=

BEGIN

RETURN sqrt((ABS(x2-x1))^2+(ABS(y2-y1))^2);

END;

same(x1,y1,x2,y2,edg):=

BEGIN

FOR jj FROM 1 TO length(edg) DO

IF (x1==edg(jj)(1) AND y1==edg(jj)(2) AND x2==edg(jj)(3) AND y2==edg(jj)(4)) OR (x2==edg(jj)(1) AND y2==edg(jj)(2) AND x1==edg(jj)(3) AND y1==edg(jj)(4))

THEN

RETURN true;

END;

END;

RETURN false;

END;

doesexist(lst,item):=

BEGIN

FOR uu FROM 1 TO length(lst) DO

IF lst(uu)==item THEN

RETURN true;

END;

END;

RETURN false;

END;

mstRob(points,edges):=

BEGIN

verticesi:={points(1)};

tree:={};

WHILE length(verticesi) != length(points) DO

  lnn:=10000;

  temp1:={};

  temp2:={};

  FOR j FROM 1 TO length(verticesi) DO

    FOR pp FROM 1 TO length(points) DO

      IF doesexist(verticesi,points(pp))==false THEN

        IF dist(points(pp)(1),points(pp)(2),verticesi(j)(1),verticesi(j)(2))<lnn AND same(points(pp)(1),points(pp)(2),verticesi(j)(1),verticesi(j)(2),edges) THEN

          lnn:=dist(points(pp)(1),points(pp)(2),verticesi(j)(1),verticesi(j)(2));

          temp1:=points(pp);

          temp2:=verticesi(j);

        END;

      END;

    END;

  END;

  verticesi:=append(verticesi,temp1);

  tree:=append(tree,{temp1(1),temp1(2),temp2(1),temp2(2)});

END;

RETURN tree;

END;

#end

I apologize, but I submitted the wrong attachment. See the two new attachments. It seems that the algorithm I presented does not take into account the minimum sum, but limits itself to joining the points in ascending order.
I can't figure out where I went wrong, translating the algorithm written in LUA (see https://github.com/asundheim/lua_minimum...er/mst.lua).

It is more natural to use complex number to represent points: (x,y) ≡ x + y*i
List of valid edges has each edge with 2 points. (to make a segment)

With your example, we have:

Code:

pts := [ /* total 12 points */

(3,4.5),(1,5),(5,6.5),(2.5,5),(7,3),(1.5,1),

(4,2.5),(2,2),(8,4),(3,9),(4.5,0.5),(5,2)];

edgs := [ /* total 12*11/2 = 66 edges */

[(3,4.5),(1,5)],

[(3,4.5),(5,6.5)],

[(3,4.5),(2.5,5)],

[(3,4.5),(7,3)],

[(3,4.5),(1.5,1)],

[(3,4.5),(4,2.5)],

[(3,4.5),(2,2)],

[(3,4.5),(8,4)],

[(3,4.5),(3,9)],

[(3,4.5),(4.5,0.5)],

[(3,4.5),(5,2)],

[(1,5),(5,6.5)],

[(1,5),(2.5,5)],

[(1,5),(7,3)],

[(1,5),(1.5,1)],

[(1,5),(4,2.5)],

[(1,5),(2,2)],

[(1,5),(8,4)],

[(1,5),(3,9)],

[(1,5),(4.5,0.5)],

[(1,5),(5,2)],

[(5,6.5),(2.5,5)],

[(5,6.5),(7,3)],

[(5,6.5),(1.5,1)],

[(5,6.5),(4,2.5)],

[(5,6.5),(2,2)],

[(5,6.5),(8,4)],

[(5,6.5),(3,9)],

[(5,6.5),(4.5,0.5)],

[(5,6.5),(5,2)],

[(2.5,5),(7,3)],

[(2.5,5),(1.5,1)],

[(2.5,5),(4,2.5)],

[(2.5,5),(2,2)],

[(2.5,5),(8,4)],

[(2.5,5),(3,9)],

[(2.5,5),(4.5,0.5)],

[(2.5,5),(5,2)],

[(7,3),(1.5,1)],

[(7,3),(4,2.5)],

[(7,3),(2,2)],

[(7,3),(8,4)],

[(7,3),(3,9)],

[(7,3),(4.5,0.5)],

[(7,3),(5,2)],

[(1.5,1),(4,2.5)],

[(1.5,1),(2,2)],

[(1.5,1),(8,4)],

[(1.5,1),(3,9)],

[(1.5,1),(4.5,0.5)],

[(1.5,1),(5,2)],

[(4,2.5),(2,2)],

[(4,2.5),(8,4)],

[(4,2.5),(3,9)],

[(4,2.5),(4.5,0.5)],

[(4,2.5),(5,2)],

[(2,2),(8,4)],

[(2,2),(3,9)],

[(2,2),(4.5,0.5)],

[(2,2),(5,2)],

[(8,4),(3,9)],

[(8,4),(4.5,0.5)],

[(8,4),(5,2)],

[(3,9),(4.5,0.5)],

[(3,9),(5,2)],

[(4.5,0.5),(5,2)]];

This is the PPL code for mstTree(points, edges)

Code:

#cas

mstValid(edges, edge) := 

  contains(edges,edge) OR 

  contains(edges,reverse(edge));

mstTree(points,edges):=

BEGIN

LOCAL vertices, tree, j, k, p, v, d1, d2, pointsd, n, temp;

vertices, tree := {points[1]}, {}

n, pointsd := 1, SIZE(points);  // n = SIZE(vertices)

WHILE n < pointsd DO

  d1, temp := inf, (inf, inf)

  FOR j FROM 1 TO n DO

    v := vertices[j];

    FOR k FROM 1 TO pointsd DO

      p := points[k];

      IF contains(vertices, p) THEN continue END;      

      IF (d2 := abs(p-v)) >= d1 THEN continue END;

      IF NOT mstValid(edges, [p,v]) THEN continue END;

      d1, temp := d2, [p,v]

    END

  END

  tree[n] := temp;

  vertices[n+=1] := temp[1];

END

return tree;

END;

#end

Input pts, and edgs into HP Prime emulator, we get

CAS> r := mstTree(pts,edgs)

\(\left(\begin{array}{cc}
2.5+5.0*i & 3.0+4.5*i \\
1.0+5.0*i & 2.5+5.0*i \\
4.0+2.5*i & 3.0+4.5*i \\
5.0+2.0*i & 4.0+2.5*i \\
4.5+0.5*i & 5.0+2.0*i \\
2.0+2.0*i & 4.0+2.5*i \\
1.5+i & 2.0+2.0*i \\
7.0+3.0*i & 5.0+2.0*i \\
8.0+4.0*i & 7.0+3.0*i \\
5.0+6.5*i & 3.0+4.5*i \\
3.0+9.0*i & 5.0+6.5*i
\end{array}\right) \)

This matched the results of Lua code.
When plotted line segments, it matched your supplied image.

CAS> apply((x)->segment(x[0],x[1]),r)

Perhaps we should let the code make the valid edges ?
This way mstValid(edges, [p,v]) can be removed, with boost in speed.

(10-28-2021 04:39 PM)Albert Chan Wrote: [ -> ]It is more natural to use complex number to represent points: (x,y) ≡ x + y*i
List of valid edges has each edge with 2 points. (to make a segment)

With your example, we have:

Code:

pts := [ /* total 12 points */ (3,4.5),(1,5),(5,6.5),(2.5,5),(7,3),(1.5,1), (4,2.5),(2,2),(8,4),(3,9),(4.5,0.5),(5,2)]; edgs := [ /* total 12*11/2 = 66 edges */ [(3,4.5),(1,5)], [(3,4.5),(5,6.5)], [(3,4.5),(2.5,5)], [(3,4.5),(7,3)], [(3,4.5),(1.5,1)], [(3,4.5),(4,2.5)], [(3,4.5),(2,2)], [(3,4.5),(8,4)], [(3,4.5),(3,9)], [(3,4.5),(4.5,0.5)], [(3,4.5),(5,2)], [(1,5),(5,6.5)], [(1,5),(2.5,5)], [(1,5),(7,3)], [(1,5),(1.5,1)], [(1,5),(4,2.5)], [(1,5),(2,2)], [(1,5),(8,4)], [(1,5),(3,9)], [(1,5),(4.5,0.5)], [(1,5),(5,2)], [(5,6.5),(2.5,5)], [(5,6.5),(7,3)], [(5,6.5),(1.5,1)], [(5,6.5),(4,2.5)], [(5,6.5),(2,2)], [(5,6.5),(8,4)], [(5,6.5),(3,9)], [(5,6.5),(4.5,0.5)], [(5,6.5),(5,2)], [(2.5,5),(7,3)], [(2.5,5),(1.5,1)], [(2.5,5),(4,2.5)], [(2.5,5),(2,2)], [(2.5,5),(8,4)], [(2.5,5),(3,9)], [(2.5,5),(4.5,0.5)], [(2.5,5),(5,2)], [(7,3),(1.5,1)], [(7,3),(4,2.5)], [(7,3),(2,2)], [(7,3),(8,4)], [(7,3),(3,9)], [(7,3),(4.5,0.5)], [(7,3),(5,2)], [(1.5,1),(4,2.5)], [(1.5,1),(2,2)], [(1.5,1),(8,4)], [(1.5,1),(3,9)], [(1.5,1),(4.5,0.5)], [(1.5,1),(5,2)], [(4,2.5),(2,2)], [(4,2.5),(8,4)], [(4,2.5),(3,9)], [(4,2.5),(4.5,0.5)], [(4,2.5),(5,2)], [(2,2),(8,4)], [(2,2),(3,9)], [(2,2),(4.5,0.5)], [(2,2),(5,2)], [(8,4),(3,9)], [(8,4),(4.5,0.5)], [(8,4),(5,2)], [(3,9),(4.5,0.5)], [(3,9),(5,2)], [(4.5,0.5),(5,2)]];

This is the PPL code for mstTree(points, edges)

Code:

#cas mstValid(edges, edge) := contains(edges,edge) OR contains(edges,reverse(edge)); mstTree(points,edges):= BEGIN LOCAL vertices, tree, j, k, p, v, d1, d2, pointsd, n, temp; vertices, tree := {points[1]}, {} n, pointsd := 1, SIZE(points); // n = SIZE(vertices) WHILE n < pointsd DO d1, temp := inf, (inf, inf) FOR j FROM 1 TO n DO v := vertices[j]; FOR k FROM 1 TO pointsd DO p := points[k]; IF contains(vertices, p) THEN continue END; IF (d2 := abs(p-v)) >= d1 THEN continue END; IF NOT mstValid(edges, [p,v]) THEN continue END; d1, temp := d2, [p,v] END END tree[n] := temp; vertices[n+=1] := temp[1]; END return tree; END; #end

Input pts, and edgs into HP Prime emulator, we get

CAS> r := mstTree(pts,edgs)

\(\left(\begin{array}{cc}
2.5+5.0*i & 3.0+4.5*i \\
1.0+5.0*i & 2.5+5.0*i \\
4.0+2.5*i & 3.0+4.5*i \\
5.0+2.0*i & 4.0+2.5*i \\
4.5+0.5*i & 5.0+2.0*i \\
2.0+2.0*i & 4.0+2.5*i \\
1.5+i & 2.0+2.0*i \\
7.0+3.0*i & 5.0+2.0*i \\
8.0+4.0*i & 7.0+3.0*i \\
5.0+6.5*i & 3.0+4.5*i \\
3.0+9.0*i & 5.0+6.5*i
\end{array}\right) \)

This matched the results of Lua code.
When plotted line segments, it matched your supplied image.

CAS> apply((x)->segment(x[0],x[1]),r)

Perhaps we should let the code make the valid edges ?
This way mstValid(edges, [p,v]) can be removed, with boost in speed.

Congratulations, your solution is truly "elegant". However, I modified my code (see below), modifying the "same" statement, and the result obtained is like the one you proposed just now. I will save your code to compare it with my code.
Thanks a lot, Roberto.
Right code:

Code:

same(x1,y1,x2,y2,edg):=

BEGIN

FOR jj FROM 1 TO length(edg) DO

IF (x1≠edg(jj)(1) AND y1≠edg(jj)(2) AND x2≠edg(jj)(3) AND y2≠edg(jj)(4)) OR (x2≠edg(jj)(1) AND y2≠edg(jj)(2) AND x1≠edg(jj)(3) AND y1≠edg(jj)(4)) <----

THEN

RETURN true;

END;

END;

RETURN false;

END;

Wrong code:

Code:

same(x1,y1,x2,y2,edg):=

BEGIN

FOR jj FROM 1 TO length(edg) DO

IF (x1==edg(jj)(1) AND y1==edg(jj)(2) AND x2==edg(jj)(3) AND y2==edg(jj)(4)) OR (x2==edg(jj)(1) AND y2==edg(jj)(2) AND x1==edg(jj)(3) AND y1==edg(jj)(4)) <----

THEN

RETURN true;

END;

END;

RETURN false;

END;

This version shuffled the list of points to vertices, and reduce loop counts by half.

Code:

#cas

mstValid(edges, edge) := 

  contains(edges,edge) OR 

  contains(edges,reverse(edge));

mstTree(points,edges):=

BEGIN

LOCAL tree, j, k, p, v, d1, d2, pointsd, n, temp;

tree := {}

n, pointsd := 1, SIZE(points);  // n = SIZE(vertices)

WHILE n < pointsd DO

  d1, temp := inf, (inf, inf)

  FOR j FROM 1 TO n DO

    v := points[j];             // filled vertices

    FOR k FROM n+1 TO pointsd DO

      p := points[k];

      IF (d2 := abs(p-v)) >= d1 THEN continue END;

      IF NOT mstValid(edges, [p,v]) THEN continue END;

      d1, temp := d2, [k,p,v]

    END

  END

  [k,p,v] := temp

  tree[n] := [p,v];

  points[k] := points[n+=1];    // make space

  points[n] := p;               // save vertice

END

return tree;

END;

#end

(10-28-2021 05:31 PM)Albert Chan Wrote: [ -> ]This version shuffled the list of points to vertices, and reduce loop counts by half.

Code:

#cas mstValid(edges, edge) := contains(edges,edge) OR contains(edges,reverse(edge)); mstTree(points,edges):= BEGIN LOCAL tree, j, k, p, v, d1, d2, pointsd, n, temp; tree := {} n, pointsd := 1, SIZE(points); // n = SIZE(vertices) WHILE n < pointsd DO d1, temp := inf, (inf, inf) FOR j FROM 1 TO n DO v := points[j]; // filled vertices FOR k FROM n+1 TO pointsd DO p := points[k]; IF (d2 := abs(p-v)) >= d1 THEN continue END; IF NOT mstValid(edges, [p,v]) THEN continue END; d1, temp := d2, [k,p,v] END END [k,p,v] := temp tree[n] := [p,v]; points[k] := points[n+=1]; // make space points[n] := p; // save vertice END return tree; END; #end

Very good!
I ask a question: what if the points in the sample are not represented by a pair of numbers, but by three or four numbers, as in multivariate statistics? For instance:
points: = {{15,17,24,14}, {17,15,32,26}, {15,14,29,23}, {13,12,10,16}, {20,17,26,28}, {15,21,26,21}}.

edges:={{15,17,24,14,17,15,32,26,15,17,24,14,15,14,29,23,15,17,24,14,13,12,10,16,15,17,24,14,20,17,26,28,15,17,24,14,15,21,26,21},{17,15,32,26,15,14,29,23,17,15,32,26,13,12,10,16,17,15,32,26,20,17,26,28,17,15,32,26,15,21,26,21},{15,14,29,23,13,12,10,16,15,14,29,23,20,17,26,28,15,14,29,23,15,21,26,21},{13,12,10,16,20,17,26,28,13,12,10,16,15,21,26,21},{20,17,26,28,15,21,26,21}}.

With these points and edges, "my" code results in:

[[17,15,15,17],[15,14,17,15],[13,12,15,14],[20,17,17,15],[15,21,15,17]].

With these points and edges, what does your code provide, Albert?

(10-28-2021 05:04 PM)robmio Wrote: [ -> ]However, I modified my code (see below), modifying the "same" statement ...

No ! OP post were correct ! (except for minor leaking of variables, jj in same())

It was your inputs that were wrong.
Each point has 2 numbers, Edge have 2 points, thus 4 numbers.
Remember, each edge contains 4 values.

pts0 := [ [3,4.5], [1,5], ..., ,[5,2]]; /* 12 points */
edgs0 := [ [3,4.5 , 1,5], [3,4.5 , 5,6.5], ..., [4.5,0.5 , 5,2]]; /* 66 edges */

Running OP code, on my laptop, it finished in 0.273 sec, return these line segments.
(for comparison, my complex-number as points 2nd version, it finished in 0.006 sec)

CAS> mstRob(pts0, edgs0)
{{2.5,5,3,4.5},
{1,5,2.5,5},
{4,2.5,3,4.5},
{5,2,4,2.5},
{4.5,0.5,5,2},
{2,2,4,2.5},
{1.5,1,2,2},
{7,3,5,2},
{8,4,7,3},
{5,6.5,3,4.5},
{3,9,5,6.5}}

Version 3, extending to multi-dimensional points.

List of possible edges may be huge, this version switched the logic, and disallowed "bad" edges.
Example, if edge [p1,p2] (and of course, [p2, p1]) is not allowed: badEdges = [[p1,p2], [p2,p1]]
Note that badEdges list length always even.

For comparing "distance", it uses norm2, squared, of difference of 2 points.

Without doing square roots, this version is even faster than my previous version.
Running same example, 12 points, 66 edges (badEdges = []) , it finished in 0.004 sec.

Code:

#cas  

mstTree(points, badEdges):=

BEGIN

LOCAL tree, j, k, p, v, d1, d2, pointsd, n, temp;

tree := {}

n, pointsd := 1, len(points);   // n = len(vertices)

WHILE n < pointsd DO

  d1, temp := inf, (inf, inf)

  FOR j FROM 1 TO n DO

    v := points[j];             // filled vertices

    FOR k FROM n+1 TO pointsd DO

      p := points[k];

      d2 = p .- v;

      IF (d2 *= d2) >= d1 THEN continue END;

      IF member([k,j], badEdges) THEN continue END;

      d1, temp := d2, [k,p,v];

    END

  END

  [k,p,v] := temp;

  tree[n] := [p,v];

  points[k] := points[n+=1];    // make space

  points[n] := p;               // save vertice

END

return tree;

END;

#end

With your example, 6 points, 6*5/2=15 edges (again, badEdges = []), this is the result:

CAS> pts := [[15,17,24,14], [17,15,32,26], [15,14,29,23], [13,12,10,16], [20,17,26,28], [15,21,26,21]]
CAS> mstTree(pts, [])

{[[15,21,26,21] , [15,17,24,14]],
[[15,14,29,23] , [15,21,26,21]],
[[17,15,32,26] , [15,14,29,23]],
[[20,17,26,28] , [17,15,32,26]],
[[13,12,10,16] , [15,17,24,14]]}

sum of line segments length = √69 + √62 + √23 + √53 + √229 ≈ 43.3893

(10-28-2021 09:07 PM)Albert Chan Wrote: [ -> ]Version 3, extending to multi-dimensional points.

List of possible edges may be huge, this version switched the logic, and disallowed "bad" edges.
Example, if edge [p1,p2] (and of course, [p2, p1]) is not allowed: badEdges = [[p1,p2], [p2,p1]]
Note that badEdges list length always even.

For comparing "distance", it uses norm2, squared, of difference of 2 points.

Without doing square roots, this version is even faster than my previous version.
Running same example, 12 points, 66 edges (badEdges = []) , it finished in 0.004 sec.

Code:

#cas mstTree(points, badEdges):= BEGIN LOCAL tree, j, k, p, v, d1, d2, pointsd, n, temp; tree := {} n, pointsd := 1, len(points); // n = len(vertices) WHILE n < pointsd DO d1, temp := inf, (inf, inf) FOR j FROM 1 TO n DO v := points[j]; // filled vertices FOR k FROM n+1 TO pointsd DO p := points[k]; d2 = p .- v; IF (d2 *= d2) >= d1 THEN continue END; IF member([k,j], badEdges) THEN continue END; d1, temp := d2, [k,p,v]; END END [k,p,v] := temp; tree[n] := [p,v]; points[k] := points[n+=1]; // make space points[n] := p; // save vertice END return tree; END; #end

With your example, 6 points, 6*5/2=15 edges (again, badEdges = []), this is the result:

CAS> pts := [[15,17,24,14], [17,15,32,26], [15,14,29,23], [13,12,10,16], [20,17,26,28], [15,21,26,21]]
CAS> mstTree(pts, [])

{[[15,21,26,21] , [15,17,24,14]],
[[15,14,29,23] , [15,21,26,21]],
[[17,15,32,26] , [15,14,29,23]],
[[20,17,26,28] , [17,15,32,26]],
[[13,12,10,16] , [15,17,24,14]]}

sum of line segments length = √69 + √62 + √23 + √53 + √229 ≈ 43.3893

I understood: my mistake was to propose the edges to the code in this way:

{{3,4.5,1,5,3,4.5,5,6.5,3,4.5,2.5,5,3,4.5,7,3,3,4.5,1.5,1,3,4.5,4,2.5,3,4.5,2,2,3,4.5,8,4,3,4.5,3,9,3,4.5,4.5,0.5,3,4.5,5,2},{1,5,5,6.5,1,5,2.5,5,1,5,7,3,1,5,1.5,1,1,5,4,2.5,1,5,2,2,1,5,8,4,1,5,3,9,1,5,4.5,0.5,1,5,5,2},{5,6.5,2.5,5,5,6.5,7,3,5,6.5,1.5,1,5,6.5,4,2.5,5,6.5,2,2,5,6.5,8,4,5,6.5,3,9,5,6.5,4.5,0.5,5,6.5,5,2},{2.5,5,7,3,2.5,5,1.5,1,2.5,5,4,2.5,2.5,5,2,2,2.5,5,8,4,2.5,5,3,9,2.5,5,4.5,0.5,2.5,5,5,2},{7,3,1.5,1,7,3,4,2.5,7,3,2,2,7,3,8,4,7,3,3,9,7,3,4.5,0.5,7,3,5,2},{1.5,1,4,2.5,1.5,1,2,2,1.5,1,8,4,1.5,1,3,9,1.5,1,4.5,0.5,1.5,1,5,2},{4,2.5,2,2,4,2.5,8,4,4,2.5,3,9,4,2.5,4.5,0.5,4,2.5,5,2},{2,2,8,4,2,2,3,9,2,2,4.5,0.5,2,2,5,2},{8,4,3,9,8,4,4.5,0.5,8,4,5,2},{3,9,4.5,0.5,3,9,5,2},{4.5,0.5,5,2}}

without making groups of 4 numbers into 4 numbers.
Thank you for having brilliantly solved the "multivariate" version as well. Are you a computer scientist?

Best wishes, Roberto.

(10-29-2021 05:53 AM)robmio Wrote: [ -> ]Thank you for having brilliantly solved the "multivariate" version as well. Are you a computer scientist?

No. I had never heard of minimum spanning tree Big Grin

(10-28-2021 09:07 PM)Albert Chan Wrote: [ -> ]Version 3, extending to multi-dimensional points
...

CAS> pts := [[15,17,24,14], [17,15,32,26], [15,14,29,23], [13,12,10,16], [20,17,26,28], [15,21,26,21]]
CAS> mstTree(pts, [])

{[[15,21,26,21] , [15,17,24,14]],
[[15,14,29,23] , [15,21,26,21]],
[[17,15,32,26] , [15,14,29,23]],
[[20,17,26,28] , [17,15,32,26]],
[[13,12,10,16] , [15,17,24,14]]}

sum of line segments length = √69 + √62 + √23 + √53 + √229 ≈ 43.3893

Instead of badEdges = [[p1,p2],[p2,p1], ...], we could simplify to [[1,2],[2,1], ...]
Same for tree = [[p(j1),p(k1)], [p(j2),p(k2)], ...], simplified to [[j1,k1],[j2,k2], ...]

But, this required pts not be shuffled (otherwise indexes meant nothing)
Version 4 solve the problems with another level of indirection (*)

Code:

#cas

mstTree(points, badEdges):=

BEGIN

LOCAL tree, j, k, v, d1, d2, pointsd, n, temp;

LOCAL m, mj, mk;                // mapped indexes

n, pointsd := 1, len(points);   // n = len(vertices)

tree, m := {}, range(1,pointsd+1);

WHILE n < pointsd DO

  d1 := temp := inf;

  FOR j FROM 1 TO n DO

    v := points[mj := m[j]];    // filled vertices

    FOR k FROM n+1 TO pointsd DO

      d2 = v .- points[mk := m[k]];

      IF (d2 *= d2) >= d1 THEN continue END;            

      IF member([mk,mj], badEdges) THEN continue END;

      d1, temp := d2, [k, [mk,mj]];

    END

  END

  k, tree[n] := temp;

  m[k] := m[n += 1];            // swap m[k], m[n+1]

  m[n] := temp[2][1];

END

return tree;

END;

#end

Redo the same example, with this new version

CAS> pts := [[15,17,24,14], [17,15,32,26], [15,14,29,23], [13,12,10,16], [20,17,26,28], [15,21,26,21]]
CAS> tree := mstTree(pts, [])

[[6,1],[3,6],[2,3],[5,2],[4,1]]

CAS> norm2(x) := sqrt(dot(x,x))
CAS> map(tree, x -> norm2(pts[x[1]] - pts[x[2]])) → [√69,√62,√23,√53,√229]
CAS> sum(float(Ans)) → 43.3893190999

What if points #2, #3 line segment not allowed ?

CAS> tree := mstTree(pts, [[2,3],[3,2]])

[[6,1],[3,6],[5,3],[2,5],[4,1]]

CAS> map(tree, x -> norm2(pts[x[1]] - pts[x[2]])) → [√69,√62,2*√17,√53,√229]
CAS> sum(float(Ans)) → 46.8396988279

(*) [Image: quote-all-problems-in-computer-science-c...4-0454.jpg]

[Image: quote-all-problems-in-computer-science-c...4-0454.jpg]

Perfect, the substitution (in the output) of the coordinates with the numbers (referring to each node), will allow me to program the Friedman-Rafsky test.

Again thank you for solving the problem brilliantly.
Best regards, Roberto.

I address the question to Albert Chan: using the last algorithm you wrote, for a matrix of dimensions {100,4} the "physical" calculator PRIME_G2 takes about 13 seconds to calculate the "MST". Does your algorithm already include the implementation with the Fibonacci heap?

Version 3 assumed all edges vaild, unless specified otherwise.
I think the speed is due to not spend time to validate edges.

By shuffling mapped indexes, it keep the points to examine low.
Points to examine reduced to about 1/3 (and, no code to validate vertices)

Note: test for member in list has O(n)

---

Points examined, for oringinal Lua code, n = total points, k = vertices examined

sum(n*k, k=1..n-1)
= n * n*(n-1)/2
= n^3/2 - n^2/2

Points examined, with shuffled indexes, splitting off vertices and non-vertices.

sum((n-k)*k, k=1..n-1)
= sum((n-1)*k - k*(k-1), k=1..n-1) // see Funny Factorials and Slick Sums
= (n-1) * n*(n-1)/2 - n*(n-1)*(n-2)/3
= n*(n-1)*(n+1)/6
= n^3/6 - n/6

(11-02-2021 11:57 AM)Albert Chan Wrote: [ -> ]Version 3 assumed all edges vaild, unless specified otherwise.
I think the speed is due to not spend time to validate edges.

By shuffling mapped indexes, it keep the points to examine low.
Points to examine reduced to about 1/3 (and, no code to validate vertices)

Note: test for member in list has O(n)

---

Points examined, for oringinal Lua code, n = total points, k = vertices examined

sum(n*k, k=1..n-1)
= n * n*(n-1)/2
= n^3/2 - n^2/2

Points examined, with shuffled indexes, splitting off vertices and non-vertices.

sum((n-k)*k, k=1..n-1)
= sum((n-1)*k - k*(k-1), k=1..n-1) // see Funny Factorials and Slick Sums
= (n-1) * n*(n-1)/2 - n*(n-1)*(n-2)/3
= n*(n-1)*(n+1)/6
= n^3/6 - n/6

Good. I'll try to translate your algorithm into Python, to see if it gets even faster. Thanks for your help.

(10-29-2021 10:52 AM)Albert Chan Wrote: [ -> ]Instead of badEdges = [[p1,p2],[p2,p1], ...], we could simplify to [[1,2],[2,1], ...]

I thought it is trivial to supply 1 sided badEdges, then "doubled" them.
For python, joining two list is simply '+'

>>> bad = [[1,2],[3,6],[6,10]]
>>> bad2 = [ [k[1],k[0]] for k in bad]
>>> bad + bad2
[[1, 2], [3, 6], [6, 10], [2, 1], [6, 3], [10, 6]]

For HP Prime, close equivalent for '+' is extend.
But, instead of extending to end of list, it extend to items inside list.

Cas> bad := [[1,2],[3,6],[6,10]]
Cas> bad2 := map(bad, reverse) → [[2,1],[6,3],[10,6]]
Cas> extend(bad, bad2) → [[1,2,2,1],[3,6,6,3],[6,10,10,6]]

To get what we wanted, we need to use transpose (3 times !)
This work, but ugly. Better solution welcome.

Cas> transpose(extend(transpose(bad), transpose(bad2))
[[1,2],[3,6],[6,10],[2,1],[6,3],[10,6]]

Version 5:

Code:

#cas

mstTree(points, badEdges):=

BEGIN

LOCAL tree, j, k, v, d1, d2, pointsd, n, temp;

LOCAL m, mj, mk;                // mapped indexes

n, pointsd := 1, len(points);   // n = len(vertices)

tree, m := {}, range(1,pointsd+1);

WHILE n < pointsd DO

  d1 := temp := inf;

  FOR j FROM 1 TO n DO

    v := points[mj := m[j]];    // filled vertices

    FOR k FROM n+1 TO pointsd DO

      d2 = v .- points[mk := m[k]];

      IF (d2 *= d2) >= d1 THEN continue END;            

      IF member([mj,mk], badEdges) THEN continue END;

      d1, temp := d2, [k, [mj,mk]];

    END

  END

  k, tree[n] =< temp;

  m[k] =< m[n += 1];            // swap m[k], m[n+1]

  m[n] =< temp[2][2];

END

return tree;

END;

#end

Since it is more natural to read from left to right, I flip the order.
Also, HP Prime always create a new list, even for just updating element.
We might as well update list element by reference, '=<' instead of ':='

OP example:

Cas> pts := [[3,4.5],[1,5],[5,6.5],[2.5,5],[7,3],[1.5,1],[4,2.5],[2,2],[8,4],[3,9],[4.5,0.5],[5,2]]
Cas> mstTree(pts, [])
[[1,4],[4,2],[1,7],[7,12],[12,11],[7,8],[8,6],[12,5],[5,9],[1,3],[3,10]]

It is easier to sketch the tree.

3 10
|
1 4 2
|
7 12 11
| |
| 5 9
|
8 6

(11-04-2021 03:47 PM)Albert Chan Wrote: [ -> ]
(10-29-2021 10:52 AM)Albert Chan Wrote: [ -> ]Instead of badEdges = [[p1,p2],[p2,p1], ...], we could simplify to [[1,2],[2,1], ...]

I thought it is trivial to supply 1 sided badEdges, then "doubled" them.
For python, joining two list is simply '+'

>>> bad = [[1,2],[3,6],[6,10]]
>>> bad2 = [ [k[1],k[0]] for k in bad]
>>> bad + bad2
[[1, 2], [3, 6], [6, 10], [2, 1], [6, 3], [10, 6]]

For HP Prime, close equivalent for '+' is extend.
But, instead of extending to end of list, it extend to items inside list.

Cas> bad := [[1,2],[3,6],[6,10]]
Cas> bad2 := map(bad, revserse) → [[2,1],[6,3],[10,6]]
Cas> extend(bad, bad2) → [[1,2,2,1],[3,6,6,3],[6,10,10,6]]

To get what we wanted, we need to use transpose (3 times !)
This work, but ugly. Better solution welcome.

Cas> transpose(extend(transpose(bad), transpose(bad2))
[[1,2],[3,6],[6,10],[2,1],[6,3],[10,6]]

Version 5:

Code:

#cas mstTree(points, badEdges):= BEGIN LOCAL tree, j, k, v, d1, d2, pointsd, n, temp; LOCAL m, mj, mk; // mapped indexes n, pointsd := 1, len(points); // n = len(vertices) tree, m := {}, range(1,pointsd+1); WHILE n < pointsd DO d1 := temp := inf; FOR j FROM 1 TO n DO v := points[mj := m[j]]; // filled vertices FOR k FROM n+1 TO pointsd DO d2 = v .- points[mk := m[k]]; IF (d2 *= d2) >= d1 THEN continue END; IF member([mj,mk], badEdges) THEN continue END; d1, temp := d2, [k, [mj,mk]]; END END k, tree[n] =< temp; m[k] =< m[n += 1]; // swap m[k], m[n+1] m[n] =< temp[2][2]; END return tree; END; #end

Since it is more natural to read from left to right, I flip the order.
Also, HP Prime always create a new list, even for just updating element.
We might as well update list element by reference, '=<' instead of ':='

OP example:

Cas> pts := [[3,4.5],[1,5],[5,6.5],[2.5,5],[7,3],[1.5,1],[4,2.5],[2,2],[8,4],[3,9],[4.5,0.5],[5,2]]
Cas> mstTree(pts, [])
[[1,4],[4,2],[1,7],[7,12],[12,11],[7,8],[8,6],[12,5],[5,9],[1,3],[3,10]]

It is easier to sketch the tree.

3 10
|
1 4 2
|
7 12 11
| |
| 5 9
|
8 6

Indeed, it is easier to view the entire tree

It seems XCas already has "safe" extend built-in, called semi-augment

semi_augment(A,B) := when(len(A)!=len(B), extend(A,B), extend(append(A,head(B)),tail(B)))

Cas> semi_augment([[1,2],[3,6],[6,10]] , [[2,1],[6,3],[10,6]])

[[1,2], [3,6], [6,10], [2,1], [6,3], [10,6]] // 3+3 = 6 items

Cas> semi_augment([[1,2]] , [[3,4,5]])

[[1,2], [3,4,5]] // 1 + 1 = 2 items

Had we use extend instead, we get these:

Cas> extend([[1,2],[3,6],[6,10]] , [[2,1],[6,3],[10,6]])

[[1,2,2,1], [3,6,6,3], [6,10,10,6]]

Cas> extend([[1,2]] , [[3,4,5]])

[[1,2,3,4,5]]