1:- module(path_count5, [
    2	]).    3
    4:- use_module(library(apply)).    5:- use_module(library(apply_macros)).    6:- use_module(library(clpfd)).    7:- use_module(library(statistics)).    8:- use_module(zdd('zdd-array')).    9:- use_module(util(math)).   10:- use_module(util(meta2)).   11:- use_module(pac(basic)).   12:- use_module(pac(meta)).   13:- use_module(util(misc)).   14:- use_module(pac('expand-pac')). % For the kind block.
   15:- use_module(zdd('zdd-misc')).   16:- use_module(pac(op)).   17:- use_module(zdd(zdd)).   18
   19:- set_prolog_flag(stack_limit, 10_200_147_483_648).   20
   21 :- op(1060, xfy, ~).		% equivalence
   22 :- op(1060, xfy, #).		% exor
   23 :- op(1060, xfy, <->).		% equivalence
   24 :- op(1050, yfx, <-).   25 :- op(1060, xfy, <=> ).	% equivalence
   26 :- op(1040, xfy, \/).		% OR
   27 :- op(1030, xfy, /\).		% AND
   28 :- op(1020, fy, \).		% NOT
   29 :- op(700, xfx, :=).		% Assignment
   30 :- op(1000, xfy, &).   31
   32% for pac query.
   33 :- pac:op(1000, xfy, &).   34 :- pac:op(700, xfx, :=).   35
   36term_expansion --> pac:expand_pac.
   37
   38% some tiny.
   39x(p(X, _), X).
   40y(p(_, Y), Y).
   41link_symbol(_->_).
   42change_symbol(A-B, A->B).
   43%
   44touch(A-_, _-A):-!.
   45touch(_-A, A-_).
   46
   47% ?- connect_pairs(a-b, [b-c, d-a], W).
   48% ?- connect_pairs(a-b, [a-c, a-b], W).  % false
   49
   50connect_pairs(X-Y, [U-X, Y-V], U-V):-!.
   51connect_pairs(X-Y, [Y-V, U-X], U-V).
   52
   53% ?- both_ends_on_frontier(p(3, 2)-p(3,5), p(2, 1)-p(3, 1)).
   54% ?- both_ends_on_frontier(p(2, 2)-p(3,5), p(2, 1)-p(3, 1)).  %false
   55% ?- both_ends_on_frontier(p(2, 2)-p(2,1), p(2, 1)-p(3, 1)).  %false
   56both_ends_on_frontier(p(S, _)-p(S, _), p(J, _)-p(K, _)):- S is max(J, K).
   57
   58% ?- compare_hybrid(C, c-d, c->d).
   59% ?- compare_hybrid(C, c->d, c->d).
   60% ?- compare_hybrid(C, c->d, c-d).
   61% ?- compare_hybrid(C, c-d, c-d).
   62% ?- predsort(compare_hybrid, [a->b, b->a], R).
   63% ?- predsort(compare_hybrid, [a-b, b-a], R).
   64% ?- predsort(compare_hybrid, [b->a, a->b, a-b, a->b, c-b, a-b, b-a], R).
   65% ?- zdd((set_compare(compare_hybrid), zdd_compare(C, a-b, c-d))).
   66% ?- zdd((set_compare(compare_hybrid), zdd_compare(C, c-d, c->d))).
   67% ?- zdd((set_compare(compare_hybrid), zdd_compare(C, c->d, c-d))).
   68% ?- zdd((set_compare(compare_hybrid), zdd_insert_atoms([a-b, c-d], 1, X))).
   69% ?- zdd((set_compare(compare_hybrid), zdd_insert_atoms([a->b, c-d], 1, X), psa(X))).
   70
   71compare_hybrid(C, X, Y):- compare(C0, X, Y),
   72	(	C0 = (=)  -> C = (=)
   73	;	functor(X, Fx, _),
   74		functor(Y, Fy, _),
   75		(	Fx == Fy
   76		-> 	(	Fx = (-)
   77			->  ( C0 = (<) -> C = (>) ; C = (<) )
   78			;	C = C0
   79			)
   80		;	(	Fx = (-) -> C = (<) ; C = (>) )
   81		)
   82	).
   83
   84:- meta_predicate rect(:).   85
   86rect(X):- zdd((set_compare(compare_hybrid), X)).
   87
   88% ?- rect((X<<{[p(0,0)-p(2,1)], [p(1,1)-p(2,3), p(2,2)-p(2,2)]},
   89%	prune_frontier(X, 2, p(0,0), Y), psa(Y))).
   90% ?- rect((X<<{[p(0,0)-p(2,1)], [p(1,1)-p(1,1), p(2,2)-p(2,2)]},
   91%	prune_frontier(X, 2, p(0,0), Y), psa(Y))).
   92% ?- rect((X<<{[p(0,0)-p(2,1)], [p(1,1)-p(1,1), p(2,2)-p(2,3)]},
   93%	prune_frontier(X, 2, p(0,0), Y), psa(Y))).
   94
   95prune_frontier(F, I, P, F0):- @(prune_frontier(F, I, P, F0)).
   96%
   97prune_frontier(F, _, _, F, _):- F<2, !.
   98prune_frontier(F, I, P, F0, S):- cofact(F, t(A, L, R), S),
   99	prune_frontier(L, I, P, L0, S),
  100	(	A = (U-V)
  101	->	(	(	U = V
  102			;	x(U, I),
  103				x(V, I)
  104			;	U = P,
  105				x(V, I)
  106			)
  107		->	prune_frontier(R, I, P, R0, S)
  108		;	R0 = 0
  109		)
  110	;	R0 = R
  111	),
  112	cofact(F0, t(A, L0, R0), S).
  113
  114
  115		/**********************
  116		*     bridge block    *
  117		**********************/
  118
  119fold_bridge_block(Es, X, Y, Z):- @(fold_bridge_block(Es, X, Y, 0, Z)).
  120
  121fold_bridge_block(Es, X, Y, Z0, Z):- @(fold_bridge_block(Es, X, Y, Z0, Z)).
  122
  123%
  124fold_bridge_block(_, _, [], _, Z, Z, _):- !.
  125fold_bridge_block(E, As, X, Y, Z0, Z, S):- cofact(X, t(A, L, R), S),
  126	(  link_symbol(A) ->   Z = Z0
  127	;  bridge_block(E, As, L, Y, Z0, Z1, S),
  128	   (	touch(E, A)
  129	   ->	bridge_block(A, E, As, R, Y, Z1, Z, S)
  130	   ;	bridge_block(E, [A|As], R, Y, Z1, Z, S)
  131	   )
  132	).
  133
  134
  135% ?- zdd((
  136%	E = p(1,0)-p(0,0),
  137%	X << {[p(0,0)-p(0,0)]},
  138%	Y << {[p(1,0)-p(1,0), p(1,1)-p(1,1), p(3,1)->p(3,2)]},
  139%	bridge_block(E, X, Y, Z),
  140%	psa(Z))).
  141
  142% ?- zdd((
  143%	E = p(1,0)-p(0,0), E0=p(0,0)-p(1,0),
  144%	X << {[p(0,0)-p(0,0), p(0,1)-p(0,1)]},
  145%	Y << {[p(1,0)-p(1,0), p(1,1)-p(1,1)]},
  146%	bridge_block(E, X, Y, 0, Z0),
  147%	zdd_join(Z0, X, X0),
  148%	bridge_block(E0, X0, Y, Z0, Z1),
  149%	psa(Z0), psa(Z1))).
  150
  151bridge_block(E, X, Y, Z):-  % assuming block order X < Y
  152	@(bridge_block(E, [], X, Y, 0, Z)).
  153
  154bridge_block(E, X, Y, Z0, Z):- @(bridge_block(E, [], X, Y, Z0, Z)).
  155
  156%
  157bridge_block(_, _, X, _, Z, Z, _):- X<2, !.
  158bridge_block(E, As, X, Y, Z0, Z, S):- cofact(X, t(A, L, R), S),
  159	(  link_symbol(A) ->   Z = Z0
  160	;  bridge_block(E, As, L, Y, Z0, Z1, S),
  161	   (	touch(E, A)
  162	   ->	bridge_block(A, E, As, R, Y, Z1, Z, S)
  163	   ;	bridge_block(E, [A|As], R, Y, Z1, Z, S)
  164	   )
  165	).
  166%
  167bridge_block(_, _, _, _, Y, Z, Z, _):- Y<2, !.
  168bridge_block(A, E, As, X, Y, Z0, Z, S):- cofact(Y, t(B, L, R), S),
  169	(	link_symbol(B) ->  Z = Z0
  170	;	bridge_block(A, E, As, X, L, Z0, Z1, S),
  171		(	touch(E, B)
  172		->	connect_pairs(E, [A, B], C),
  173			zdd_merge(X, R, U, S),
  174			change_symbol(E, E0),
  175			zdd_insert_atoms([E0, C|As], U, V, S),
  176			zdd_join(V, Z1, Z, S)
  177		;	bridge_block(A, E, [B|As], X, R, Z1, Z, S)
  178		)
  179	).
  180
  181
  182easy_bridge_block(Es, X, Y, Z):-  @(easy_bridge_block(Es, X, Y, Z)).
  183
  184easy_bridge_block(Es, X, Y, Z, S):-  s_zdd_merge(X, Y, Z0, S),
  185	mqp_update_links(Es, Z0, Z, S).
  186
  187		/*************
  188		*     end    *
  189		*************/
  190
  191		/****************
  192		*     bridge    *
  193		****************/
  194
  195
  196% ?- zdd((count_path(p(0,0), p(0, 1), C))).
  197% ?- zdd((count_path(p(0,0), p(1, 0), C))).
  198% ?- zdd((count_path(p(0,0), p(1,1), C))).
  199% ?- zdd((count_path(p(0,0), p(2,2), C))).
  200% ?- zdd((count_path(p(0,0), p(3,3), C))).
  201% ?- zdd((count_path(p(1,1), p(2,1), C))).
  202% ?- zdd((count_path(p(1,2), p(3,2), C))).
  203% ?- zdd((count_path(p(1,2),p(10,2), C))).
  204% ?- zdd((count_path(p(1,2),p(12,2), C))).
  205% ?- zdd((count_path(p(1,2),p(13,2), C))).
  206% ?- zdd((count_path(p(1,2),p(14,2), C))).
  207
  208% ?- zdd((count_path(p(1,2),p(100,2), C))).
  209% ?- time(zdd((count_path_line(p(1,2),p(1000, 2), C)))).
  210%@ % 4,691,349 inferences, 0.982 CPU in 1.154 seconds
  211% ?- time(zdd((count_path_line(p(1,2),p(10000,2), C)))).
  212%@ % 61,745,569 inferences, 59.803 CPU in 60.290 seconds
  213
  214count_path(X, Y, C):- @(count_path(X, Y, C)).
  215%
  216count_path(X, Y, C, S):-
  217	qp_block(X, Y, R, S),
  218	memo_final(R, X, Y, S),
  219	qp_list(Q, [X-X], S),
  220	path_tree(Q, Tree, S),
  221	card(Tree, C, S).
  222%
  223memo_final([], _, _, _).
  224memo_final([A|R], X, Y, S):-
  225	qp_list(A, List, S),
  226	(	qp_final(List, X, Y)
  227	->  memo(qp_final(A)-true, S)
  228	;  true
  229	),
  230	memo_final(R, X, Y, S).
  231%
  232qp_block(X, Y, R):- @(qp_block(X, Y, R)).
  233%
  234qp_block(X, X, [Q], S):-!,
  235	zdd_insert_atoms([X-X], 1, Q, S),
  236	memo(qp_suc(Q)-List, S),
  237	( var(List) -> List=[]; true).
  238qp_block(p(I, J), p(I, K), R, S):-!,
  239	M is (J+K)//2,
  240	qp_block(p(I, J), p(I, M), R1, S),
  241	M0 is M+1,
  242	qp_block(p(I, M0), p(I, K), R2, S),
  243	qp_joint([p(I, M)-p(I, M0), p(I, M0)-p(I, M)], R1, R2, R3, S),
  244	prune_block_local(R3, p(I, J), p(I, K), R, S).
  245qp_block(p(I, J), p(K, L), R, S):-
  246	M is (I+K)//2,
  247	qp_block(p(I, J), p(M, L), R1, S),
  248	M0 is M+1,
  249	qp_block(p(M0,J), p(K, L), R2, S),
  250	findall(p(M, H)-p(M0, H), between(J, L, H), U),
  251	findall(p(M0, H0)-p(M, H0), between(J, L, H0), V),
  252	append(U, V, W),
  253	sort(W, Es),
  254	qp_joint(Es, R1, R2, R3, S),
  255	prune_block_local(R3, p(I, J), p(K, L), R, S).
  256%
  257prune_block_local(X, _, _, X, _).
  258
  259% ?- zdd(rect_bundle(p(0,0), p(0,0), Bs, Cs)), length(Bs, Len).
  260% ?- zdd(rect_bundle(p(0,0), p(1,1), Bs, Cs)), length(Bs, Len).
  261% ?- zdd(rect_bundle(p(0,0), p(1,0), Bs, Cs)), length(Bs, Len).
  262% ?- time(zdd(rect_bundle(p(0,0), p(10,10), Bs, Cs))),
  263%	length(Bs, Len), length(Cs, LenC).
  264
  265rect_bundle(LL, UR, Bridges, Column_qps):-
  266	bridge_bundle(LL, UR, Bridges),
  267	column_qp_bundle(LL, UR, Column_qps).
  268%
  269column_qp_bundle(X, Y, Z):- @(column_qp_bundle(X, Y, Z)).
  270%
  271column_qp_bundle(p(I, J), p(K, L), Z, S):-
  272	qp_along_line(J, L, List),
  273	convert_to_col(List, Q, I, S),
  274	U is K - I,
  275	column_qp_bundle(U, [Q], Z0, S),
  276	reverse(Z0, Z).
  277
  278%
  279column_qp_bundle(0, Z, Z, _):- !.
  280column_qp_bundle(K, [Q|Z], Z0, S):-
  281	map_h_shift(Q, Q0, S),
  282	K0 is K-1,
  283	column_qp_bundle(K0, [Q0, Q|Z], Z0, S).
  284
  285%
  286zdd_repeat(0, _, X, X, _):-!.
  287zdd_repeat(N, F, X, Y, S):- N>0,  call(F, X, U, S),
  288	N0 is N-1,
  289	zdd_repeat(N0, F, U, Y, S).
  290
  291
  292% ?- bridge_bundle(p(0,0), p(0,0), S).
  293% ?- bridge_bundle(p(0,0), p(1,0), S).
  294% ?- bridge_bundle(p(0,0), p(1,1), S).
  295% ?- bridge_bundle(p(0,0), p(3,1), S), maplist(writeln, S).
  296% ?- time((bridge_bundle(p(0,0), p(100,2), S), length(S, Len))).
  297
  298bridge_bundle(p(I, _), p(I, _), []):-!.
  299bridge_bundle(p(I, J), p(K, L), Bundle):-
  300	findall(Bridge, (	between(I, K, U),
  301						U0 is U-1,
  302						U0 >= I,
  303						findall( H,
  304								(  between(J, L, V),
  305								   (  H = p(U0, V)-p(U, V)
  306								   ;  H = p(U, V)-p(U0, V)
  307								   )
  308								),
  309								Bridge)
  310					),
  311			Bundle).
  312
  313% ?- bridge_shift([[p(1,1)-p(1,1)]], X).
  314bridge_shift(X, Y):-maplist(maplist(h_shift_link), X, Y).
  315
  316%
  317h_shift_link(p(X, Y)-p(Z, U), p(X0, Y)-p(Z0, U)):-
  318	X0 is X+1,
  319	Z0 is Z+1.
  320%
  321v_shift_link(p(X, Y)-p(Z, U), p(X, Y0)-p(Z, U0)):-
  322	Y0 is Y+1,
  323	U0 is U+1.
  324
  325% ?- zdd(convert_to_row([[1-2], [3-4]], Row, 0)).
  326convert_to_row(X, Y, Z):- @(convert_to_row(X, Y, Z)).
  327%
  328convert_to_row([], [], _, _).
  329convert_to_row([X|Y], [X0|Y0], H, S):-
  330	convert_to_row_(X, X1, H),
  331	zdd_insert_atoms(X1, 1, X0, S),
  332	convert_to_row(Y, Y0, H, S).
  333
  334% ?- convert_to_row_([1-2, 3-4], R, 5).
  335convert_to_row_([], [], _).
  336convert_to_row_([I-J|X], [p(I, H)-p(J, H)|X0], H):-
  337	convert_to_row_(X, X0, H).
  338
  339% ?- zdd(convert_to_col([[1-2], [3-4]], Col, 0)).
  340convert_to_col(X, Y, Z):- @(convert_to_col(X, Y, Z)).
  341%
  342convert_to_col([], [], _, _).
  343convert_to_col([X|Y], [X0|Y0], H, S):-
  344	convert_to_col_(X, X1, H),
  345	zdd_insert_atoms(X1, 1, X0, S),
  346	convert_to_col(Y, Y0, H, S).
  347
  348% ?- convert_to_col_([1-2, 3-4], R, 5).
  349convert_to_col_([], [], _).
  350convert_to_col_([I-J|X], [p(H, I)-p(H, J)|X0], H):-
  351	convert_to_col_(X, X0, H).
  352
  353% ?- zdd((convert_to_col([[1-2], [3-4]], Col, 0),
  354%	map_h_shift(Col, Col0))).
  355map_h_shift(X, Y):- @(map_h_shift(X, Y)).
  356%
  357map_h_shift([], [], _).
  358map_h_shift([X|Y], [X0|Y0], S):- h_shift(X, X0, S),
  359	map_h_shift(Y, Y0, S).
  360
  361% ?- zdd((convert_to_row([[1-2], [3-4]], Row, 0),
  362%	map_v_shift(Row, Row0))).
  363map_v_shift(X, Y):- @(map_v_shift(X, Y)).
  364%
  365map_v_shift([], [], _).
  366map_v_shift([X|Y], [X0|Y0], S):- v_shift(X, X0, S),
  367	map_v_shift(Y, Y0, S).
  368
  369% ?- zdd((qp_along_line(1, 2, X), convert_to_col(X, Y, 3),
  370%	maplist(h_shift, Y, Y0),
  371%	forall(member(M, Y0), (qp_list(M, L), writeln(L))))).
  372%
  373h_shift(X, Y):- @(h_shift(X, Y)).
  374%
  375h_shift(X, X, _):- X < 2, !.
  376h_shift(X, Y, S):- cofact(X, t(p(I, J)-p(I, K), L, R), S),
  377	I0 is I + 1,
  378	h_shift(L, L0, S),
  379	h_shift(R, R0, S),
  380	cofact(Y, t(p(I0, J)-p(I0, K), L0, R0), S).
  381
  382% ?- zdd((qp_along_line(1, 2, X), convert_to_row(X, Y, 3),
  383%	maplist(v_shift, Y, Y0),
  384%	forall(member(M, Y0), (qp_list(M, L), writeln(L))))).
  385
  386v_shift(X, Y):- @(v_shift(X, Y)).
  387%
  388v_shift(X, X, _):- X < 2, !.
  389v_shift(X, Y, S):- cofact(X, t(p(I, J)-p(K, J), L, R), S),
  390	v_shift(L, L0, S),
  391	v_shift(R, R0, S),
  392	J0 is J + 1,
  393	cofact(Y, t(p(I, J0)-p(K, J0), L0, R0), S).
  394
  395
  396% ?- spy(solve_rect).
  397% ?- zdd((solve_rect(p(0,0),p(1,0), Z))).
  398% ?- zdd((solve_rect(p(0,0),p(1,1), Z))).
  399% ?- zdd((solve_rect(p(0,0),p(2,2), Z))).
  400% ?- zdd((solve_rect(p(0,0),p(3,3), Z))).
  401% ?- call_with_time_limit(300, zdd((solve_rect(p(0,0),p(4,4), Z)))).
  402%@ n(3,41)
  403%@ n(2,1681)
  404%@ n(1,68921)
  405%@ n(0,2825761)
  406%@ ERROR: Unhandled exception: Time limit exceeded
  407
  408
  409% ?- zdd((solve_rect(p(0,0),p(4,4), Z))).
  410%@ n(3)
  411%@ n(2)
  412%@ n(1)
  413%@ Z = [4996364, 4996368, 4996372, 4996376, 4996380, 4996384, 4996388, 4996392, 4996396|...] .
  414
  415
  416% ?- zdd((solve_rect(p(0,0),p(15,15), Z))).
  417
  418
  419solve_rect(X, Y, Z):- @(solve_rect(X, Y, Z)).
  420%
  421solve_rect(p(I, J), p(I, L), Z, S):-!, initial_column(I, J, L, Z, S).
  422solve_rect(p(I, J), p(K, L), Z, S):-
  423	initial_col_bridge(I, J, L, Col, B, S),
  424	N is (K-I)-1,
  425	map_h_shift(Col, Col0, S),
  426	repeat_bridge(N, Col, B, Col0, Z, S).
  427
  428%
  429repeat_bridge(0, X, B, C, Y, S):-!, length(X, Len), writeln(n(0,Len)),
  430	qp_joint(B, X, C, Y, S).
  431repeat_bridge(N, X, B, C, Y, S):- length(X, Len), writeln(n(N, Len)),
  432	qp_joint(B, X, C, Y0, S),
  433	bridge_shift(B, B0),
  434	map_h_shift(C, C0, S),
  435	N0 is N-1,
  436	repeat_bridge(N0, Y0, B0, C0, Y, S).
  437
  438
  439% ?- zdd((initial_col_bridge(0, 0, 3, C, B))).
  440% ?- time(zdd((initial_col_bridge(0, 0, 15, C, B)))).
  441initial_col_bridge(I, L, H, C, B):-	@(initial_col_bridge(I, L, H, C, B)).
  442%
  443initial_col_bridge(I, L, H, C, B, S):- initial_column(I, L, H, C, S),
  444	initial_bridge(I, L, H, B).
  445
  446
  447% ?- zdd((initial_column(1, 2, 4, Col),
  448%   forall(member(M, Col), (qp_list(M, L), writeln(L))))).
  449initial_column(I,Low,Hi, Col):- @(initial_column(I, Low, Hi, Col)).
  450%
  451initial_column(I,Low,Hi, Col, S):-
  452	qp_along_line(Low, Hi, X),
  453	convert_to_col(X, Col, I, S).
  454
  455% ?- initial_bridge(0, 1, 4, Bridge), maplist(writeln, Bridge).
  456initial_bridge(I, Low, Hi, Bridge):- J is I + 1,
  457	findall([p(I, V)-p(J, V), p(J, V)-p(I, V)],
  458			between(Low, Hi, V),
  459			Bridge).
  460
  461
  462		/****************
  463		*     prune     *
  464		****************/
  465
  466% ?- zdd((qp_along_line(1, 13, X), length(X, L), convert_to_col(X, Y, 3),
  467%	prune_by_column(Y, 3, p(3,1)-p(3,1), Y0), length(Y, Ly))).
  468
  469prune_by_column(F, X, P, F0):- @(prune_by_column(F, X, P, F0)).
  470%
  471prune_by_column([], _, _, [], _).
  472prune_by_column([Q|F], X, P, [Q|F0], S):- admissible_qp(Q, X, P, S), !,
  473	prune_by_column(F, X, P, F0, S).
  474prune_by_column([_|F], X, P, F0, S):-
  475	prune_by_column(F, X, P, F0, S).
  476%
  477admissible_qp(Q, X, P):- @(admissible_qp(Q, X, P)).
  478%
  479admissible_qp(Q, X, P, S):- qp_list(Q, List, S),
  480	admissible_qp_list(List, X, P).
  481
  482% ?- admissible_qp_list([p(0,0)-p(2,1), p(3,1)-p(3,1)], 2, p(0,0)).
  483admissible_qp_list([], _, _):-!.
  484admissible_qp_list([P-p(X, _)|List], X, P):-!,
  485	admissible_qp_list(List, X, P).
  486admissible_qp_list([p(X, _)-p(X, _)|List], X, P):- !,
  487	admissible_qp_list(List, X, P).
  488admissible_qp_list([p(Z, _)-_|_], X, _):- Z>X, !.
  489admissible_qp_list([P-P|List], X, P):- admissible_qp_list(List, X, P).
  490
  491% ?- rect((qp_along_line(1,1, Ls), maplist_qp(Ls, Q))).
  492% ?- rect((qp_along_line(1,10, Ls), maplist_qp(Ls, Q))).
  493% ?- time(rect((qp_along_line(1, 15, Ls), maplist_qp(Ls, Q)))).
  494%@ % 99,784,517 inferences, 8.136 CPU in 8.371 seconds (97% CPU, 12264944 Lips)
  495
  496%
  497maplist_qp(X, Y):- @(maplist_qp(X, Y)).
  498%
  499maplist_qp([], [], _).
  500maplist_qp([L|Ls], [Q|Qs], S):-
  501	qp_list(Q, L, S),
  502	maplist_qp(Ls, Qs, S).
  503
  504% simple line quasi path
  505qp_line(I, I, Q, S):-!, zdd_singleton(I, Q, S).
  506qp_line(I, J, Q, S):-
  507	M is (I+J)//2,
  508	M0 is M+1,
  509	qp_line(I, M, R, S),
  510	qp_line(M0, J, R0, S),
  511	bridge_line(M, M0, R, R0, Q0, S),
  512	bridge_line(M, M0, Q0, R0, Q1, S),
  513	zdd_join(R, R0, R1, S),
  514	zdd_join(R1, Q1, Q, S).
  515%
  516bridge_line(I, J, X, Y, Z):- @(bridge_line(I, J, X, Y, 0, Z)).
  517
  518bridge_line(_, _, X, 0, Z, Z, _):- X < 2, !.
  519bridge_line(_, _, 0, X, Z, Z, _):- X < 2, !.
  520bridge_line(I, J, X, Y, Z0, Z, S):- cofact(X, t(J0-K, L, R), S),
  521	( J0 < I
  522	->	Z = Z0
  523	;  	bridge_line(I, J, L, Y, Z0, Z1, S),
  524		(	J0 = J
  525		->	zdd_insert(I-K, R, R0, S),
  526			zdd_join(R0, Z1, Z, S)
  527		;	bridge_line(I, J, R, Y, 0, Z2, S),
  528			zdd_insert(J0-K, Z2, Z3, S),
  529			zdd_join(Z3, Z1, Z, S)
  530		)
  531	).
  532
  533
  534% ?- qp_along_line(1,1,X).
  535% ?- I = 6, qp_along_line(0,I,X), length(X, L), N is 2^(I+1).
  536% ?- I = 7, qp_along_line(0,I,X), length(X, L), N is 2^(I+1).
  537% ?- I = 20, qp_along_line(0,I,X), length(X, L), N is 2^(I+1).
  538%@ I = 20,
  539%@ X = [[0-20], [0-19, 20-20], [0-16, 17-20], [0-16, 17-19, 20-20], [0-17, 18-20], [0-17, 18-19, ... - ...], [0-16, ... - ...|...], [... - ...|...], [...|...]|...],
  540%@ L = 54608393,
  541%@ N = 2097152 .
  542% ?- qp_along_line(0,3,X), maplist(writeln, X).
  543
  544% simple line quasi path
  545qp_along_line(I, I, [[I-I]]).
  546qp_along_line(I, J, Q):-
  547	M is (I+J)//2,
  548	M0 is M+1,
  549	qp_along_line(I, M, R),
  550	qp_along_line(M0, J, R0),
  551	simple_prod_union(R, R0, R1, []),
  552	simple_bridge_links([M-M0, M0-M], R, R0, Q, R1).
  553
  554
  555%
  556map_qp_list([], X, X, _).
  557map_qp_list([L|Ls], X, Y, S):- qp_list(Q, L, S),
  558	zdd_join(Q, X, X0, S),
  559	map_pq_list(Ls, X0, Y, S).
  560
  561% ?- simple_bridge(1-2, [1-1], [2-2], R).
  562% ?- simple_bridge(2-1, [1-1], [2-2], R).
  563% ?- simple_bridge(2-3, [1-1], [2-2], R).
  564simple_bridge(A-B, L, M, N):-	% assuming qp  X < qp Y
  565	(	select(B-V, M, M0),
  566		select(U-A, L, L0)
  567	;	select(B-V, L, L0),
  568		select(U-A, M, M0)
  569	),
  570	!,
  571	(	U @> V
  572	-> 	ord_union(L0, [U-V|M0], N)
  573	; 	ord_union([U-V], M0, M1),
  574		ord_union(L0, M1, N)
  575	).
  576simple_bridge(_, _, _, []).
  577
  578%
  579simple_bridge_links([], _, _, Z, Z).
  580simple_bridge_links([E|Es], X, Y, Z, U):-
  581	simple_bridge_link(X, E, Y, Z, Z0),
  582	simple_bridge_links(Es, X, Y, Z0, U).
  583%
  584simple_bridge_link([], _, _, Z, Z).
  585simple_bridge_link([X|Xs], E, Y, Z, U):-
  586	simple_bridge_basic(Y, X, E, Z, Z0),
  587	simple_bridge_link(Xs, E, Y, Z0, U).
  588%
  589simple_bridge_basic([], _, _, Z, Z).
  590simple_bridge_basic([Y|Ys], X, E, Z, U):- simple_bridge(E, X, Y, V),
  591	(	V = []
  592	->	Z0 = Z
  593	;	Z  = [V|Z0]
  594	),
  595	simple_bridge_basic(Ys, X, E, Z0, U).
  596
  597% ?- Q = [a-a], Q0=[b-b], Q1 = [c-c],
  598%	simple_prod_union([Q, Q0], [Q1], Z).
  599
  600%
  601simple_prod_union(X, Y, Z):- simple_prod_union(X, Y, Z, []).
  602%
  603simple_prod_union([], _, U, U).
  604simple_prod_union([L|Y], Z, U, V):-
  605	simple_union(Z, L, U, U0),
  606	simple_prod_union(Y, Z, U0, V).
  607%
  608simple_union([], _, U, U).
  609simple_union([M|As], L, [N|U], V):-
  610	ord_union(L, M, N),
  611	simple_union(As, L, U, V).

?- zdd(( qp_list(X, [a-a]), qp_list(Y, [b-b]), qp_joint([a-b, b-a], [X], [Y], Z), maplist(pred(([U]:-qp_list(U, List), writeln(List))), Z))). ?- zdd(( qp_list(X, [a-a]), qp_list(Y, [b-b]), qp_joint([a-b, b-a], [X], [Y], Z), zdd_join(X, Y, A), memo(qp_suc(A)-L), writeln(qp_suc(A)-L), maplist(pred(([U]:-qp_list(U, List), writeln(List))), Z))). ?- zdd(( qp_list(X, [a-a]), qp_list(Y, [b-b]), qp_joint([a-b, b-a], [X], [Y], Z), zdd_join(X, Y, A), memo(qp_suc(A)-L), writeln(qp_suc(A)-L), maplist(pred(([U]:-qp_list(U, List), writeln(List))), Z))).

  631% qp_joint(X, Y, Z, U):-@(qp_joint(X, Y, Z, U)).
  632%
  633qp_joint(X, Y, Z, U, S):- map_prod_sum(Y, Z, V, [], S),
  634	qp_bridge_links(X, Y, Z, U, V, S).
  635
  636% ?- zdd((qp_list(X, [a-a]), qp_list(Y, [b-b]),
  637%	qp_bridge(a-b, X, Y, Z), qp_list(Z, U))).
  638% ?- zdd((qp_list(X, [a-a]), qp_list(Y, [b-b]),
  639%	qp_bridge(b-a, X, Y, Z), qp_list(Z, U))).
  640
  641qp_bridge(E, X, Y, Z):- @(qp_bridge(E, X, Y, Z)).
  642%
  643qp_bridge(A-B, X, Y, Z, S):- % assuming qp  X < qp Y
  644	qp_list(X, L, S),
  645	qp_list(Y, M, S),
  646	(	select(U-A, L, L0),
  647		select(B-V, M, M0)
  648	;	select(B-V, L, L0),
  649		select(U-A, M, M0)
  650	),
  651	!,
  652	ord_union(L0, M0, N0),
  653	ord_union([U-V], N0, N),
  654	qp_list(Z, N, S).
  655qp_bridge(_, _, _, 0, _).
  656
  657%
  658qp_bridge_links(X, Y, Z, U, V):- @(qp_bridge_links(X, Y, Z, U, V)).
  659%
  660qp_bridge_links([], _, _, Z, Z, _).
  661qp_bridge_links([E|Es], X, Y, Z, U, S):-
  662	qp_bridge_link(X, E, Y, Z, Z0, S),
  663	qp_bridge_links(Es, X, Y, Z0, U, S).
  664%
  665qp_bridge_link([], _, _, Z, Z, _).
  666qp_bridge_link([X|Xs], E, Y, Z, U, S):-
  667	qp_bridge_basic(Y, X, E, Z, Z0, S),
  668	qp_bridge_link(Xs, E, Y, Z0, U, S).
  669%
  670qp_bridge_basic([], _, _, Z, Z, _).
  671qp_bridge_basic([Y|Ys], X, E, Z, U, S):- qp_bridge(E, X, Y, V, S),
  672	(	V = 0
  673	->	Z0 = Z
  674	;	Z  = [V|Z0],
  675		select_goto(E, X, Y, V, S)
  676	),
  677	qp_bridge_basic(Ys, X, E, Z0, U, S).
  678%
  679select_goto(E, X, Y, Z):- @(select_goto(E, X, Y, Z)).
  680%
  681select_goto(A-B, X, Y, Z, S):-
  682	G = (A-B)-Z,
  683	( A@<B -> W = X ; W = Y ),
  684	getmemo(qp_suc(W)-U, V, S),
  685	(	var(V) -> V = [] ; true ),
  686	add_new(G, V, U).
  687
  688		/***********************
  689		*     end of bridge    *
  690		***********************/
  691
  692% ?-  belong(3, 3, 4).
  693% ?-  belong(3, 4, 1).
  694belong(X, I, J):- I @=< J, !, I @=< X, X @=< J.
  695belong(X, I, J):- J @=< X, X @=< I.
  696%
  697sub_interval(I, J, X, Y):- belong(I, X, Y), belong(J, X, Y).
  698%
  699disjoint_interval(I, J, X, Y):- \+ belong(I, X, Y), \+ belong(J, X, Y).
  700
  701% ?- compatible_interval(1, 2, 3, 4).
  702% ?- compatible_interval(1, 3, 2, 4). % false
  703% ?- compatible_interval(2, 4, 5, 7).
  704% ?- compatible_interval(7, 5, 2, 4).
  705
  706compatible_interval(I, J, X, Y):- disjoint_interval(I, J, X, Y), !.
  707compatible_interval(I, J, X, Y):- sub_interval(I, J, X, Y), !.
  708compatible_interval(I, J, X, Y):- sub_interval(X, Y, I, J).
  709
  710% ?- rect_dg(rect(1,1), X, Y).
  711% ?- rect_dg(rect(0,1), X, Y).
  712% ?- rect_dg(rect(0,0), X, Y).
  713
  714rect_dg(rect(W, H), Nodes, Links):-
  715	rect_nodes(W, H, Nodes),
  716	rect_links(W, H, Nodes, Links).
  717% ?- rect_nodes(4,5,Ns), length(Ns, L).
  718
  719rect_nodes(W, H, Nodes):-
  720	findall(p(I, J),
  721			(	between(0, W, I),
  722				between(0, H, J)
  723			),
  724			Nodes0),
  725	sort(Nodes0, Nodes).
  726
  727%
  728rect_links(W, H, Ns, Links):-
  729	findall(P-Q,
  730			(	member(P, Ns),
  731				P = p(I, J),
  732				(	Q = p(I, J1),
  733					(	J1 is J-1,
  734						J1 >= 0
  735					;	J1 is J+1,
  736						J1 =< H
  737					)
  738				;   Q = p(I1, J),
  739					(	I1 is I-1,
  740						I1 >= 0
  741					;	I1 is I+1,
  742						I1 =< W
  743					)
  744				)
  745			),
  746			Links0),
  747	sort(Links0, Links).
  748
  749		/********************
  750		*   frontier based  *
  751		********************/
  752
  753% ?- zdd(path_count(rect(0,0), C)).
  754% ?- zdd(path_count(rect(1,0), C)).
  755% ?- zdd(path_count(rect(0,1), C)).
  756% ?- zdd(path_count(rect(1,1), C)).
  757% ?- zdd(path_count(rect(2, 1), C)).
  758% ?- zdd(path_count(rect(3, 1), C)).
  759% ?- zdd(path_count(rect(1, 3), C)).
  760% ?- zdd(path_count(rect(4, 1), C)).
  761% ?- zdd(path_count(rect(2, 3), C)).
  762% ?- zdd(path_count(rect(4,3), C)).
  763%@ C = 976 .
  764% ?- zdd(path_count(rect(4,4), C)).
  765%@ C= 8512
  766% ?- call_with_time_limit(360, zdd(path_count(rect(5, 5), C))).
  767
  768path_count(Rect, C):- @(path_count(Rect, C)).
  769
  770%
  771path_count(rect(0,0), 1, _):-!.
  772path_count(Rect, C, S):-
  773	rect_dg(Rect, Ns, Links),
  774	findall(A-A, member(A, Ns), Ids),
  775	Rect = rect(W, H),
  776	set_key(final, p(W, H), S),
  777	qp_list(Q, Ids, S),
  778	memo(qp_suc(Q)-[], S),
  779	frontier(Links, [Q], FinalQs, S),
  780	fold_qp_memo_final(FinalQs, S),
  781	path_tree(Q, Tree, S),
  782	card(Tree, C, S).
  783
  784%
  785frontier([], Qs, Qs, _).
  786frontier([E|Es], Qs, Qs0, S):-
  787	frontier_link(E, Qs, Qs1, Qs, S),
  788	prune_frontier(Qs1, E, Qs2, S),
  789	frontier(Es, Qs2, Qs0, S).
  790
  791%
  792frontier_link(_, [], Qs, Qs, _):-!.
  793frontier_link(E, [X|Qs], Qs0, Qs1, S):- X<2, !,
  794 	frontier_link(E, Qs, Qs0, Qs1, S).
  795frontier_link(E, [Q|Qs], Next, Qs1, S):-
  796	qp_update(E, Q, Q0, S),
  797	(	Q0 < 2
  798	->	Next = Qs0
  799	;	qp_memo_final(Q0, S),
  800		Next = [Q0|Qs0]
  801	),
  802	qp_update_succs(E-Q0, Q, S),
  803	frontier_link(E, Qs, Qs0, Qs1, S).
  804
  805% ?- zdd((qp_list(Q, [a]), qp_list(Q0, [b]), qp_list(Q1, [c]),
  806%	map_prod_sum([Q, Q0], [Q1], Z),
  807%	maplist(pred([A, M]:-qp_list(A, M)), Z, W))).
  808%
  809map_prod_sum(X, Y, Z):- map_prod_sum(X, Y, Z, []).
  810%
  811map_prod_sum(X, Y, Z, U):- @(map_prod_sum(X, Y, Z, U)).
  812
  813map_prod_sum([], _, U, U, _).
  814map_prod_sum([Q|Y], Z, U, V, S):- qp_list(Q, L, S),
  815	map_qp_sum(Z, L, U, U0, S),
  816	map_prod_sum(Y, Z, U0, V, S).
  817%
  818map_qp_sum([], _, U, U, _).
  819map_qp_sum([Q|As], L, [Q0|U], V, S):-qp_list(Q, M, S),
  820	ord_union(L, M, N),
  821	qp_list(Q0, N, S),
  822	map_qp_sum(As, L, U, V, S).
  823
  824
  825% ?- zdd((path_count(rect(1,0), C))).
  826% ?- zdd((path_count(rect(0,1), C))).
  827% ?- zdd((path_count(rect(1,1), C))).
  828% ?- zdd((path_count(rect(3,1), C))).
  829% ?- time(zdd((path_count(rect(1,8), C)))).   % take time.
  830%@ % 86,018,490 inferences, 10.604 CPU in 11.323
  831% ?- time(zdd((path_count(rect(8,1), C)))).
  832%@ % 9,254,518 inferences, 3.089 CPU in 4.850
  833% ?- zdd((path_count(rect(2,3), C))).
  834% ?- zdd((path_count(rect(3,2), C))).
  835% ?- time(zdd((path_count(rect(3,3), C)))).
  836%@ % 50,584,360 inferences, 7.362 CPU in 9.480
  837%@ C = 184 .
  838% ?- time(zdd((path_count(rect(4,4), C)))).
  839%@ % 9,474,429,287 inferences, 1256.779 CPU in 1433.734 seconds (88% CPU, 7538660 Lips)
  840%@ C = 8512 .
  841
  842qp_update(E, X, Y):- @(qp_update(E, X, Y)).
  843%
  844qp_update(_, X, 0, _):- X < 2.
  845qp_update(E, X, Y, S):- qp_list(X, U, S),
  846	memo(qp_update(E, X)-Y, S),
  847	(	nonvar(Y)  ->	true
  848	;	select_touch(E, [A, B], U, U0)
  849	->	connect_pairs(E, [A, B], W),
  850		ord_union([W], U0, V),
  851		qp_list(Y, V, S)
  852	;	Y = 0
  853	).
  854qp_update(_, _, 0,  _).
  855
  856% ?- zdd((qp_list(X, [a-a, b-b]), memo(qp_suc(X)-[]))).
  857qp_update_succs(G, Q):- @(qp_update_succs(G, Q)).
  858
  859% ?- zdd(qp_update_succs(a, 1)).
  860qp_update_succs(G, Q, S):- getmemo(qp_suc(Q)-U, V, S),
  861	G = E-Q0,
  862	(	var(V) -> U = [E-Q0]
  863	;	add_new(E-Q0, V, U)
  864	).
  865
  866qp_admissible(L):-
  867	forall((select(A-B, L, L0), select(C-D, L0, _)),
  868		   compatible_interval(A, B, C, D)).
  869
  870% ?- compatible_interval(1, 2, 3, 4).
  871% ?- compatible_interval(1, 3, 2, 4). % false
  872% ?- compatible_interval(2, 4, 5, 7).
  873% ?- compatible_interval(7, 5, 2, 4).
  874
  875
  876% compatible_interval(I, J, X, Y):- disjoint_interval(I, J, X, Y), !.
  877% compatible_interval(I, J, X, Y):- sub_interval(I, J, X, Y), !.
  878% compatible_interval(I, J, X, Y):- sub_interval(X, Y, I, J).
  879%@ true.
  880%@ path-count4.pl qcompiled.
  881
  882% ?- select_touch(a-b, [A, B], [a-a, b-b], X).
  883% ?- select_touch(a-b, [A, B], [a-a, b-b], X),
  884%	connect_pairs(a-b, [A, B], Y).
  885%@ false.
  886
  887% ?- select_touch(a-b, [A, B], [a-a, b-c], X),
  888%	connect_pairs(a-b, [A, B], Y).
  889% ?- select_touch(a-b, [A, B], [a-a, c-b], X).
  890
  891select_touch(_, [], V, V).
  892select_touch(E, [A|As], [A|U], V):- touch(E, A), !,
  893	select_touch(E, As, U, V).
  894select_touch(E, As, [A|U], [A|V]):-select_touch(E, As, U, V).
  895
  896% ?- add_new(a,[b,c,a], X).
  897add_new(X, Y, Y):- memberchk(X, Y), !.
  898add_new(X, Y, [X|Y]).
  899
  900% ?- zdd((qp_list(Q, [p(0,0)-p(1,0)]), qp_list(Q, L))).
  901%@ L = [p(0, 0)-p(1, 0)].
  902%
  903qp_list(Q, List):- @(qp_list(Q, List)).
  904%
  905qp_list(0, 0, _):-!.
  906qp_list(1, [], _):-!.   

 ???

  907qp_list(Q, List, S):- var(Q), !,
  908	zdd_insert_atoms(List, 1, Q, S),
  909	memo(qp_suc(Q)-U, S),
  910	(	var(U) -> U = []
  911	;	true			% Already registered.
  912	).
  913qp_list(Q, [A|As], S):- cofact(Q, t(A, 0, R), S),
  914	qp_list(R, As, S).
  915
  916%
  917qp_memo_final(X):- @(qp_memo_final(X)).
  918%
  919qp_memo_final(Q, S):-  memoq(qp_final(Q)-true, S), !.
  920qp_memo_final(Q, S):- qp_list(Q, List, S),
  921	get_key(final, F, S),
  922	(	qp_final(List, F)
  923	->	memo(qp_final(Q)-true, S)
  924	;	true
  925	).
  926%
  927fold_qp_memo_final([], _).
  928fold_qp_memo_final([Q|Qs], S):- qp_memo_final(Q, S),
  929	fold_qp_memo_final(Qs, S).
  930
  931%
  932path_tree(Q, T):- @(path_tree(Q, T)).
  933%
  934
  935path_tree(X, 0,  _):- X<1, !.
  936path_tree(Q, 1, S):- memoq(qp_final(Q)-true, S), !.
  937path_tree(Q, T, S):- memo(path_tree(Q)-T, S),
  938	(	nonvar(T) -> true
  939	;	memo(qp_suc(Q) - L, S),		% writeln(qp_suc(Q)-L),
  940		path_tree_list(L, 0, T, S)
  941	).
  942%
  943path_tree_list([], C, C, _).
  944path_tree_list([E-Q|Qs], C0, C, S):-
  945	path_tree(Q, A, S),
  946	zdd_insert(E, A, A0, S),
  947	zdd_join(C0, A0, C1, S),
  948	path_tree_list(Qs, C1, C, S).
  949
  950% ?- qp_final([p(0,0)-p(2,3),p(0,3)-p(0,3),p(1,3)-p(1,3)], p(2,3)).
  951
  952qp_final([X-Y|R], X, Y):-!, qp_all_id(R).
  953qp_final([A-A|R], X, Y):- qp_final(R, X, Y).
  954
  955qp_final([p(0,0)-F|R], F):-!, qp_all_id(R).
  956qp_final([A-A|R], F):- qp_final(R, F).
  957
  958%
  959qp_all_id([]).
  960qp_all_id([A-A|R]):- qp_all_id(R).
  961
  962%
  963dag_card(Q, C):- @(dag_card(Q, C)).
  964%
  965dag_card(Q, C, S):- open_state(M),
  966	dag_card(Q, C, S, M),
  967	close_state(M).
  968%
  969dag_card(Q, 0,  _, _):- Q<2, !.
  970dag_card(Q, 1, S, _):- memoq(qp_final(Q)-true, S), !.
  971dag_card(Q, C, S, M):- memo(dag_card(Q)-C, M),
  972	(	nonvar(C) -> true
  973	;	memo(qp_suc(Q) - L, S),
  974		dag_card_list(L, 0, C, S, M)
  975	).
  976%
  977dag_card_list([], C, C, _, _).
  978dag_card_list([_-Q|Qs], C0, C, S, M):-
  979	dag_card(Q, A, S, M),
  980	C1 is C0 + A,
  981	dag_card_list(Qs, C1, C, S, M).
  982
  983
  984		/***************
  985		*     debug    *
  986		***************/
  987%
  988snap_qp_update(E, X, Y, S):- qp_update(E, X, Y, S),
  989	qp_list(X, X0, S),
  990	qp_list(Y, Y0, S),
  991	format("~w\n~w\n~w\n\n", [qp_update(E, X, Y), X0, Y0]).
  992%
  993dump_frontier(Qs, S):- maplist(pred(S,
  994	([Q]:- qp_list(Q, U, S), writeln(U))), Qs).
  995%
  996dump_qp(X, _):-  X<2, !, writeln(X).
  997dump_qp(X, S):-  qp_list(X, L, S),
  998	writeln(X=L).
  999
 1000%
 1001delta_star([], Q, Q, _).
 1002delta_star([E|Es], Q,  Q0, S):-  memo(qp_suc(Q)-L, S),
 1003	writeln(E),
 1004	memberchk(E-Q1, L),
 1005	dump_qp(Q1, S),
 1006	delta_star(Es, Q1, Q0, S).
 1007
 1008%
 1009mirror(X, Y):- @(mirror(X, Y)).
 1010%
 1011mirror(X, X, _):- X<2, !.
 1012mirror(X, Y, S):- cofact(X, t(p(I,J)-p(K,L), Lx, Rx), S),
 1013	mirror(Lx, MLx, S),
 1014	mirror(Rx, MRx, S),
 1015	zdd_insert(p(J, I)-p(L, K), MRx, MRx0, S),
 1016	zdd_join(MLx, MRx0, Y, S).
 1017
 1018%
 1019rect_path_tree(Rect, T):-@(rect_path_tree(Rect, T)).
 1020
 1021rect_path_tree(Rect, T, S):-
 1022	sort_links_by_rank(Rect, Rs, Ids),
 1023	Rect = rect(W, H),
 1024	set_key(final, p(W, H), S),
 1025	qp_list(Q, Ids, S),
 1026	repeat_frontier(Rs, [Q], FinalQs, S),
 1027	fold_qp_memo_final(FinalQs, S),
 1028	path_tree(Q, T, S).
 1029
 1030% ?- spy(remove_end_links).
 1031% ?- remove_end_links([1-[c-a]], R0, a-d).
 1032% ?- remove_end_links([1-[c-a, d-c, a-d]], R0, a-d).
 1033
 1034remove_end_links(R, R0, ST):- maplist(pred(ST,
 1035	(	[J-Es, J-Es0]:- foldl(
 1036		pred([ST, Es, Es0],
 1037			 ( [P-Q, U, V]:- ST = S-T,
 1038					(	( Q==S; P==T)
 1039					->  V = U
 1040					;	U = [P-Q|V]
 1041					))),  Es,  Es0, []))),  R, R0).
 1042
 1043%
 1044rest(Q, X, Y, S):- memo(qp_suc(Q)-L, S),
 1045	(	select(Q, X, X0)-> true
 1046	;	X0 = X
 1047	),
 1048	rest_list(L, X0, Y, S).
 1049%
 1050rest_list([], X, X, _).
 1051rest_list([_-Q|Qs], X, Y, S):-
 1052	rest(Q, X, X0, S),
 1053	rest_list(Qs, X0, Y, S).
 1054
 1055
 1056% ?- zdd((dbg_qp_update(p(0,0)-p(1,0), [p(0,0)-p(0,0),p(1,0)-p(1,0)], L))).
 1057% ?- spy(dbg_qp_update).
 1058% ?- zdd((dbg_qp_update(p(2,2)-p(2,3), [p(0,0)-p(2,2),p(0,3)-p(0,3), p(1,3)-p(1,3), p(2,3)-p(2,3)], L))).
 1059
 1060dbg_qp_update(E, X, Y):- @(dbg_qp_update(E, X, Y)).
 1061%
 1062dbg_qp_update(E, X, Y, S):- qp_list(Q, X, S),
 1063	qp_update(E, Q, Q0, S),
 1064	qp_list(Q0, Y, S).
 1065
 1066
 1067% ?- zdd((rect_path_tree(rect(3,2), T), card(T, C),
 1068%   zdd_insert_atoms([p(0,0)-p(2,3),p(0,3)-p(0,3),p(1,3)-p(1,3)], 1, W),
 1069%	mirror(T, MT), card(MT, CMT),
 1070%	rect_path_tree(rect(2,3), T0), card(T0, C0), mirror(T0, MT0),
 1071%	card(MT0, CM0),
 1072%	zdd_subtr(T, MT0, D), card(D, CD), psa(D)
 1073%)).
 1074
 1075
 1076% Test to show getmemo/2,3 works well.
 1077% ?- time((numlist(1, 10000, Ns),
 1078%	zdd((maplist( pred([A]:- getmemo(a(A)-A, _)), Ns),
 1079%		 maplist( pred(([A]:- getmemo(a(A)-B, Old), B is 2*Old)),  Ns),
 1080%		 foldl( pred(([I, B, C]:- getmemo(a(I)-J, J),
 1081%							 C is B + J)), Ns, 0, Sum))))).
 1082%@ % 938,719 inferences, 0.074 CPU in 0.075 seconds (99% CPU, 12701698 Lips)
 1083%@ Ns = [1, 2, 3, 4, 5, 6, 7, 8, 9|...],
 1084%@ Sum = 100010000.