Term_factorized/3 error for compound_name_arity zero term

kuniaki.mukai · September 7, 2026, 3:13pm

I happened to find a bug (?) on compound_name_arity zero terms in ‘term_facotized/3’.

% ?- term_factorized(f(q()), S, E).
%@ S = f(q()),
%@ E = [].

% ?- term_factorized(f(q(), q()), S, E).
%@ ERROR: Domain error: `compound_non_zero_arity' expected, found `q()'
%@ ERROR: In:
%@ ERROR:   [14] functor(q(),_7688,_7690)
%@ ERROR:   [13] terms:mk_subst([... - _7738],[_7748=_7750|_7744],t(black('',_7762,_7764,''),black(...,...,_7774,...))) at /Users/cantor/lib/swipl/library/terms.pl:309
%@ ERROR:   [11] toplevel_call(pt_mincoa:pt_mincoa: ...) at /Users/cantor/lib/swipl/boot/toplevel.pl:1678
%@ ERROR:

I wanted to compare my version of ‘term_factoirzed’/3 via minimum coalgebra
representation, which shows statistics below ( imac intel 2020, macOS 26 Swiprolog latest ).
I hope it would be marginal compared with the builtin ‘term_foctorized/3’.

% ?-N=100000, time(forall(between(1, N, _), pt_term_factorize_test)).
%@ % 30,222,332 inferences, 2.885 CPU in 2.890 seconds (100% CPU, 10475896 Lips)
%@ N = 100000.

pt_term_factorize_test :-
	random_term(R),
	pt_term_coa([R], [I], Coa),
	pt_term_factorize(I, Coa, _, _).

jan · September 7, 2026, 3:58pm

Yes Fixed. Curious to how your code compares!

kuniaki.mukai · September 7, 2026, 4:32pm

Thanks for fix.

% ?- term_factorized(f(q()), S, E).
%@ S = f(q()),
%@ E = [].
% ?- term_factorized(f(q(), q()), S, E).
%@ S = f(_A, _A),
%@ E = [_A=q()].

However, it still ssems sometning wrong. Simply I use this predicate for test query below,
and got errors. So I hope it is reproducible for large N .

term_factorize_test :-
	random_term(R),
	term_factorized(R, _, _).

% ?-N=1, time(forall(between(1, N, _), term_factorize_test)).
%@ % 87 inferences, 0.000 CPU in 0.000 seconds (83% CPU, 2485714 Lips)
%@ N = 1.
% ?-N=2, time(forall(between(1, N, _), term_factorize_test)).
%@ % 869 inferences, 0.000 CPU in 0.000 seconds (94% CPU, 6437037 Lips)
%@ N = 2.
% ?-N=10, time(forall(between(1, N, _), term_factorize_test)).
%@ % 2,451 inferences, 0.000 CPU in 0.000 seconds (97% CPU, 9111524 Lips)
%@ N = 10.
% ?-N=100, time(forall(between(1, N, _), term_factorize_test)).
%@ % 5,356 inferences, 0.001 CPU in 0.001 seconds (97% CPU, 10260536 Lips)
%@ ERROR: Domain error: `compound_non_zero_arity' expected, found `h()'
%@ ERROR: In:
%@ ERROR:   [17] throw(error(domain_error(compound_non_zero_arity,...),context(...,_1322)))
%@ ERROR:   [15] catch('<garbage_collected>','<garbage_collected>','<garbage_collected>') at /Users/cantor/lib/swipl/boot/init.pl:565
%@ ERROR:   [14] call_cleanup(prolog_statistics:catch(...,...,...),prolog_statistics:(_1400=true)) at /Users/cantor/lib/swipl/boot/init.pl:685
%@ ERROR:   [13] prolog_statistics:time('<garbage_collected>') at /Users/cantor/lib/swipl/library/statistics.pl:287
%@ ERROR:   [11] toplevel_call('<garbage_collected>') at /Users/cantor/lib/swipl/boot/toplevel.pl:1678
%@ ERROR: 
%@ ERROR: Note: some frames are missing due to last-call optimization.
%@ ERROR: Re-run your program in debug mode (:- debug.) to get more detail.
%@    Call: (16) _192=true ? creep
%@    Exit: (16) true=true ? abort

jan · September 8, 2026, 7:04am

Fixed that too. You now get occasional stack overflows. This is the result of using an rb tree for the lookup. This depends on proper standard order, but as we have seen before, this is not guaranteed for (all) cyclic terms.

kuniaki.mukai · September 8, 2026, 12:20pm

Thanks for quick fix. Now I can compare the builtin ‘term_factorized/3’
and ‘pt_term_factorized/3’ of mine. pt is for Paige-Tarjan method
for term graphs. My version uses neither rbtree,
structural unification, nor structural comparison to collect subterms
for feeding paige-tarjan method. As far as this quick test for shallow terms
is concerned my version is much stable and robust.
The problem if exists, C-codes are generated by AI, and
I am now not familiar with C-codes at all. I think
if you asked the AI, much better C-codes must be returned in a blink.

term_factorized_test :-
	random_term(R),
	term_factorized(R, _, _).

% ?- N=1000, time(forall(between(1, N, _), term_factorized_test)).
%@ % 103,882,118 inferences, 2.884 CPU in 3.104 seconds (93% CPU, 36016355 Lips)
%@ ERROR: Stack limit (1.0Gb) exceeded
%@ ERROR:   Stack sizes: local: 0.7Gb, global: 0.3Gb, trail: 5Kb
%@ ERROR:   Stack depth: 12,197,990, last-call: 50%, Choice points: 11
%@ ERROR:   In:
%@ ERROR:     [12,197,990] rbtrees:lookup(<, <compound f/3>, _1020, <compound black/4>)
%@ ERROR:     [12,197,987] terms:insert_vars(<compound f/3>, _1060, <compound t/2>)
%@ ERROR:     [12,197,986] terms:insert_arg_vars(1, <compound f/3>, <compound f/3>, <compound t/2>)
%@ ERROR:     [12,197,984] terms:insert_arg_vars(1, <compound f/3>, <compound f/3>, <compound t/2>)
%@ ERROR:     [12,197,982] terms:insert_arg_vars(1, <compound f/3>, <compound f/3>, <compound t/2>)
%@ ERROR: 
%@ ERROR: Use the --stack_limit=size[KMG] command line option or
%@ ERROR: ?- set_prolog_flag(stack_limit, 2_147_483_648). to double the limit.
%@    Call: (16) _192=true ? creep
%@    Exit: (16) true=true ? abort

pt_term_factorized_test :-
	random_term(R),
	pt_term_coa([R], [I], Coa),
	pt_term_factorized(I, Coa, _, _).

% ?-N=100000, time(forall(between(1, N, _), pt_term_factorized_test)).
%@ % 32,861,119 inferences, 1.482 CPU in 1.655 seconds (90% CPU, 22170772 Lips)
%@ N = 100000.

jan · September 8, 2026, 3:29pm

Can you share the code? This looks like a good improvement.

kuniaki.mukai · September 8, 2026, 4:19pm

where and how to send this gz file to you ?
pt_c_sources.tar.gz

For just in case, here is load_c to load the SOs.

:- setenv('C', 'devel/zdd/C').   % ~/devel/zdd/C
%
:- use_module(library(process)).

user:load_c(BaseName) :-
	getenv('C', CLoc),
	atomic_list_concat(["~", CLoc, BaseName], "/", FBN),
	expand_file_name(FBN, [FullBaseName]),
    atomic_list_concat([FullBaseName, '.c'], CFile),
    atomic_list_concat([FullBaseName, '.so'], SoFile),
    (   ( \+ exists_file(SoFile) ; time_file(CFile, TC), time_file(SoFile, TS), TC > TS )
    ->  format('% 🔄 ビルド開始: ~w...~n', [CFile]),
        catch(
            % Added -Wl,-exported_symbol,_install to force symbol visibility
            ( process_create(path('swipl-ld'),
%                             ['-shared', '-Wl,-exported_symbol,_install', '-o', FullBaseName, CFile],
							 ['-shared', '-o', FullBaseName, CFile],
                             [process(PID)]),
              process_wait(PID, Status) ),
            Err,
            ( print_message(error, Err), Status = error )
        ),
        (   Status == exit(0)
        ->  format('% ✅ ビルド成功: ~w~n', [SoFile])
        ;   format('% ❌ ビルド失敗 (Status: ~w)~n', [Status]), fail
        )
    ;   true
    ),
    user:use_foreign_library(FullBaseName).

kuniaki.mukai · September 8, 2026, 4:55pm

I have sent by private email service of this forum.
Drop extension .c of attached for expansion.

jan · September 8, 2026, 5:16pm

Thanks! Good to see some nice Unicode and Japanese at work as well . I’ll share how it works out.

jan · September 8, 2026, 5:41pm

I have a bit of trouble connecting the dots. The two C files implement extract_coa/3 and paige_tarjan/3. What is missing?

kuniaki.mukai · September 8, 2026, 6:00pm

Time is money. I prefer to send it as it is than taking time to polish for readality.
Sorry for many possible junks for you.

(pt-min-coa.pl)

:- module(pt_mincoa, [pt_term_min_coa/3, pt_term_coa/3]).

:- use_module(library(random_terms)).
:- load_c('extract_coa').
:- load_c('paige_tarjan').

% ?- test_paige_tarjan_1.

% Test 1: Acyclic term graph with duplicate subterms (nodes 2 and 3 are identical atoms)
test_paige_tarjan_1 :-
    X = nodes(f(2, 3), a, a),
    paige_tarjan(X, M, Q),
    format('M: ~w~nQ: ~w~n', [M, Q]).

% ?- test_paige_tarjan_2.
% Test 2: Symmetric cyclic rational term graph (mutual recursion)
test_paige_tarjan_2 :-
    X = nodes(f(2), f(1)),
    paige_tarjan(X, M, Q),
    format('M: ~w~nQ: ~w~n', [M, Q]).

% ?- test_paige_tarjan_3.
% Test 3: Self-referential single node cycle (rational term)
test_paige_tarjan_3 :-
    X = nodes(f(1)),
    paige_tarjan(X, M, Q),
    format('M: ~w~nQ: ~w~n', [M, Q]).

% ?- test_paige_tarjan_4.
% Test 4: Mixed types with distinct equivalence classes
test_paige_tarjan_4 :-
    X = nodes(g(2, 3), f(4), a, a),
    paige_tarjan(X, M, Q),
    format('M: ~w~nQ: ~w~n', [M, Q]).


%! random_pt_variant(R, T).
%  create random term R, mimimize to T, check it that R and T are variant.
% ?- N = 10000, forall(between(1, N, _), random_pt_variant(R, T)).
% ?- N = 100000, time(forall(between(1, N, _), random_pt_variant(R, T))).

% ?- N = 1000, forall(between(1, N, _), random_pt_variant(R, T)).
% ?- N = 10000, forall(between(1, N, _), random_pt_variant(R, T)).
% ?- N = 100000, forall(between(1, N, _), random_pt_variant(R, T)).

% ?- generate_hydras(10, _L), time(user:extract_coa(_L, _Js, _Tvec)),
%	time(user:paige_tarjan(_Tvec, Coa, _Qvec)),
%	functor(Coa,_,Csize), functor(_Qvec, _, Qsize).
%@ % 3,079 inferences, 0.002 CPU in 0.003 seconds (69% CPU, 1542585 Lips)
%@ % -1 inferences, 0.000 CPU in 0.000 seconds (78% CPU, -58824 Lips)
%@ Coa = u(h(1, 1)),
%@ Csize = 1,
%@ Qsize = 55 .

% ?- generate_hydras(100, _L), time(user:extract_coa(_L, _Js, _Tvec)),
%	time(user:paige_tarjan(_Tvec, Coa, _Qvec)),
%	functor(Coa,_,Csize), functor(_Qvec, _, Qsize).
%@ % 25,507,549 inferences, 17.221 CPU in 24.645 seconds (70% CPU, 1481191 Lips)
%@ % -1 inferences, 0.001 CPU in 0.001 seconds (90% CPU, -1218 Lips)
%@ Coa = u(h(1, 1)),
%@ Csize = 1,
%@ Qsize = 5050 .

% ?- generate_hydras(100, _L), time(user:extract_coa(_L, _Js, _Tvec)).
% ?- generate_hydras(100, _L).

%
random_pt_variant(R, T):-
	random_term(R),
	pt_term_coa([R], [I], Coa),
	coa_to_equations(Coa, V, E),
	maplist(call, E),
	arg(I, V, T),
	variant(R, T).

%! random_pt_variant_swi(R, T)
% ?- N = 10, forall(between(1, N, _), random_pt_variant_swi(R, T)).
% ?- N = 100, forall(between(1, N, _), random_pt_variant_swi(R, T)).
% ?- N = 1000, time((forall(between(1, N, _), random_pt_variant_swi(R, T)))).
random_pt_variant_swi(R, T):-
	random_term(R),
	pt_term_coa_swi([R], [I], Coa),
	coa_to_equations(Coa, V, E),
	maplist(call, E),
	arg(I, V, T),
	variant(R, T).

%
coa_to_equations(Coa, U, E):- functor(Coa, F, N),
	functor(U, F, N),
	coa_to_equations(1, Coa, U, E).

%
coa_to_equations(I, V, U, [X = Tc|E]):- arg(I, V, T), !,
	arg(I, U, X),
	convert_index_in_args(V, U, T, Tc),
	J is I + 1,
	coa_to_equations(J, V, U, E).
coa_to_equations(_, _, _, []).

%
convert_index_in_args(_, _, T, T):-
	(	var(T); atomic(T); compound_name_arity(T, _, 0) ), !.
convert_index_in_args(V, U, T, Tc):-  functor(T, F, K),
	functor(Tc, F, K),
	convert_index_in_args(1, V, U, T, Tc).
%
convert_index_in_args(I, V, U, T, Tc):-arg(I, T, J), !,
	arg(J, U,  X),
	arg(I, Tc , X),
	K is I + 1,
	convert_index_in_args(K, V, U, T, Tc).
convert_index_in_args(_, _, _, _, _).

%  hydra as test data.

hydra(0, _).
hydra(N, h(X, X)):- N>0, N0 is N-1, hydra(N0, X).

% ?- hydra(3, H, A), A=H.
hydra(0, A, A).
hydra(N, h(X, X), A):- N>0, N0 is N-1, hydra(N0, X, A).
%

% ?- generate_hydras(3, X).
generate_hydras(0, []) :- !.
generate_hydras(N, [H|Rest]) :-
    N > 0,
    hydra(N, H, A), A=H,
    N1 is N - 1,
    generate_hydras(N1, Rest).

% ?- create_hydra_list(1, 100, _L), pt_term_coa(_L, _V, Coa, Q).
% ?- time((create_hydra_list(1, 1000, _L), pt_term_coa(_L, _V, Coa, Q))).
create_hydra_list(N, M, L):- numlist(N, M, Is),
	maplist(create_hydra, Is, L).
%
create_hydra(K, H):- hydra(K, H, H).

% ?- N = 10, Len = 10000,
%	time(forall(between(1, N, _), (find_shared_nodes(Len, Cs), writeln(Cs)))).
find_shared_nodes(Len, Cs):-compound_random_terms(Len, Rs),
			   pt_term_coa(Rs, _, Mvec, Qmap),
			   non_singleton_equiv_classes(Qmap, Mvec, Cs).

% ?-forall(between(1, 100000, _), compound_random_term(_)).
compound_random_term(R) :- random_term(R, []),
	compound(R),
	compound_name_arity(R, _, K),
	K > 0,
	!.
compound_random_term(R) :- compound_random_term(R).
% ?-  compound_random_terms(1000, _Rs).
compound_random_terms(N, Rs):- length(Rs, N),
	  maplist(compound_random_term,  Rs).
%
random_term(R):- random_term(R, []).
%
random_terms(N, Rs):- length(Rs, N),  maplist(random_term,  Rs).
%
rel_to_fun(L, R):- sort(L, L0), rel_to_fun(L0, [], R).
%
rel_to_fun([], X, X).
rel_to_fun([A-B|L], [A-U|V], R):-!,
	rel_to_fun(L, [A-[B|U]|V], R).
rel_to_fun([A-B|L], U, R):-!,
	rel_to_fun(L, [A-[B]|U], R).

% ?- non_singleton_equiv_classes(a(1,2,1,2), v(a, f(b)), R).
% ?- non_singleton_equiv_classes(a(1,2,1,2), v(A, f(b)), R).
non_singleton_equiv_classes(Qmap, Mvec, Cs):- Qmap=..[_|Js],
	functor(Qmap, _, N),
	zip(Js, Is, R),
	numlist(1, N, Is),
	sort(R, R0),
	rel_to_fun(R0, Fun),
	findall(X, (member(X-V, Fun),
				V=[_,_|_],
				arg(X, Mvec, A),
				compound(A)
			   ),
			Cs).

% ?- pt_term_coa([a,A,"abc", 1], Is, Tvec, Mvec).
% ?- pt_term_coa([a,a()], Is, Tvec, Mvec).
% ?- pt_term_coa([b(),a()], Is, Tvec, Mvec).
% ?- pt_term_coa([a,b], Is, Tvec, Mvec).
% ?- X=f(X, X), Y=f(Y,Y),  pt_term_coa([X, Y], Is, Tvec, Mvec).

% ?- N = 100000, hydra(N, _H, _R), _H = _R,
%	time(term_factorized(_H, _A, _L)),
%	time(pt_term_coa([_H], _I, _M, _Q)).
%@ % 96,662 inferences, 56.576 CPU in 56.641 seconds (100% CPU, 1709 Lips)
%@ % 4 inferences, 0.004 CPU in 0.005 seconds (75% CPU, 985 Lips)
%@ N = 100000 .

% ?- N = 100, hydra(N, _H, _R), _H = _R,
%	term_factorized(_H, _A, _L), length(_L, Len), writeln(shared=Len),
%	pt_term_coa([_H], _I, _M, Q).

% ?- time((generate_hydras(100, _Xs), extract_coa(_Xs, _Js, _Tvec))),
% time(user:paige_tarjan(_Js, _Tvec, Is)).
%@ % 25,517,950 inferences, 17.169 CPU in 24.511 seconds (70% CPU, 1486251 Lips)
%@ % -1 inferences, 0.000 CPU in 0.000 seconds (71% CPU, -100000 Lips)
%@ Is = q(1, 2) .


% ?- time((generate_hydras(2, _Xs), extract_coa(_Xs, _Js, _Tvec))),
% time(user:paige_tarjan(_Js, _Tvec, Is)).

% ?- generate_hydras(100, _Xs), time(pt_term_coa(_Xs, _Js, Coa)).
%@ % 25,507,753 inferences, 17.224 CPU in 24.660 seconds (70% CPU, 1480911 Lips)
%@ Coa = u(h(1, 1)) .

% ?- pt_term_coa([a], X, Y).
% ?- current_predicate(user:extract_coa/3).

pt_term_coa(Xs, Is, Mvec) :- pt_term_coa(Xs, Is, Mvec, _).
%
pt_term_coa(Xs, Is, Mvec, Qmap) :-
    user:extract_coa(Xs, Js, Tvec),
	user:paige_tarjan(Tvec, Mvec, Qmap),
%	print_size_diff(Tvec, Mvec),
    apply_quotient_map(Js, Qmap, Is).
%
print_size_diff(X, Y):-
	functor(X, _, A),
	functor(Y, _, B),
	C is A - B,
	(	C = 0 -> true
	;	writeln(C)
	).


% ?- time((generate_hydras(1000, _Xs), pt_term_coa_swi(_Xs, Is, _Mvec))).
% ?- pt_term_coa_swi([a,A,"abc", 1], Is, Tvec, Mvec).
% ?- pt_term_coa_swi([a,a()], Is, Tvec, Mvec).
% ?- pt_term_coa_swi([b(),a()], Is, Tvec, Mvec).
% ?- pt_term_coa_swi([a,b], Is, Tvec, Mvec).
% ?- X=f(X, X), Y=f(Y,Y),  pt_term_coa_swi([X, Y], Is, Tvec, Mvec).
%
pt_term_coa_swi(Xs, Is, Mvec) :- pt_term_coa_swi(Xs, Is, Mvec, _).
%
pt_term_coa_swi(Xs, Is, Mvec, Qmap) :-
    user:extract_coa(Xs, Js, Tvec),
	paige_tarjan_p(Tvec, Mvec, Qmap),
    apply_quotient_map(Js, Qmap, Is).

%
zip([], [], []).
zip([X|Xs], [Y|Ys], [X-Y|XYs]):- zip(Xs, Ys, XYs).
%
same_member(X, [Y|_]):- same_term(X, Y),!.
same_member(X, [_|R]):- same_member(X, R).

% ?- indexed_member(I-a, [J-a]).
% ?- indexed_member(I-a, [J-b]).
indexed_member(I-X, [J-Y|_]):- same_term(X, Y), !, I = J.
indexed_member(A, [_|R]):- indexed_member(A, R).

%%--------------------------------------------------
%!	pt_minimum_coa(+Xs, -Is, -Mvec) id det.
%
% ?- X=f(X), pt_minimum_coa([X, f(X), f(f(X))], Is, Mvec).
% ?- X=f(Y,0),Y=f(X,1), pt_minimum_coa([X, Y, f(X,Y), f(Y,X)], Is, Mvec).
% ?- X=f(Y,0),Y=f(X,0), pt_minimum_coa([X, Y, f(X,Y), f(Y,X)], Is, Mvec).
% ?- pt_minimum_coa([X, Y, 1, "ab", ab], Is, Mvec).

pt_minimum_coa(Xs, Is, Mvec):-
	pt_minimum_coa_by_same_term(Xs, Is, Mvec).
%
pt_minimum_coa_by_same_term(Xs, Is, Mvec):-
	prepare_coa_as_vector(Xs, Js, Tvec),
	paige_tarjan_p(Tvec, Mvec, Qmap),
	apply_quotient_map(Js, Qmap, Is).

% ?- prepare_coa_as_vector([a], Is, Tvec).
% ?- X=f(Y,0), Y=f(X, 1), prepare_coa_as_vector([X, Y], Is, Tvec).
prepare_coa_as_vector(Xs, Is, Tvec):-
	collect_subterms(Xs, [], Ts),
	zip(_, Ts, Dict),
	term_index(Xs, Is, Dict),
	dict_to_coa(Dict, Dict, Coa),
	coa_vector(Coa, Tvec).

% ?- X = f(X), collect_subterms([X], [], Ts), zip(_, Ts, Dict),
% dict_to_coa(Dict, Dict,  R).
% ?- X = f(X,Y), Y=f(X, f(Y)), collect_subterms([X,Y], [], Ts), zip(_, Ts, Dict),
% dict_to_coa(Dict, Dict, R).
% ?- X = f(X,Y), Y=f(X, f(Y)), collect_subterms([X,Y], [], Ts), zip(_, Ts, Dict),
% dict_to_coa(Dict, Dict, R), index_coa(R).
% ?- X = f(X,Y), Y=f(X, f(Y)), collect_subterms([X,Y], [], Ts), zip(_, Ts, Dict),
% dict_to_coa(Dict, Dict, R), coa_vector(R, D).

dict_to_coa([], _, []).
dict_to_coa([I-T|D0], Dict, [I-T0|Coa]):-
	(	\+ compound(T) -> T0 = T
	;	T=..[F|As],
		term_index(As, Is,  Dict),
		T0 =..[F|Is]
	),
	dict_to_coa(D0, Dict, Coa).

%
term_index([], [], _).
term_index([A|As], [J|Is], Dict):-
	indexed_member(J-A, Dict),
	!,
	term_index(As, Is, Dict).

% ?- apply_quotient_map([4, 3], q(2, 1, 4, 1, 2), R).
apply_quotient_map([], _, []).
apply_quotient_map([I|Is], Vec, [J|Js]):- arg(I, Vec, J),
	apply_quotient_map(Is, Vec, Js).

%!	collect_subterms(+L, +X, -Y) is det.
%	L, X, Y, are  a list of terms.
%	y in Y iff y in X or y is a subterm of some element of L.
%	identity is decided by same_term/2, not ==/2.
%
% ?- collect_subterms([a,B,B, C], [], X).
% ?- A=f(A), collect_subterms([A], [], X).
% ?- B=f(C,D), C=f(B, 0), D=g(B, C), collect_subterms([a,B,C], [], X).
% ?- collect_subterms([a, b], [], X).
% ?- collect_subterms([a(1), b], [], X).
% ?- X=f(X), collect_subterms([X], [], R).

collect_subterms(L, X, Y):- collect_same_subterms(L, X, Y).

%
collect_same_subterms([], X, X):-!.
collect_same_subterms([A|As], X, Y):-
	(	same_member(A, X) -> X1 = X
	;	\+ compound(A) -> X1 = [A|X]
	;	X0 = [A|X],
		A =..[_|Bs],
		collect_same_subterms(Bs, X0, X1)
	),
	collect_same_subterms(As, X1, Y).

%
coa_vector(Coa, Vec):- zip(Is, Ls, Coa),
	length(Is, N),
	numlist(1, N, Is),
	Vec=..[v|Ls].

% ?- print_vector(g(a(b), c(d))).
print_vector(Vec):- functor(Vec,_,N),
	forall(between(1, N, J), (arg(J, Vec, A), writeln(J=A))).

%%---------------------------------------------------------
%!  paige_tarjan_p(+X, -M, -Q) is det.
%
%   paige_tarjan.c と入出力（クラスID順、代表元順）が完全に一致する同値版。

% ?- paige_tarjan_p(f(1), M, Q).
%@ M = u(1),
%@ Q = q(1).

%	?- paige_tarjan_p(#(f(1,2), f(2,1)), M, Q).
%@ M = u(f(1, 1)),
%@ Q = q(1, 1).

% ?- X = #(f(), f, g(), g, A, "bbb"),
%	 paige_tarjan_p(X, M, Q).
%@ X = #(..),
%@ M = u(f(), f, g(), g, A, "bbb"),
%@ Q = q(1, 2, 3, 4, 5, 6).

paige_tarjan_p(X, M, Q) :-
    functor(X, _, N),
    (   N =:= 0
    ->  M = u, Q = q
    ;   functor(BlockMap, q, N),
        initial_partition(X, N, InitialBlocks),
        initialize_block_map(InitialBlocks, BlockMap),
        refine_partition(InitialBlocks, X, BlockMap, FinalBlocks),
        normalize_blocks_appearance(FinalBlocks, OrderedBlocks),
        build_outputs(OrderedBlocks, X, N, M, Q)
    ).

% --- Phase 1: 初期分割（インデックス直接走査） ---
initial_partition(X, N, Partition) :-
    collect_terms(1, N, X, Items),
    partition_initial(Items, Partition).

collect_terms(I, N, _, []) :- I > N, !.
collect_terms(I, N, X, [I-Term|Rest]) :-
    arg(I, X, Term),
    I1 is I + 1,
    collect_terms(I1, N, X, Rest).

partition_initial([], []).
partition_initial([I-Term|Rest], [[I|Block]|Blocks]) :-
    extract_initial_equiv(Rest, Term, Block, Remainder),
    partition_initial(Remainder, Blocks).

extract_initial_equiv([], _, [], []).
extract_initial_equiv([J-T2|Rest], T1, [J|Matched], Unmatched) :-
    initial_equiv(T1, T2), !,
    extract_initial_equiv(Rest, T1, Matched, Unmatched).
extract_initial_equiv([Item|Rest], T1, Matched, [Item|Unmatched]) :-
    extract_initial_equiv(Rest, T1, Matched, Unmatched).

initial_equiv(A, B) :-
    (   var(A)
    ->  A == B
    ;   var(B)
    ->  fail
    ;   compound(A)
    ->  compound(B),
        compound_name_arity(A, F1, N1),
        compound_name_arity(B, F2, N2),
        F1 == F2, N1 =:= N2
    ;   \+ compound(B),
        A == B
    ).

% --- Phase 2: 反復細分化 ---
initialize_block_map([], _).
initialize_block_map([[Rep|Rest]|Blocks], BlockMap) :-
    arg(Rep, BlockMap, Rep),
    maplist(init_map_to(Rep, BlockMap), Rest),
    initialize_block_map(Blocks, BlockMap).

init_map_to(Rep, BlockMap, I) :-
    arg(I, BlockMap, Rep).

refine_partition(Blocks, X, BlockMap, FinalBlocks) :-
    split_blocks(Blocks, X, BlockMap, NewBlocks, false, Changed),
    (   Changed == true
    ->  refine_partition(NewBlocks, X, BlockMap, FinalBlocks)
    ;   FinalBlocks = Blocks
    ).

split_blocks([], _, _, [], Changed, Changed).
split_blocks([Block|Blocks], X, BlockMap, NewBlocks, FlagIn, FlagOut) :-
    split_single_block(Block, X, BlockMap, Split),
    (   Split = [_]
    ->  NewBlocks = [Block|RestNew],
        split_blocks(Blocks, X, BlockMap, RestNew, FlagIn, FlagOut)
    ;   update_split_map(Split, BlockMap),
        append(Split, RestNew, NewBlocks),
        split_blocks(Blocks, X, BlockMap, RestNew, true, FlagOut)
    ).

split_single_block(Block, X, BlockMap, Split) :-
    partition_by_transitions(Block, X, BlockMap, Split).

partition_by_transitions([], _, _, []).
partition_by_transitions([I|Rest], X, BlockMap, [[I|SubBlock]|SubBlocks]) :-
    extract_same_transition(Rest, I, X, BlockMap, SubBlock, Remainder),
    partition_by_transitions(Remainder, X, BlockMap, SubBlocks).

extract_same_transition([], _, _, _, [], []).
extract_same_transition([J|Rest], I, X, BlockMap, [J|Matched], Unmatched) :-
    same_transition(I, J, X, BlockMap), !,
    extract_same_transition(Rest, I, X, BlockMap, Matched, Unmatched).
extract_same_transition([J|Rest], I, X, BlockMap, Matched, [J|Unmatched]) :-
    extract_same_transition(Rest, I, X, BlockMap, Matched, Unmatched).

same_transition(I, J, X, BlockMap) :-
    arg(I, X, TI),
    arg(J, X, TJ),
    (   compound(TI)
    ->  compound_name_arity(TI, _, A),
        (   A =:= 0
        ->  true
        ;   check_args(1, A, TI, TJ, BlockMap)
        )
    ;   true
    ).

check_args(K, A, TI, TJ, BlockMap) :-
    (   K > A
    ->  true
    ;   arg(K, TI, CI),
        arg(K, TJ, CJ),
        arg(CI, BlockMap, TargetI),
        arg(CJ, BlockMap, TargetJ),
        TargetI == TargetJ,
        K1 is K + 1,
        check_args(K1, A, TI, TJ, BlockMap)
    ).

update_split_map([_|NewBlocks], BlockMap) :-
    maplist(destructively_update_block(BlockMap), NewBlocks).

destructively_update_block(BlockMap, [Rep|Rest]) :-
    setarg(Rep, BlockMap, Rep),
    maplist(set_map_destructively(BlockMap, Rep), Rest).

set_map_destructively(BlockMap, Rep, I) :-
    setarg(I, BlockMap, Rep).

% --- 出現順正規化（最小インデックス昇順） ---
normalize_blocks_appearance(Blocks, OrderedBlocks) :-
    maplist(tag_min_index, Blocks, Keyed),
    keysort(Keyed, Sorted),
    maplist(get_block, Sorted, OrderedBlocks).

tag_min_index(Block, Min-Block) :-
    min_list(Block, Min).

get_block(_-Block, Block).

% --- Phase 3: 解の構成 ---
build_outputs(FinalBlocks, X, N, M, Q) :-
    length(FinalBlocks, M_Size),
    functor(Q, q, N),
    functor(M, u, M_Size),
    populate_q(FinalBlocks, 1, Q),
    populate_m(FinalBlocks, 1, X, Q, M).

populate_q([], _, _).
populate_q([Block|Blocks], ClassId, Q) :-
    maplist(set_q(Q, ClassId), Block),
    NextClass is ClassId + 1,
    populate_q(Blocks, NextClass, Q).

set_q(Q, ClassId, I) :-
    arg(I, Q, ClassId).

populate_m([], _, _, _, _).
populate_m([Block|Blocks], ClassId, X, Q, M) :-
    min_list(Block, Rep),
    arg(Rep, X, Term),
    (   compound(Term)
    ->  compound_name_arity(Term, F, A),
        (   A > 0
        ->  Term =.. [F|Args],
            maplist(get_map(Q), Args, MappedArgs),
            MTerm =.. [F|MappedArgs]
        ;   MTerm = Term
        )
    ;   MTerm = Term
    ),
    arg(ClassId, M, MTerm),
    NextClass is ClassId + 1,
    populate_m(Blocks, NextClass, X, Q, M).

get_map(BlockMap, I, TargetId) :-
    arg(I, BlockMap, TargetId).

%-----------------------------------------------------
%!  find_back_edge_sources(+Root, +V, -Sources) is det.

% ?- V = u(f(1)),
%	find_back_edge_sources(1, V, L).
%@ V = u(f(1)),
%@ L = [1].

% ?- V = u(f(2), f(3), f(1), a),
% find_back_edge_sources(1, V, L).
%@ V = u(f(2), f(3), f(1), a),
%@ L = [3].

% ?- V = u(f(2), f(2), a),
%	find_back_edge_sources(1, V, L).
%@ V = u(f(2), f(2), a),
%@ L = [2].

% ?- V = u(h(1, 1)),
%	find_back_edge_sources(1, V, L).
%@ V = u(h(1, 1)),
%@ L = [1].

% ?- V = u(f(2), g(3, 4), a, f(1)),
% find_back_edge_sources(1, V, L).
%@ V = u(f(2), g(3, 4), a, f(1)),
%@ L = [4].

find_back_edge_sources(Root, V, Sources) :-
    functor(V, _, N),
    functor(ColorTable, c, N),   % var: 白, gray: 灰, black: 黒
    dfs_sources(Root, V, ColorTable, RawSources, []),
    sort(RawSources, Sources).

dfs_sources(U, V, ColorTable, Sources, Tail) :-
    arg(U, ColorTable, Color),
    (   Color == black
    ->  Sources = Tail
    ;   setarg(U, ColorTable, gray),
        arg(U, V, Term),
        (   compound(Term)
        ->  compound_name_arity(Term, _, Arity),
            traverse_edges(1, Arity, Term, U, V, ColorTable, Sources, Mid),
            setarg(U, ColorTable, black),
            Mid = Tail
        ;   setarg(U, ColorTable, black),
            Sources = Tail
        )
    ).

traverse_edges(K, Arity, _, _, _, _, Sources, Tail) :-
    K > Arity, !,
    Sources = Tail.
traverse_edges(K, Arity, Term, U, V, ColorTable, Sources, Tail) :-
    arg(K, Term, Next),
    (   integer(Next)
    ->  arg(Next, ColorTable, Color),
        (   Color == gray
        ->  % Back Edge 発見: 始点 U を収集
            Sources = [U|Mid1],
            traverse_edges_next(K, Arity, Term, U, V, ColorTable, Mid1, Tail)
        ;   var(Color)
        ->  dfs_sources(Next, V, ColorTable, Sources, Mid1),
            traverse_edges_next(K, Arity, Term, U, V, ColorTable, Mid1, Tail)
        ;   % black
            traverse_edges_next(K, Arity, Term, U, V, ColorTable, Sources, Tail)
        )
    ;   traverse_edges_next(K, Arity, Term, U, V, ColorTable, Sources, Tail)
    ).

traverse_edges_next(K, Arity, Term, U, V, ColorTable, Sources, Tail) :-
    K1 is K + 1,
    traverse_edges(K1, Arity, Term, U, V, ColorTable, Sources, Tail).


%%---------------------------------------------------------
%!  pt_term_factorize(+Root, +Coa, -Skeleton, -Equations) is det.
%
%   Root:      ルートのインデックス
%   Coa:       Paige-Tarjan の出力ベクトル u(...)
%   Skeleton:  循環が変数で切断された骨格項
%   Equations: [Var = SubTerm, ...] のバインド方程式リスト

% ?- pt_term_factorize(1, coa(f(1)), S, E).
%@ E = [S=f(S)].
% ?- pt_term_factorize(1, coa(a), S, E).


% ?-N=10, time(forall(between(1, N, _), term_factorize_test)).
%@ % 1,967 inferences, 0.000 CPU in 0.000 seconds (97% CPU, 9322275 Lips)
%@ N = 10.
% ?-N=100, time(forall(between(1, N, _), term_factorize_test)).


% ?- term_factorized(p(),_,_).
% 共有なし -> 成功（素通り）
% ?- term_factorized(f(q()), S, E).

% 共有あり（同一の q() が2箇所に出現） -> クラッシュ
% ?- term_factorized(f(q(), q()), S, E).

% ?- term_factorized(f(q()), S, E).
%@ S = f(q()),
%@ E = [].
% ?- term_factorized(f(q(), q()), S, E).
%@ S = f(_A, _A),
%@ E = [_A=q()].


% ?-N=1, time(forall(between(1, N, _), term_factorize_test)).
%@ % 87 inferences, 0.000 CPU in 0.000 seconds (83% CPU, 2485714 Lips)
%@ N = 1.
% ?-N=2, time(forall(between(1, N, _), term_factorize_test)).
%@ % 869 inferences, 0.000 CPU in 0.000 seconds (94% CPU, 6437037 Lips)
%@ N = 2.
% ?-N=10, time(forall(between(1, N, _), term_factorize_test)).
%@ % 2,451 inferences, 0.000 CPU in 0.000 seconds (97% CPU, 9111524 Lips)
%@ N = 10.
% ?-N=100, time(forall(between(1, N, _), term_factorize_test)).

term_factorize_test :-
	random_term(R),
	term_factorized(R, _, _).


% ?-N=100000, time(forall(between(1, N, _), pt_term_factorize_test)).
%@ % 32,908,725 inferences, 4.378 CPU in 5.050 seconds (87% CPU, 7516297 Lips)
%@ N = 100000.

pt_term_factorize_test :-
	random_term(R),
	pt_term_coa([R], [I], Coa),
	pt_term_factorize(I, Coa, _, _).

pt_term_factorize(Root, Coa, Skeleton, Equations) :-
    find_back_edge_sources(Root, Coa, BackSources),
    functor(Coa, _, N),
    functor(Vars, v, N),
    build_factorized(1, N, Coa, Vars, BackSources, Equations),
    arg(Root, Vars, Skeleton).

build_factorized(I, N, _, _, _, []) :-
    I > N, !.
build_factorized(I, N, Coa, Vars, BackSources, Equations) :-
    arg(I, Coa, Term),
    arg(I, Vars, VarI),
    (   memberchk(I, BackSources)
    ->  % Back Edge 始点: VarI は自由変数のまま残し、方程式 VarI = SubTerm を外に出す
        instantiate_subterm(Term, Vars, SubTerm),
        Equations = [VarI = SubTerm | RestEqs]
    ;   % 通常ノード: VarI に展開後の部分項を直接ユニフィケーション（インライン結合）
        instantiate_subterm(Term, Vars, VarI),
        Equations = RestEqs
    ),
    I1 is I + 1,
    build_factorized(I1, N, Coa, Vars, BackSources, RestEqs).

instantiate_subterm(Term, Vars, Instantiated) :-
    (   compound(Term)
    ->  compound_name_arity(Term, F, Arity),
        functor(Instantiated, F, Arity),
        fill_args(1, Arity, Term, Vars, Instantiated)
    ;   % 定数・アトム
        Instantiated = Term
    ).

fill_args(K, Arity, _, _, _) :-
    K > Arity, !.
fill_args(K, Arity, Term, Vars, Instantiated) :-
    arg(K, Term, Next),
    (   integer(Next)
    ->  arg(Next, Vars, TargetVar),
        arg(K, Instantiated, TargetVar)
    ;   arg(K, Term, Leaf),
        arg(K, Instantiated, Leaf)
    ),
    K1 is K + 1,
    fill_args(K1, Arity, Term, Vars, Instantiated).

jan · September 9, 2026, 6:00pm

Thanks. I had Claude working a lot on this. The current GIT version implements the following new stuff based on your prototype:

term_minimal(+Term, -Minimal)
Creates a minimal version of Term, i.e., a term that is ==/2 identical, but using maximal sharing of sub terms.
term_factorized/3 now uses your stuff. That makes it cycle safe and, on random terms, about 6 times faster than the old implementation.
variant_sha1/2 and variant_hash/2 now handle cyclic terms, producing the same hash or X = [a|X] and X = [a,a|X]. This makes “nb sets” cycle-safe.

These are the public predicates. src/pl-bisim.c also defines a couple of $-predicates. I’m open to make more API accessible if there is a good reason.

Claude suggests this could also fix the currently inconsistent ordering of cyclic terms. I have not studied this. If this provides an answer it might be worth implementing as this allows sort/2, binary trees, etc., to reliably handle arbitrary cycling terms.

kuniaki.mukai · September 10, 2026, 12:35am

I am glad to hear that. Gemini (I use) thinks [] as atom, and did not know that ‘[X|Y]’ is equalt to [|](X, Y)', etc. So extract_coa/3` has stupid bugs, though now (yesterday) it was fixed.

Now I need to get physical address of terms in order to get O(N) time efficiency to build
the COA (flat term equations sytem) for paige-tarjan famous O(m log n) mehod for term graphs.
Currently gemini found physical address access is not allowed in API. So I have to dig into the git source for to get in C-source.

Claude is known to be wise and you are best human to dictate it to create SWI-Prolog
libraries. On the other hand, my current interest is very very narrow.
It is a small library around Page-Tarjan in a way as far as I can ask Gemini to write C-source for that.

I am satisfied with the paige_tarjan.c in efficency, which I sent you, and
now consider to ask Gemini to create Paige-Tarjan on general directed graphs. Gemini seems
familiar with Paige-Tarjan. It will surely create C-source in a blink after I write a prompt.
Although my concern is narrow, I feel AI is useful to add easily C-leve efficiency to Prolog programming.

jan · September 10, 2026, 9:32am

Claude did some serious rewrites. Using physical addresses in the internal low-level API is possible, so it uses those. The public API does not expose these as they are not constant: GC and stack shifts can change the addresses due to pretty much any API call. In the low level interface we can do so, doing the allocation with specialized routines that return a stack overflow rather than enlarging the stacks. If that happens we revert everything to the start state, enlarge the stacks and try again. Pretty much all internal term manipulation works that way.

It also discovered the paige_tarjan was essentially O(N^2) and modified it to become nearly O(N).

The current implementation is a lot faster than yours as a result. As is, the low level APIs are not exposed. I’m happy to consider exposing these as well. Just let me know what you want.

kuniaki.mukai · September 10, 2026, 10:20am

I need it for only to build “pre colagegra” vector, to
pass it “paige-tarjan” metho, for instance. The physical address is used as unique index for a node in the input prolog term, that’s all.
GC and outer process are irrelevant, as far as they do nothing while building the vector, though slightly not sure for now, on free variables in the input terms.

O(N) sounds big inovation. In general it is known as O(M log_2 N)
where M is the number of edges, N the number of nodes. Congratulation !

kuniaki.mukai · September 10, 2026, 10:43am

Hi Jan,

Regarding the low-level API, exposing a minimal interface to extract the unreduced Kripke structure (automaton) from rational terms would be ideal. Specifically, we need to traverse cyclic terms and register nodes based strictly on physical cell identity rather than semantic equality (which performs Bisimulation reduction prematurely).

Having access to either of the following would solve the problem cleanly:

A cell-identity test or hash primitive:
- An identity check function:
  
  int PL_is_identical_cell(term_t t1, term_t t2);
  
  (returns whether two terms dereference to the exact same heap cell / Word*)
- Or a physical-identity term hash:
  
  uintptr_t PL_physical_cell_id(term_t t);
  
  (returning an ID or dereferenced address valid within a non-shifting scope)
A scoped, non-shifting traversal block (if direct Word* access is exposed):
- A way to freeze stack reallocation / GC during the graph extraction phase, or access to the underlying deRefprimitive so we can inspect functor/args without triggering allocation:
  
  Word PL_deref(term_t t);

If your rewritten internal implementation of build_coa / Paige-Tarjan already extracts the raw transition vector and runs the near-O(N) minimization internally, having that exported as a high-level foreign predicate (or seeing how you structure the low-level loop) would also be fantastic.

Please let me know how you prefer to expose this. Thanks again for your support.

Best regards,

Kuniaki

jan · September 10, 2026, 10:51am

That is PL_same_term(). The problem is you can’t hash on that, so you’ll have to test every term to be “same”. The internal version replaces the functor with a tagged integer that denotes the location in the array. That way, we can do a simple linear scan through the term to build the initial array. No sensible public API can do so because it makes the stacks “invalid”. That is fine as long as we can guarantee nothing that requires valid stacks happens before we restore the stack.

This cannot meaningfully exist. It would be way to fragile. The API is hard enough as it is, considering the frequent issues people have with it.

This has more or less the same problem

I think that makes sense. I’ll look into it.

kuniaki.mukai · September 10, 2026, 11:11am

Hi Jan,

Thank you for the crystal-clear explanation. The in-place marking technique (temporarily overwriting the functor with a tagged array index and restoring it afterwards) makes complete sense. It is the classic, most efficient way to achieve a true O(N) traversal without hashing overhead, and it fully explains why exposing raw cell pointers to the public API is both fragile and impossible.

Having the entire pipeline—from the raw graph extraction via in-place traversal to the optimized Paige-Tarjan minimization—packaged and exported as an internal/high-level foreign predicate is definitely the cleanest and most robust approach.

I look forward to seeing it when you have time to look into it. Thank you again for your incredible insight and support.

Best regards,

Kuniaki

jan · September 10, 2026, 5:12pm

Give the current source version a try. See

?- apropos(automaton).

to see the predicates. If you do so in swipl-win, Ctrl-click will open the docs for the predicate.

Let is know what yiu think about the API and performance!

kuniaki.mukai · September 10, 2026, 5:36pm

Thanks for the update. I’ll test the new automaton predicates in the latest source right away.

In parallel, I did a quick benchmark extracting the raw graph from Prolog using ‘$term_id’/2 and library(rbtrees) (O(N log N)) followed by Paige-Tarjan minimization:

?- N = 1000,
   time((create_hydra_list(1, N, _Xs),
         pt_term_raw_coa(_Xs, _Js, _Coa),
         functor(_Coa, _, Size),
         length(_Js, JLen))).
% 209,404,134 inferences, 13.010 CPU in 13.583 seconds (96% CPU, 16095667 Lips)
% N = JLen, JLen = 1000,
% Size = 1 .

Separating physical cells with ‘$term_id’/2 prevented pre-collapse and correctly reduced the 1,000 Hydra instances to Size = 1 in ~13s without O(N^2) explosion.

I will now benchmark your C-level implementation against this.

Topic		Replies	Views
Fast Term Factorization via Paige-Tarjan: Bypassing ==/2 and compare/2 General	16	299	July 15, 2026
Term_factorized, term_minimal, term_automaton, automaton_minimal General	14	215	September 20, 2026
Cyclic terms unification `X=f(f(X), X), Y=f(Y, f(Y)), X = Y. ` General	70	1535	August 29, 2025
Instability of standard order Request For Comments	16	247	October 5, 2026
How to compare/3 without surprises on non-ground terms? Help!	81	2777	March 31, 2023

Term_factorized/3 error for compound_name_arity zero term

Related topics