1module triang;
2
3% Redistribution and use in source and binary forms, with or without
4% modification, are permitted provided that the following conditions are met:
5%
6%    * Redistributions of source code must retain the relevant copyright
7%      notice, this list of conditions and the following disclaimer.
8%    * Redistributions in binary form must reproduce the above copyright
9%      notice, this list of conditions and the following disclaimer in the
10%      documentation and/or other materials provided with the distribution.
11%
12% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
13% AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
14% THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
15% PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNERS OR
16% CONTRIBUTORS
17% BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
18% CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
19% SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
20% INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
21% CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
22% ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
23% POSSIBILITY OF SUCH DAMAGE.
24%
25
26
27COMMENT
28
29                ##########################################
30                ##                                      ##
31                ##    Solving zerodimensional systems   ##
32                ##          Triangular systems          ##
33                ##                                      ##
34                ##########################################
35
36
37Zerosolve returns lists of dpmats in prefix form, that consist of
38triangular systems in the sense of Lazard, provided the input is
39radical. For the corresponding definitions and concepts see
40
41[Lazard] D. Lazard: Solving zero dimensional algebraic systems.
42                J. Symb. Comp. 13 (1992), 117 - 131.
43and
44
45[EFGB] H.-G. Graebe: Triangular systems and factorized Groebner
46                bases. Report Nr. 7 (1995), Inst. f. Informatik,
47                        Univ. Leipzig.
48
49
50The triangularization of zerodim. ideal bases is done by Moeller's
51approach, see
52
53[Moeller] H.-M. Moeller : On decomposing systems of polynomial
54        equations with finitely many solutions.
55        J. AAECC 4 (1993), 217 - 230.
56
57We present three implementations :
58        -- the pure lex gb (zerosolve)
59        -- the "slow turn to pure lex" (zerosolve1)
60and
61        -- the mix with [FGLM] (zerosolve2)
62
63END COMMENT;
64
65symbolic procedure triang!=trsort(a,b);
66  mo_dlexcomp(dp_lmon a,dp_lmon b);
67
68symbolic procedure triang!=makedpmat x;
69  makelist for each p in x collect dp_2a p;
70
71% =================================================================
72% The pure lex approach.
73
74symbolic operator zerosolve;
75symbolic procedure zerosolve m;
76  if !*mode='algebraic then makelist zerosolve!* dpmat_from_a m
77  else zerosolve!* m;
78
79symbolic procedure zerosolve!* m;
80% Solve a zerodimensional dpmat ideal m, first groebfactor it and then
81% triangularize it. Returns a list of dpmats in prefix form.
82 if (dpmat_cols m>0) or (dim!* m>0) then
83        rederr"ZEROSOLVE only for zerodimensional ideals"
84 else if not !*noetherian or ring_degrees cali!=basering then
85        rederr"ZEROSOLVE only for pure lex. term orders"
86 else for each x in groebfactor!*(m,nil) join triang_triang car x;
87
88symbolic procedure triang_triang m;
89% m must be a zerodim. ideal gbasis (recommended to be radical)
90% wrt. a pure lex term order.
91% Returns a list l of dpmats in triangular form.
92 if (dpmat_cols m>0) or (dim!* m>0) then
93        rederr"Triangularization only for zerodimensional ideals"
94 else if not !*noetherian or ring_degrees cali!=basering then
95        rederr"Triangularization only for pure lex. term orders"
96 else for each x in triang!=triang(m,ring_names cali!=basering) collect
97                triang!=makedpmat x;
98
99symbolic procedure triang!=triang(a,vars);
100% triang!=triang(A,vars)={f1.x for x in triang!=triang(B,cdr vars)}
101%                       \union triang!=triang(A:<B>,vars)
102% where A={f1,...,fr}, B={f2~,...fr~}, see [Moeller].
103% Returns a list of polynomial lists.
104  if dpmat_unitideal!? a then nil
105  else begin scalar x,f1,m1,m2,b;
106    x:=car vars;
107    m1:=sort(for each x in dpmat_list a collect bas_dpoly x,
108                function triang!=trsort);
109    if length m1 = length vars then return {m1};
110    f1:=car m1;
111    m2:=for each y in cdr m1 collect bas_make(1,dp_xlt(y,x));
112    b:=interreduce!* dpmat_make(length m2,0,m2,nil,nil);
113    return append(
114    for each u in triang!=triang(b,cdr vars) collect (f1 . u),
115                triang!=triang(matstabquot!*(a,b),vars));
116    end;
117
118% =================================================================
119% Triangularization wrt. an arbitrary term order
120
121symbolic operator zerosolve1;
122symbolic procedure zerosolve1 m;
123  if !*mode='algebraic then makelist zerosolve1!* dpmat_from_a m
124  else zerosolve1!* m;
125
126symbolic procedure zerosolve1!* m;
127   for each x in groebfactor!*(m,nil) join triang_triang1 car x;
128
129symbolic procedure triang_triang1 m;
130% m must be a zerodim. ideal gbasis (recommended to be radical)
131% Returns a list l of dpmats in triangular form.
132 if (dpmat_cols m>0) or (dim!* m>0) then
133        rederr"Triangularization only for zerodimensional ideals"
134 else if not !*noetherian then
135        rederr"Triangularization only for noetherian term orders"
136 else for each x in triang!=triang1(m,ring_names cali!=basering) collect
137                triang!=makedpmat x;
138
139symbolic procedure triang!=triang1(a,vars);
140% triang!=triang(A,vars)={f1.x for x in triang!=triang1(B,cdr vars)}
141%                       \union triang!=triang1(A:<B>,vars)
142% where A={f1,...,fr}, B={f2~,...fr~}, see [Moeller].
143% Returns a list of polynomial lists.
144  if dpmat_unitideal!? a then nil
145  else if length vars = 1 then {{bas_dpoly first dpmat_list a}}
146  else (begin scalar u,x,f1,m1,m2,b,vars1,res;
147    x:=car vars; vars1:=ring_names cali!=basering;
148    setring!* ring_define(vars1,eliminationorder!*(vars1,{x}),
149                        'revlex,ring_ecart cali!=basering);
150    a:=groebfactor!*(dpmat_neworder(a,nil),nil);
151    % Constraints in dimension zero may be skipped :
152    a:=for each x in a collect car x;
153    for each u in a do
154    << m1:=sort(for each x in dpmat_list u collect bas_dpoly x,
155                function triang!=trsort);
156       f1:=car m1;
157       m2:=for each y in cdr m1 collect bas_make(1,dp_xlt(y,x));
158       b:=interreduce!* dpmat_make(length m2,0,m2,nil,nil);
159       res:=nconc(append(
160       for each v in triang!=triang1(b,cdr vars) collect (f1 . v),
161                triang!=triang1a(matstabquot!*(u,b),vars)),res);
162    >>;
163    return res;
164    end) where cali!=basering=cali!=basering;
165
166symbolic procedure triang!=triang1a(a,vars);
167% triang!=triang(A,vars)={f1.x for x in triang!=triang1(B,cdr vars)}
168%                       \union triang!=triang1(A:<B>,vars)
169% where A is already a gr basis wrt. the elimination order.
170% Returns a list of polynomial lists.
171  if dpmat_unitideal!? a then nil
172  else if length vars = 1 then {{bas_dpoly first dpmat_list a}}
173  else begin scalar u,x,f1,m1,m2,b;
174    x:=car vars;
175    m1:=sort(for each x in dpmat_list a collect bas_dpoly x,
176                function triang!=trsort);
177    f1:=car m1;
178    m2:=for each y in cdr m1 collect bas_make(1,dp_xlt(y,x));
179    b:=interreduce!* dpmat_make(length m2,0,m2,nil,nil);
180    return append(
181       for each u in triang!=triang1(b,cdr vars) collect (f1 . u),
182                triang!=triang1a(matstabquot!*(a,b),vars));
183    end;
184
185% =================================================================
186% Triangularization wrt. an arbitrary term order and FGLM approach.
187
188symbolic operator zerosolve2;
189symbolic procedure zerosolve2 m;
190  if !*mode='algebraic then makelist zerosolve2!* dpmat_from_a m
191  else zerosolve2!* m;
192
193symbolic procedure zerosolve2!* m;
194% Solve a zerodimensional dpmat ideal m, first groebfactoring it and
195% secondly triangularizing it.
196   for each x in groebfactor!*(m,nil) join triang_triang2 car x;
197
198symbolic procedure triang_triang2 m;
199% m must be a zerodim. ideal gbasis (recommended to be radical)
200% Returns a list l of dpmats in triangular form.
201 if (dpmat_cols m>0) or (dim!* m>0) then
202        rederr"Triangularization only for zerodimensional ideals"
203 else if not !*noetherian then
204        rederr"Triangularization only for noetherian term orders"
205 else for each x in triang!=triang2(m,ring_names cali!=basering)
206        collect triang!=makedpmat x;
207
208symbolic procedure triang!=triang2(a,vars);
209% triang!=triang(A,vars)={f1.x for x in triang!=triang2(B,cdr vars)}
210%                       \union triang!=triang2(A:<B>,vars)
211% where A={f1,...,fr}, B={f2~,...fr~}, see [Moeller].
212% Returns a list of polynomial lists.
213  if dpmat_unitideal!? a then nil
214  else if length vars = 1 then {{bas_dpoly first dpmat_list a}}
215  else (begin scalar u,x,f1,m1,m2,b,vars1,vars2,extravars,res,c1;
216    x:=car vars; vars1:=ring_names cali!=basering;
217    extravars:=dpmat_from_a('list . (vars2:=setdiff(vars1,vars)));
218    % We need this to make A truely zerodimensional.
219    c1:=ring_define(vars1,eliminationorder!*(vars1,{x}),
220                        'revlex,ring_ecart cali!=basering);
221    a:=matsum!* {extravars,a};
222    u:=change_termorder!*(a,c1);
223    a:=groebfactor!*(dpmat_sieve(u,vars2,nil),nil);
224    % Constraints in dimension zero may be skipped :
225    a:=for each x in a collect car x;
226    for each u in a do
227    << m1:=sort(for each x in dpmat_list u collect bas_dpoly x,
228                function triang!=trsort);
229       f1:=car m1;
230       m2:=for each y in cdr m1 collect bas_make(1,dp_xlt(y,x));
231       b:=interreduce!* dpmat_make(length m2,0,m2,nil,nil);
232       res:=nconc(append(
233       for each v in triang!=triang2(b,cdr vars) collect (f1 . v),
234                triang!=triang2a(matstabquot!*(u,b),vars)),res);
235    >>;
236    return res;
237    end) where cali!=basering=cali!=basering;
238
239symbolic procedure triang!=triang2a(a,vars);
240% triang!=triang(A,vars)={f1.x for x in triang!=triang2(B,cdr vars)}
241%                       \union triang!=triang2(A:<B>,vars)
242% where A is already a gr basis wrt. the elimination order.
243% Returns a list of polynomial lists.
244  if dpmat_unitideal!? a then nil
245  else if length vars = 1 then {{bas_dpoly first dpmat_list a}}
246  else begin scalar u,x,f1,m1,m2,b;
247    x:=car vars;
248    m1:=sort(for each x in dpmat_list a collect bas_dpoly x,
249                function triang!=trsort);
250    f1:=car m1;
251    m2:=for each y in cdr m1 collect bas_make(1,dp_xlt(y,x));
252    b:=interreduce!* dpmat_make(length m2,0,m2,nil,nil);
253    return append(
254       for each u in triang!=triang2(b,cdr vars) collect (f1 . u),
255                triang!=triang2a(matstabquot!*(a,b),vars));
256    end;
257
258endmodule; % triang
259
260end;
261