Artifact 1e77d65bc79d06a399dc5f034092a6c8b48f975f987d7f92c18ea7bcae85d061:
- Executable file
r38/packages/assist/dummycnt.red
— part of check-in
[f2fda60abd]
at
2011-09-02 18:13:33
on branch master
— Some historical releases purely for archival purposes
git-svn-id: https://svn.code.sf.net/p/reduce-algebra/code/trunk/historical@1375 2bfe0521-f11c-4a00-b80e-6202646ff360 (user: arthurcnorman@users.sourceforge.net, size: 48766) [annotate] [blame] [check-ins using] [more...]
module dummycnt; fluid '(g_dvnames g_dvbase g_sc_ve g_init_stree g_skip_to_level !*distribute); %%%%%%%%%%%%%%%%%%%%%%% MISCELANEOUS ROUTINES %%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure ad_splitname(u); if idp u then begin scalar uu, nn; uu := reverse(explode u); while uu and (charnump!: car uu) do << nn := car uu . nn; uu := cdr uu; >>; if uu then uu := intern(compress(reverse(uu))); if nn then nn := compress(nn); return (uu.nn); end; symbolic procedure anticom_assoc(u,v); begin scalar next_cell; if null v then return nil else if u = caar v then return (1 . car v) else << next_cell := anticom_assoc(u, cdr v); if null next_cell then return nil; if oddp(length(cdar v)) then rplaca(next_cell, - car next_cell); return next_cell; >>; end; %%%%%%%%%%%%%%%%%%%%%%% ORDERING FUNCTIONS %%%%%%%%%%%%%%%%%%%% symbolic procedure ad_signsort(l, fn); begin scalar tosort, sorted, insertl, dig; integer thesign; tosort := copy l; thesign := 1; sorted := nil; while tosort do if null sorted then << sorted := car tosort . sorted; tosort := cdr tosort; >> else if car tosort = car sorted then << thesign := 0; sorted := tosort := nil; >> else if apply(fn, {car sorted, car tosort}) then << sorted := car tosort . sorted; tosort := cdr tosort; >> else << thesign := - thesign; insertl := sorted; dig := t; while dig do if null cdr insertl then dig := nil else if cadr insertl = car tosort then << insertl := {nil}; dig := nil; thesign := 0; sorted := tosort := nil; >> else if not apply(fn, {cadr insertl, car tosort}) then << insertl := cdr insertl; thesign := - thesign; >> else dig := nil; if tosort then << rplacd(insertl, (car tosort) . cdr insertl); tosort := cdr tosort; >>; >>; return (thesign . reverse sorted); end; symbolic procedure cdr_sort(lst, fn); begin scalar tosort, sorted, insertl; tosort := lst; while tosort do << if (null sorted) or apply(fn, {cdar sorted, cdar tosort}) then << sorted := car tosort . sorted; tosort := cdr tosort; >> else << insertl := sorted; while (cdr insertl) and not(apply(fn, {cdadr insertl, cdar tosort})) do insertl := cdr insertl; rplacd(insertl, (car tosort) . cdr insertl); tosort := cdr tosort >> >>; return reverse sorted; end; symbolic procedure cdr_signsort(l, fn); begin scalar tosort, sorted, insertl, dig; integer thesign; tosort := copy l; thesign := 1; sorted := nil; while tosort do if null sorted then << sorted := car tosort . sorted; tosort := cdr tosort; >> else if cdar tosort = cdar sorted then << thesign := 0; sorted := tosort := nil; >> else if apply(fn, {cdar sorted, cdar tosort}) then << sorted := car tosort . sorted; tosort := cdr tosort; >> else << thesign := - thesign; insertl := sorted; dig := t; while dig do if null cdr insertl then dig := nil else if cdadr insertl = cdar tosort then << dig := nil; thesign := 0; sorted := tosort := nil; >> else if not apply(fn, {cdadr insertl, cdar tosort}) then << insertl := cdr insertl; thesign := - thesign; >> else dig := nil; if tosort then << rplacd(insertl, (car tosort) . cdr insertl); tosort := cdr tosort >>; >>; return (thesign . reverse sorted); end; symbolic procedure num_signsort(l); ad_signsort(l, function(lambda(x,y); x <= y)); symbolic procedure cons_ordp(u,v, fn); if null u then t else if null v then nil else if pairp u then if pairp v then if car u = car v then cons_ordp(cdr u, cdr v, fn) else cons_ordp(car u, car v, fn) else nil else if pairp v then t else apply2(fn,u,v); symbolic procedure atom_compare(u,v); if numberp u then numberp v and not(u < v) else if idp v then orderp(u,v) else numberp v; symbolic procedure idcons_ordp(u,v); cons_ordp(u, v, function atom_compare); symbolic procedure skp_ordp(u,v); cons_ordp(car u, car v, function atom_compare); symbolic procedure numlist_ordp(u,v); cons_ordp(u,v,function(lambda(x,y); x <= y)); symbolic procedure ad_numsort(l); sort(l,function(lambda(x,y); x <= y)); %%%%%%%%%%%%%%%%%%%%%%% ACCESS ROUTINES %%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure sc_kern(ind); caddr venth(g_sc_ve, ind); symbolic procedure sc_rep(ind); cadr venth(g_sc_ve, ind); symbolic procedure sc_desc(ind); car venth(g_sc_ve, ind); symbolic procedure dummyp(var); begin scalar varsplit; integer count, res; if not idp var then return nil; count := 1; while count <= upbve(g_dvnames) do << if var = venth(g_dvnames, count) then << res := count; count := upbve(g_dvnames) + 1 >> else count := count + 1; >>; if res eq 0 then << varsplit := ad_splitname(var); if (car varsplit eq g_dvbase) then return cdr varsplit >> else return res; end; symbolic procedure dv_ind2var(ind); if ind <= upbve(g_dvnames) then venth(g_dvnames, ind) else mkid(g_dvbase, ind); %%%%%%%%%%%%%%%%%%%%%% SYMMETRY CELLS %%%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure sc_repkern(s_cell, n); if car s_cell eq '!* then % nil symmetric cell begin scalar kern, rest, next_rest; integer head, rep; rest := cdr s_cell; rep := 0; while rest do << head := car rest; kern := {head} . kern; rest := cdr rest; next_rest := nil; rep := rep*2 + 1; for each elt in rest do << if elt eq head then rep := rep * 2 + 1 else << rep := rep * 2; next_rest := elt . next_rest >> >>; rest := reverse next_rest; >>; return {rep, pa_list2vect(reverse kern, n)}; end else begin scalar count, replist, rep, kern; integer last_count; s_cell := cdr s_cell; % s_cell supposed sorted for each elt in s_cell do if (count := assoc(elt, replist)) then rplacd (count, cdr count + 1) else replist := (elt . 1) . replist; replist := sort(replist, function(lambda(x,y); cdr x <= cdr y)); last_count := 0; for each elt in replist do if (cdr elt neq last_count) then << rep := (cdr elt . 1) . rep; kern := {car elt} . kern; last_count := cdr elt; >> else << rplacd(car rep, cdar rep + 1); rplaca(kern, car elt . car kern) >>; return {rep , pa_list2vect(kern, n)}; end; %%%%%%%%%%%%%%%%%%%%% PARTITIONS COMP %%%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure pa_list2vect(pa, n); begin scalar ve, reps; integer abs; ve := mkve(n); for each cell in pa do << reps := eval('min . cell) . reps; for each elt in cell do putve(ve, elt, car reps); >>; for count := 1:n do << if null venth(ve, count) then << if abs = 0 then << abs := count; reps := abs . reps >>; putve(ve, count, abs) >> >>; return ((reverse reps) . ve); end; symbolic procedure pa_part2list(p); begin scalar ve; integer len, rep; len := upbve(cdr p); ve := mkve(len); for count := len step -1 until 1 do << rep := venth(cdr p, count); putve(ve, rep, count . venth(ve, rep)); >>; return for each count in car p join copy venth(ve,count); end; symbolic procedure pa_vect2list(pa); begin scalar ve; integer count, rep; ve := mkve(upbve(cdr pa)); for count := 1 : upbve(cdr pa) do << rep := venth(cdr pa, count); putve(ve, rep, count . venth(ve,rep)); >>; return for each rep in (car pa) collect ordn(venth(ve, rep)); end; symbolic procedure pa_coinc_split(p1, p2); begin scalar ve1, ve2, cursplit, split_alist, split_info, coinc, split; integer count, plength; plength := upbve(cdr p1); ve1 := mkve(plength); ve2 := mkve(plength); split := mkve(plength); count := 0; for each rep in car p1 do << count := count + 1; putve(ve1, rep, count) >>; count := 0; for each rep in car p2 do << count := count + 1; putve(ve2, rep, count) >>; for count := 1 : plength do << cursplit := (venth(ve1, venth(cdr p1, count)) . venth(ve2, venth(cdr p2, count))); if (split_info := assoc(cursplit, split_alist)) then << rplacd(cdr split_info, cddr split_info + 1); putve(split, count, cadr split_info) >> else << split_info := cursplit . (count . 1); split_alist := split_info . split_alist; putve(split, count, count) >> >>; split_alist := sort(split_alist, function(lambda x,y; if caar x < caar y then t else if caar y < caar x then nil else cdar x leq cdar y)); split := (for each cell in split_alist collect cadr cell) . split; coinc := for each cell in split_alist collect (car cell) . (cddr cell); return coinc . split; end; %%%%%%%%%%%%%%%%%%%%% SYMMETRY TREES %%%%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure st_flatten(stree); if numberp(cadr stree) then cdr stree else for each elt in cdr stree join copy st_flatten(elt); symbolic procedure st_extract_symcells(stree, maxind); begin scalar ve, symcells; integer count; if null stree then return (nil . mkve(0)); symcells := st_extract_symcells1(st_consolidate(stree),nil,1); stree := car symcells; if not listp stree then % stree is a single symcell stree := {'!* , stree}; symcells := cadr symcells; ve := mkve(length(symcells)); count := upbve(ve); while symcells do << putve(ve, count, car symcells . sc_repkern(car symcells, maxind)); symcells := cdr symcells; count := count - 1 >>; return(st_consolidate(stree) . ve); end; symbolic procedure st_extract_symcells1(stree, symcells, count); begin scalar res, new_stree; if not listp cadr stree then % stree is a symcell return { count , stree . symcells, count + 1} else << new_stree := car stree . for each inner_stree in cdr stree collect << res := st_extract_symcells1(inner_stree, symcells, count); symcells := cadr res; count := caddr res; if numberp car res then {'!*, car res} else car res >>; return ({ new_stree, symcells, count }) >>; end; symbolic procedure st_signchange(ve1, ve2); car st_signchange1(g_init_stree, vect2list ve1) * car st_signchange1(g_init_stree, vect2list ve2); symbolic procedure st_signchange1(stree, eltlist); begin scalar levlist, elt_levlist, subsign; integer the_sign; the_sign := 1; levlist := for each child in cdr stree collect if numberp child then child else << subsign := st_signchange1(child, eltlist); the_sign := the_sign * car subsign; cdr subsign >>; if not cdr levlist then return (the_sign . car levlist); elt_levlist := eltlist; if member(car eltlist, levlist) then elt_levlist := 0 . elt_levlist else while not member(cadr elt_levlist, levlist) do elt_levlist := cdr elt_levlist; %% cdr elt_levlist starts with the elements of levlist %% Compute the sign change if car stree eq '!- and not permp(levlist, cdr elt_levlist) then the_sign := - the_sign; %% remove from elt_levlist (and thus from eltlist) %% the elements of levlist except the last (which will be the %% ref). rplacd(elt_levlist, pnth(cdr elt_levlist, length(levlist))); return (the_sign . cadr elt_levlist); end; symbolic procedure st_sorttree(stree, ve, fn); cdr st_sorttree1(stree, ve, fn); symbolic procedure st_sorttree1(stree, ve, fn); begin scalar schild, vallist, sorted, thesign, tosort; thesign := 1; if numberp cadr stree then << if car stree eq '!* then << vallist := for each elt in cdr stree collect venth(ve,elt); return (vallist . (1 . stree)) >>; tosort := for each elt in cdr stree collect elt . venth(ve,elt); >> else << if (car stree) eq '!* then << for each child in cdr stree do if thesign neq 0 then << schild := st_sorttree1(child, ve, fn); thesign := thesign * cadr schild; vallist := (car schild) . vallist; sorted := (cddr schild) . sorted; >>; if thesign = 0 then return (nil . 0 . nil) else << sorted := reverse sorted; vallist := reverse vallist; return (vallist . (thesign . ('!* . sorted))); >> >>; for each child in cdr stree do if thesign neq 0 then << schild := st_sorttree1(child, ve, fn); thesign := thesign * cadr schild; tosort := ((cddr schild) . (car schild)) . tosort; >>; >>; if thesign = 0 then return (nil . (0 . nil)); if car stree = '!+ then tosort := cdr_sort(tosort, fn) else << tosort := cdr_signsort(tosort, fn); if car tosort = 0 then return (nil . (0 . nil)) else thesign := thesign * car tosort; tosort := cdr tosort; >>; % fill up return structures while tosort do << sorted := (caar tosort) . sorted; vallist := (cdar tosort) . vallist; tosort := cdr tosort; >>; sorted := (car stree) . reverse sorted; vallist := reverse(vallist); return (vallist . (thesign . sorted)); end; symbolic procedure st_ad_numsorttree(stree); begin scalar sorted; sorted := st_ad_numsorttree1(stree); return car sorted . cadr sorted; end; symbolic procedure st_ad_numsorttree1(stree); begin scalar subtree, contents, tosort; integer thesign; if numberp stree then return {1, stree, stree}; thesign := 1; if car stree eq '!* then << stree := '!* . for each elt in cdr stree collect << subtree := st_ad_numsorttree1(elt); thesign := thesign * car subtree; contents := cddr subtree . contents; cadr subtree >>; contents := ad_numsort(for each elt in contents join elt); return thesign . (stree . contents); >>; tosort := for each elt in cdr stree collect << subtree := st_ad_numsorttree1(elt); thesign := thesign * car subtree; cdr subtree >>; if car stree eq '!+ then << tosort := cdr_sort(tosort, function numlist_ordp); tosort := for each elt in tosort collect << contents := (cdr elt) . contents; car elt >>; contents := ad_numsort(for each elt in reverse contents join elt); return (thesign . (('!+ . tosort) . contents)); >>; if car stree eq '!- then << tosort := cdr_signsort(tosort, function numlist_ordp); thesign := car tosort; tosort := for each elt in cdr tosort collect << contents := (cdr elt) . contents; car elt >>; contents := ad_numsort(for each elt in reverse contents join elt); return (thesign . (('!- . tosort) . contents)); >>; end; symbolic procedure st_consolidate(stree); begin scalar join_cells, children, tmp; if null stree then return nil; if numberp cadr stree then return stree; join_cells := t; for each child in reverse(cdr stree) do << tmp := st_consolidate(child); if tmp then << if cddr tmp then join_cells := nil else tmp := {'!*, cadr tmp}; children := tmp . children; >>; >>; if children then << if null cdr children then return car children; if join_cells then children := for each elt in children collect cadr elt; return (car stree) . children >> else return nil; end; %%%%%%%%%%%%%%%%%%%%%% SKELETONS %%%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure dv_cambhead(camb); begin if listp camb then << if member(car camb, {'expt, 'minus}) then return dv_cambhead(cadr camb); if listp camb then return car camb; >>; end; symbolic procedure dv_skelhead(skelpair); dv_cambhead car(skelpair); symbolic procedure dv_skelsplit(camb); begin scalar skel, stree, subskels; integer count, ind, maxind, thesign; thesign := 1; if not listp camb then if (ind := dummyp(camb)) then return {1, ind, ('!~dv . {'!*, ind})} else return {1, 0, (camb . nil)}; stree := get(car camb, 'symtree); if not stree then << stree := for count := 1 : length(cdr camb) collect count; if flagp(car camb, 'symmetric) then stree := '!+ . stree else if flagp(car camb, 'antisymmetric) then stree := '!- . stree else stree := '!* . stree >>; subskels := mkve(length(cdr camb)); count := 0; for each arg in cdr camb do << count := count + 1; if listp arg then putve(subskels, count, (arg . nil)) else if (ind := dummyp(arg)) then << maxind := max(maxind, ind); putve(subskels, count, ('!~dv . {'!*, ind})) >> else putve(subskels, count, (arg . nil)); >>; stree := st_sorttree(stree, subskels, function skp_ordp); if stree and (car stree = 0) then return nil; thesign := car stree; skel := dv_skelsplit1(cdr stree, subskels); stree := st_consolidate(cdr skel); skel := (car camb) . car skel; return {thesign, maxind, skel . stree}; end; symbolic procedure dv_skelsplit1(stree, skelve); begin scalar cell_stree, child_list, cur_cell, dv_stree, part, skel, ve; integer count, len; if numberp cadr stree then << ve := skelve; child_list := cdr stree; skel := for each elt in cdr stree collect car venth(ve,elt); >> else << len := length(cdr stree); ve := mkve(len); count := len; for each child in reverse(cdr stree) do << putve(ve, count, dv_skelsplit1(child, skelve)); skel := car(venth(ve,count)) . skel; child_list := count . child_list; count := count - 1; >>; skel := for each elt in skel join copy elt; >>; %% if root of stree is * node, then %% no partition of children is necessary if car stree eq '!* then << for each elt in reverse(child_list) do if cdr venth(ve, elt) then dv_stree := cdr venth(ve, elt) . dv_stree; if length(dv_stree) = 1 then dv_stree := car dv_stree else if dv_stree then dv_stree := '!* . dv_stree; return (skel . dv_stree); >>; %% regroup children with equal skeletons for each elt in child_list do if null cur_cell then % new skeleton cur_cell := car venth(ve, elt) . {cdr venth(ve, elt)} else if (car venth(ve, elt)) = (car cur_cell) then rplacd(cur_cell, (cdr venth(ve,elt)) . cdr cur_cell) else << part := cur_cell . part; cur_cell := car venth(ve, elt) . {cdr venth(ve, elt)}; >>; part := cur_cell . part; %% prepend contribution of each cell to dv_stree %% note that cells of part are in reverse order, %% as are elements of each cell for each cell in part do if cdr cell then << cell_stree := car stree . reverse(cdr cell); dv_stree := cell_stree . dv_stree >>; %% now set type of dv_stree, if it has more than one element if length(dv_stree) neq 1 then dv_stree := '!* . dv_stree else dv_stree := car dv_stree; return skel . dv_stree; end; symbolic procedure nodum_varp u; % u is a list or an atom (index) or !~dv or !~dva % returns true if it is neither a list nor a dummy var % nor !~dv or !~dva. if listp u then t else if flagp(u,'dummy) or car ad_splitname u = g_dvbase or u member {'!~dv,'!~dva} then nil else t; symbolic procedure list_is_all_free u; % u is a list of indices % returns nil if there is at least one dummy index % or if one of them is !~dv or !~dva. if null u then t else if nodum_varp car u then list_is_all_free cdr u else nil; symbolic procedure dv_skelprod(sklist, maxind); % This is the corrected function for commuting % operators which do not depend on dummy variables. begin scalar skform, stree, symcells, skel, apair, anticom_alist, com_alist, noncom_alist, acom_odd, acom_even, idvect, varskel; integer the_sign, count; %% sort skeletons according to lexicograpical order of dv_skelhead, %% placing commuting factors before anticommuting factors the_sign := 1; for each skelpair in sklist do << skel := car skelpair; varskel:=if listp skel then if car skel neq 'expt then cdr skel; % else % if car skel neq 'expt then cdr skel; if flagp(dv_skelhead skelpair , 'anticom) then << if (apair := anticom_assoc(skel, anticom_alist)) then << if member(cdr skelpair, cddr apair) then the_sign := 0 else the_sign := the_sign * car apair; rplacd(cdr apair, (cdr skelpair) . (cddr apair)) >> else anticom_alist := (skel . {cdr skelpair}) . anticom_alist; >> else if flagp(dv_skelhead skelpair, 'noncom) then noncom_alist := (skel . {cdr skelpair}) . noncom_alist % we do not need the "else if" for commuting operators % if no dummy variable is involved: % else if null list_is_all_free varskel or atom skel then % if(apair := assoc(skel, com_alist)) then else if ( (null list_is_all_free varskel or atom skel) and (apair := assoc(skel, com_alist)) ) then rplacd (apair, (cdr skelpair) . (cdr apair)) % else nil else com_alist := (skel . {cdr skelpair}) . com_alist; >>; if the_sign = 0 then return nil; %% restore order of factors for each anticom cell anticom_alist := for each elt in anticom_alist collect (car elt) . reverse(cdr elt); %% sort com_alist com_alist := sort(com_alist, function(lambda(x,y); idcons_ordp(car x, car y))); %% sort anticom_alist, taking care of sign changes %% isolate even prod of anticoms from odd prod of anticoms for each elt in anticom_alist do if evenp(length(cdr elt)) then acom_even := elt . acom_even else acom_odd := elt . acom_odd; acom_even := sort(acom_even, function(lambda(x,y); idcons_ordp(car x, car y))); anticom_alist := ad_signsort(acom_odd, function(lambda(x,y); idcons_ordp(car x, car y))); the_sign := the_sign * car anticom_alist; anticom_alist := merge_list1(acom_even, cdr anticom_alist, function idcons_ordp); skform := append(com_alist, anticom_alist); skform := append(skform, reverse noncom_alist); if maxind = 0 then << if the_sign = -1 then skform := ((-1) . {nil}) . skform; return skform . nil; >>; %% build complete symtree, %% omiting skels which do not depend on dummy variables for each elt in reverse noncom_alist do stree := cadr elt . stree; for each elt in reverse anticom_alist do if length(cdr elt) > 1 then stree := ('!- . cdr elt) . stree else if (cdr elt) then stree := cadr elt . stree; for each elt in reverse com_alist do if length(cdr elt) > 1 then stree := ('!+ . cdr elt) . stree else if (cdr elt) then stree := cadr elt . stree; if length(stree) > 1 then stree := '!* . stree else stree := car stree; stree := st_consolidate(stree); idvect := mkve(maxind); for count := 1 : maxind do putve(idvect, count, count); stree := st_sorttree(stree, idvect, function numlist_ordp); %% the sign change for sorting the symmetry tree does not influence %% the sign of the expression. Indeed, the symtree used to fill up %% the blanks in the expression is the symtree stored with the %% skeleton, which is not sorted. Note however that if the sign here %% is 0, then the expression is null. % the_sign := the_sign * car stree; if car stree = 0 then return nil; if the_sign = -1 then skform := ((-1) . {nil}) . skform; symcells := st_extract_symcells(cdr stree, maxind); return skform . symcells; end; symbolic procedure dv_skel2factor1(skel_kern, dvars); begin scalar dvar,scr; if null skel_kern then return nil; return if listp skel_kern then <<scr:=dv_skel2factor1(car skel_kern, dvars); scr:=scr . dv_skel2factor1(cdr skel_kern, dvars) >> else if skel_kern eq '!~dv then << dvar := car dvars; if cdr dvars then << rplaca(dvars, cadr dvars); rplacd(dvars, cddr dvars); >>; dvar >> else skel_kern; end; %%%%%%%%%%%%%%% PARTITION SYMMETRY TREES %%%%%%%%%%%%%%%%%% symbolic procedure pst_termnodep(pst); null cdr venth(cdr pst, 1); symbolic procedure pst_mkpst(stree); pst_equitable(nil . pst_mkpst1(stree)); symbolic procedure st_multtermnodep(stree); begin scalar res, subtrees; if car stree neq '!* then return nil; subtrees := cdr stree; res := t; while subtrees do << if numberp cadar subtrees then subtrees := cdr subtrees else << subtrees := nil; res := nil; >> >>; return res; end; symbolic procedure pst_mkpst1(stree); begin scalar subtrees, s_cells, ve, pst, cell; integer count, lastcount; if null stree then return nil; ve := mkve(length(cdr stree)); subtrees := cdr stree; count := 1; if numberp(car subtrees) then % terminal node with single cell while subtrees do << putve(ve, count, ({car subtrees} . nil)); count := count + 1; subtrees := cdr subtrees; >> % check if valid as pst terminal node with several cells else if st_multtermnodep(stree) then << ve := mkve(for each cell in subtrees sum (length(cdr cell))); lastcount := 0; for each s_cell in subtrees do << cell := cdr s_cell; if car s_cell eq '!* then for count := 1 : length(cell) do pst := {count + lastcount} . pst else pst := ( for count := 1 : length(cell) collect (count + lastcount) ) . pst; count := lastcount + 1; lastcount := lastcount + length(cell); for each elt in cell do << putve(ve,count, {{elt}}); count := count + 1; >>; >>; return (reverse pst . ve); >> else while subtrees do << pst := pst_mkpst1(car subtrees); s_cells := nil; for count2 := 1 : upbve(cdr pst) do s_cells := append(car venth(cdr pst, count2), s_cells); putve(ve, count, (s_cells . pst)); count := count + 1; subtrees := cdr subtrees; >>; if ((car stree) eq '!*) then % discrete partition pst := ((for count := 1 : upbve(ve) collect {count}) . ve) else % single cell partition pst := ({(for count := 1 : upbve(ve) join {count})} . ve); return pst; end; symbolic procedure pst_subpst(pst, ind); venth(cdr pst, ind); symbolic procedure pst_reduce(pst); begin scalar isolated, f_cell, rpst, tmp, npart, nsubs; integer ind, count; if null pst then return (nil . nil); if null cdr pst then return pst; f_cell := caar pst; while length(f_cell) eq 1 do << ind := car f_cell; % index of pst_subpst if pst_termnodep(pst) then << isolated := append(isolated, {caar venth(cdr pst, ind)}); %% remove first cell from pst, and set f_cell if cdar pst then % pst is not fully reduced << %% remove first cell rplaca(pst, cdar pst); %% update pst representation npart := for each cell in car pst collect for each elt in cell collect if (elt > ind) then elt - 1 else elt; nsubs := mkve(upbve(cdr pst)-1); for count := 1 : upbve(nsubs) do if count geq ind then putve(nsubs, count, venth(cdr pst, count+1)) else putve(nsubs, count, venth(cdr pst, count)); rplaca(pst, npart); rplacd(pst, nsubs); f_cell := caar pst; >> else % pst fully reduced f_cell := pst := nil; >> else << rpst := pst_reduce(cdr pst_subpst(pst,ind)); if car rpst then %% new isolates << %% add new isolates to isolated isolated := append(isolated, car rpst); if cdr rpst then %% first subtree in pst was not discrete, update subtree spec << tmp := pst_subpst(pst,ind); rplaca(tmp, setdiff(car tmp, car rpst)); rplacd(tmp, cdr rpst); f_cell := nil; >> else % first subtree in pst was discrete, so remove it << if cdar pst then % pst not fully reduced << rplaca(pst, cdar pst); npart := for each cell in car pst collect for each elt in cell collect if (elt > ind) then elt - 1 else elt; nsubs := mkve(upbve(cdr pst)-1); for count := 1 : upbve(nsubs) do if count geq ind then putve(nsubs, count, venth(cdr pst, count+1)) else putve(nsubs, count, venth(cdr pst, count)); rplaca(pst, npart); rplacd(pst, nsubs); f_cell := caar pst; >> else f_cell := pst := nil; >>; >> else %% car rpst is nil, so no more isolated d-elts << f_cell := nil; >>; >> >>; return (isolated . pst); end; symbolic procedure pst_isolable(rpst); begin scalar ve, f_cell; %% verify if fully reduced. if null cdr rpst then return nil; %% f_cell is list of elts in first cell in rpst. %% ve is vector of descriptions of elts in f_cell f_cell := caadr rpst; ve := cddr rpst; %% if the elts in f_cell are d-elts, then return the list of d-elts if null cdr venth(ve, car f_cell) then return for each ind in f_cell collect caar venth(ve, ind); return for each ind in f_cell join copy pst_isolable(nil . cdr venth(ve, ind)); end; symbolic procedure pst_isolate(s_cell, rpst); begin scalar redisol; redisol := pst_reduce(pst_isolate1(s_cell, cdr rpst)); rplaca(redisol, append(car rpst, car redisol)); return redisol; end; symbolic procedure pst_isolate1(s_cell, pst); begin scalar fcell, tmp, spst; integer ind; %% fcell is the list of elts in the first cell of rpst %% ve is the vector of descriptions of elts in fcell fcell := caar pst; %% find out which elt of fcell needs to be set aside, if any tmp := fcell; while (ind = 0) do << if null tmp then ind := -1; ind := car tmp; tmp := cdr tmp; if not member(s_cell, car (spst := pst_subpst(pst, ind))) then ind := 0 >>; %% if no elt should be set aside, then s_cell is not isolable if (ind = -1) then return nil; %% effectively isolate, splitting first cell if necessary if (length(fcell) > 1) then << tmp := delete(ind, fcell) . cdar pst; tmp := {ind} . tmp; rplaca(pst, tmp) >>; %% if the set aside elt is not a mere dummy variable, then isolate %% s_cell in the partition it represents. if not pst_termnodep(pst) then << spst := car spst . pst_isolate1(s_cell, cdr spst); putve(cdr pst, ind, spst) >>; return pst; end; symbolic procedure pst_equitable(rpst); begin scalar nrpst, reduced, isol; if null cdr rpst then return rpst; isol := car rpst; nrpst := pst_reduce(cdr rpst); rplaca(nrpst, append(isol, car nrpst)); repeat << isol := car nrpst; nrpst := isol . pst_equitable1(isol, cdr nrpst); reduced := pst_reduce(cdr nrpst); if car reduced then nrpst := (append(isol, car reduced) . cdr reduced); reduced := car reduced >> until not reduced; return nrpst; end; symbolic procedure pst_equitable1(isolated, pst); begin scalar isol, ve, alpha, beta, p1, equit, cell, psi; integer len, k, n_delems; if null pst then return nil; %% make partition to equitate, merging isolated and car pst isol := isolated; len := length(isolated); ve := mkve(upbve(cdr pst) + len); for count := 1 : upbve(cdr pst) do putve(ve, count, car venth(cdr pst, count)); alpha := car pst; for count := upbve(cdr pst) + 1 : upbve(ve) do << putve(ve, count, {car isol}); isol := cdr isol; alpha := {count} . alpha; >>; p1 := fullcopy alpha; len := length(p1); n_delems := upbve(ve); while (alpha and len < n_delems) do << beta := car alpha; alpha := cdr alpha; equit := nil; len := 0; while(p1) do << cell := car p1; p1 := cdr p1; psi := if cdr cell then pst_partition(cell, beta, ve) else {cell}; k := length(psi); equit := append(equit, psi); len := len + k; if k geq 2 then alpha := append(cdr psi, alpha); >>; p1 := equit; >>; equit := pnth(p1,length(isolated)+1); %%% make every child of pst equitable w.r.t. isolated if not pst_termnodep(pst) then for count := 1 : upbve(cdr pst) do << p1 := venth(cdr pst, count); putve(cdr pst, count, (car p1 . pst_equitable1(isolated, cdr p1))); >>; return (equit . cdr pst); end; symbolic procedure pst_d1(d1,d2, ve); for each e1 in venth(ve,d1) collect ordn for each e2 in venth(ve, d2) collect ordn car pa_coinc_split(sc_kern(e1), sc_kern(e2)); symbolic procedure pst_d(d1, d2, ve); if listp d1 then if listp d2 then ordn for each e1 in d1 collect ordn for each e2 in d2 collect pst_d(e1, e2, ve) else ordn for each e1 in d1 collect pst_d(e1, d2, ve) else if listp d2 then ordn for each e2 in d2 collect pst_d(d1, e2, ve) else pst_d1(d1, d2, ve); symbolic procedure pst_partition(s1, s2, ve); begin scalar elt_d, elt_apair, pst_alist; for each elt in s1 do << elt_d := pst_d(elt, s2, ve); if (elt_apair := assoc(elt_d, pst_alist)) then rplacd(elt_apair, elt . cdr elt_apair) else pst_alist := (elt_d . {elt}) . pst_alist; >>; % sort regrouped elts according to distance to s2 pst_alist := sort(pst_alist, function( lambda(x,y); numlist_ordp(car x, car y))); return for each elt in pst_alist collect reverse(cdr elt); end; %%%%%%%%%%%%%%%%%%%%%%%% BACKTRACKING %%%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure dv_next_choice(sc, partial_perm, rpst, comp_info); begin scalar next_perm, extensions, nrpst, new_aut; integer npoints, len, next_img ; npoints := upbve(car sc); g_skip_to_level := len := upbve(partial_perm) + 1; sc_setbase(sc, partial_perm); extensions := pst_isolable(rpst); repeat << extensions := idsort( intersection (extensions, candidate_extensions(sc, partial_perm))); if extensions then << next_img := car extensions; extensions := cdr extensions; nrpst := pst_equitable(pst_isolate(next_img, fullcopy(rpst))); next_perm := list2vect!*(car nrpst,'symbolic); comp_info := dv_compare(next_perm, comp_info, len, npoints); if (car comp_info = 0) then if (upbve(next_perm) = npoints) then << new_aut := pe_mult(pe_inv(venth(cadr comp_info, 1)), next_perm); process_new_automorphism(sc, new_aut); >> else comp_info := dv_next_choice(sc, next_perm, nrpst, comp_info) else if (car comp_info = 1) then if (upbve(next_perm) < npoints) then comp_info := dv_next_choice(sc, next_perm, nrpst, comp_info) else rplaca(comp_info, 0); rplacd(cdr comp_info, cdr cddr comp_info); >> >> until (null extensions) or (len > g_skip_to_level); return comp_info; end; symbolic procedure can_rep_cell(comp_info, level); venth(venth(cadr comp_info, 2), level); symbolic procedure last_part_kern(comp_info); car cddr comp_info; symbolic procedure dv_compare(next_perm, comp_info, len, npoints); begin scalar part_kern, part_rep, can_rep, curlev, res; if car comp_info = 1 then return dv_fill_comp_info(next_perm, comp_info, len, npoints, nil, nil); if len = 1 then << part_kern := sc_kern(venth(next_perm, 1)); part_rep := {sc_rep(venth(next_perm,1))}; >> else << part_kern := last_part_kern(comp_info); part_kern := pa_coinc_split(part_kern, sc_kern(venth(next_perm, len))); part_rep := (sc_rep(venth(next_perm, len)) . car part_kern); part_kern := cdr part_kern; >>; can_rep := can_rep_cell(comp_info, len); curlev := len; res := 0; repeat << if equal(can_rep, part_rep) then << res := 0; if (curlev < upbve(next_perm)) then << curlev := curlev + 1; part_kern := pa_coinc_split(part_kern, sc_kern(venth(next_perm, curlev))); part_rep := (sc_rep(venth(next_perm,curlev)) . car part_kern); part_kern := cdr part_kern; can_rep := can_rep_cell(comp_info, curlev); >> >> else if numlist_ordp(can_rep, part_rep) then << res := 1; rplaca(comp_info, 1); comp_info := dv_fill_comp_info(next_perm, comp_info, curlev, npoints, part_rep, part_kern); >> else << res := 2; % grow partial permutation kernel stack rplacd(cdr comp_info, nil . (cddr comp_info)); rplaca(comp_info, 2); >> >> until (res neq 0) or (curlev = upbve(next_perm)); if res = 0 then << % update partial permutation stack rplacd(cdr comp_info, part_kern . cddr comp_info); if (curlev = npoints) and dv_new_aut_hook(next_perm, comp_info) then << g_skip_to_level := 0; rplaca(comp_info, 2); >>; >>; return comp_info; end; symbolic procedure dv_fill_comp_info(pe, comp_info, len, npoints, part_rep, part_kern); begin scalar part_rep; integer level; if len = 1 then << part_kern := sc_kern(venth(pe, 1)); part_rep := {sc_rep(venth(pe,1))}; >> else if null part_kern then << part_kern := last_part_kern(comp_info); part_kern := pa_coinc_split(part_kern, sc_kern(venth(pe, len))); part_rep := (sc_rep(venth(pe, len)) . car part_kern); part_kern := cdr part_kern; >>; putve(venth(cadr comp_info, 2), len, part_rep); level := len + 1; while(level <= upbve(pe)) do << part_kern := pa_coinc_split(part_kern, sc_kern(venth(pe, level))); part_rep := (sc_rep(venth(pe, level)) . car part_kern); part_kern := cdr part_kern; putve(venth(cadr comp_info, 2), level, part_rep); level := level + 1 >>; rplacd(cdr comp_info, part_kern . (cddr comp_info)); if level = npoints+1 then if null venth(cadr comp_info, 1) and dv_null_first_kern(part_kern) then << g_skip_to_level := 0; rplaca(comp_info, 2); >> else << putve(cadr comp_info, 1, fullcopy(pe)); putve(cadr comp_info, 3, part_kern); >>; return comp_info; end; symbolic procedure dv_null_first_kern(kern); begin scalar l_kern, cell, nullexp, acell; integer count, count2; nullexp := nil; l_kern := pa_vect2list kern; for each cell in l_kern do if cdr cell and not nullexp then << count := 0; for count2 := 1 : upbve(g_sc_ve) do if (car (acell := car venth(g_sc_ve, count2)) eq '!-) and member (car cell, acell) then count := count + 1; if oddp count then nullexp := t; >>; return nullexp; end; symbolic procedure dv_new_aut_hook(pe, comp_info); begin scalar tmp1, tmp2, ve; integer count, thesign; thesign := st_signchange(venth(cadr comp_info,1), pe); tmp1 := pa_part2list(venth(cadr comp_info, 3)); tmp2 := pa_part2list(caddr comp_info); ve := mkve(length(tmp1)); count := 1; while tmp1 do << putve(ve, car tmp1, car tmp2); tmp1 := cdr tmp1; tmp2 := cdr tmp2; count := count + 1; >>; for count := 1 : upbve(g_sc_ve) do << tmp1 := car venth(g_sc_ve, count); if car tmp1 eq '!- then << tmp1 := cdr tmp1; tmp2 := for each elt in tmp1 collect venth(ve,elt); % tmp2 is the image of tmp1. Since all cells in g_sc_ve are % ordered in increased numerical order thesign := thesign * car num_signsort(tmp2); >> >>; if thesign = -1 then return t; return nil; end; %%%%%%%%%%%%%%%%%%%%%%%%%% TOP LEVEL %%%%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure dv_canon_monomial(sf); begin scalar tmp, sklist, camb, skel, skprod, aut_sc, can_kern, new_dvnames, pst, comp_info, factlst, res, fact, sorted_factlst; integer count, expnt, thesign, maxind; %% get skeleton pairs for each one of the factors thesign := 1; while not domainp sf do << tmp := lpow(sf); sf := lc(sf); %% suppose exponents are integers expnt := cdr tmp; camb := car tmp; if expnt neq 1 and flagp(dv_cambhead(camb),'anticom) then << skel := nil; sf := nil; >> else skel := dv_skelsplit(camb); if null skel then sf := nil else << if car skel < 0 then << skel := cdr skel; if oddp(expnt) then thesign := - thesign else rplacd(cdr skel, subst('!-, '!+, cdr skel)); >> else skel := cdr skel; if (car skel > maxind) then maxind := car skel; skel := cadr skel; if expnt neq 1 then rplaca(skel, {'expt, car skel, expnt}); sklist := skel . sklist; >>; >>; if null sf then return nil; sklist := reverse((sf . nil) . sklist); %% regroup factors with identical skeletons skprod := dv_skelprod(sklist, maxind); if null skprod then return nil; sklist := car skprod; if maxind > 0 then << g_sc_ve := cddr skprod; g_init_stree := cadr skprod; aut_sc := sc_create(upbve(g_sc_ve)); comp_info := mkve(3); putve(comp_info, 2, mkve(upbve(g_sc_ve))); comp_info := {1, comp_info, nil}; pst := pst_mkpst(g_init_stree); tmp := list2vect!*(car pst,'symbolic); g_skip_to_level := 1; if car pst then comp_info := dv_compare(tmp, comp_info, 1, upbve(g_sc_ve)); if cdr pst then comp_info := dv_next_choice(aut_sc, tmp, pst, comp_info); if g_skip_to_level = 0 then return nil; can_kern := pa_part2list(venth(cadr comp_info, 3)); count := 0; new_dvnames := nil; for each elt in can_kern do << count := count + 1; if elt neq count then new_dvnames := (elt . count) . new_dvnames; >>; >>; for each cell in sklist do << factlst := nil; skel := car cell; if cadr cell then << for each stree in cdr cell do << fact := dv_skel2factor( (skel . stree), new_dvnames); if car fact = -1 then thesign := - thesign; factlst := (cdr fact) . factlst; >>; factlst := reverse factlst; if flagp(dv_cambhead skel, 'anticom) then << sorted_factlst := ad_signsort(factlst, 'idcons_ordp); thesign := thesign * car sorted_factlst; sorted_factlst := cdr sorted_factlst; >> else sorted_factlst := sort(factlst, 'idcons_ordp); res := append(res, sorted_factlst); >> else res := append(res, {skel}); >>; %% transform res, list of factors, into standard form if thesign = -1 then skprod := {'minus, 'times . res} else if thesign = 1 then skprod := 'times . res else skprod := 0; return !*a2f skprod; end; symbolic procedure dv_skel2factor(skelpair, newnames); begin scalar stree, dvars; if null cdr skelpair then return car skelpair; stree := sublis(newnames, cdr skelpair); stree := st_ad_numsorttree(stree); dvars := for each elt in st_flatten(cdr stree) collect dv_ind2var elt; return (car stree . dv_skel2factor1(car skelpair, dvars)); end; %%%%%%%%%%%%%%%%%%%%%%% USER INTERFACE %%%%%%%%%%%%%%%%%%%%%%%%% symbolic procedure canonical sq; begin scalar sf, denom, res, !*distribute; res := nil; sq := simp!* car sq; denom := denr sq; on distribute; sf := distri_pol numr sq; %% process each monomial in sf while not domainp(sf) do << res := addf(res, dv_canon_monomial(lt sf .+ nil)); sf := red sf; >>; res := addf(res,sf); %% simplify the whole thing, and return return simp!*( {'!*sq, res ./ denom, nil} ) end; put ('canonical, 'simpfn, 'canonical); flag('(symtree),'opfn); symbolic procedure symtree (name, s); << put (name, 'symtree, alg_to_symb s); >>; symbolic procedure remsym u; % ALLOWS TO ELIMINATE THE DECLARED SYMMETRIES. for each j in u do if flagp(j,'symmetric) then remflag(list j,'symmetric) else if flagp(j,'antisymmetric) then remflag(list j,'antisymmetric) else remprop(j,'symtree); symbolic procedure dummy_names u; <<if g_dvbase then msgpri("The created dummy base",g_dvbase, "must be cleared",nil,t); g_dvnames := list2vect!*(u,'symbolic); flag(u, 'dummy); t>>; rlistat '(dummy_names); symbolic procedure show_dummy_names; if g_dvnames then symb_to_alg vect2list g_dvnames else symb_to_alg list('list); symbolic procedure dummy_base u; if g_dvnames then msgpri("Named variables",symb_to_alg vect2list g_dvnames, "must be eliminated",nil,t) else g_dvbase := u; symbolic procedure clear_dummy_base; << g_dvbase := nil;t>>; symbolic procedure clear_dummy_names; << g_dvnames := nil;t>>; flag ('(show_dummy_names clear_dummy_names dummy_base clear_dummy_base), 'opfn); deflist( '((clear_dummy_base endstat) (clear_dummy_names endstat)),'stat); symbolic procedure anticom u; << for each x in u do <<flag(list x, 'anticom); flag(list x, 'noncom)>>; t>>; symbolic procedure remanticom u; % ALLOWS TO ELIMINATE THE DECLARED anticom flag. % Operators becomes COMMUTATIVE operators. << for each x in u do <<remflag(x,'noncom); remflag(x,'anticom)>>; t>>; deflist('((anticom rlis) (remanticom rlis)),'stat); endmodule; end;