source: git/Tst/BuchDL/Sol_Exe_4.tst @ 8f1d129

spielwiese
Last change on this file since 8f1d129 was 9558c5f, checked in by Hans Schönemann <hannes@…>, 18 years ago
*hannes/lossen: very long examples from DL-Book git-svn-id: file:///usr/local/Singular/svn/trunk@8762 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 6.8 KB
Line 
1LIB "tst.lib";
2tst_init();
3
4
5//======================  Exercise 4.1 =============================
6ring R = 0, (t,w,x,y,z), dp;
7ideal I = w2xy+w2xz+w2z2, tx2y+x2yz+x2z2, twy2+ty2z+y2z2,
8          t2wx+t2wz+t2z2;
9if(not(defined(primdecGTZ))){ LIB "primdec.lib"; }
10list L = primdecGTZ(I);
11
12int s = size(L);   // number of components
13int i,j;
14for (i=1; i<=s; i++)
15{
16  L[i][1]=std(L[i][1]);
17  L[i][2]=std(L[i][2]);
18  // test for primeness:
19  if ( size( reduce( lead(L[i][2]), std(lead(L[i][1])) ) ) > 0 )
20  {
21    print(string(i)+"th component is not prime and of dimension "
22      +string(dim(L[i][2]))+".");
23  }
24  else
25  {
26    print(string(i)+"th component is prime and of dimension "
27      +string(dim(L[i][2]))+".");
28  }
29}
30
31for (i=1; i<=s; i++)
32{
33  for (j=1; j<=s; j++)
34  {
35    if (( dim(L[j][2]) > dim(L[i][2]) ) and
36       ( size(reduce(L[j][2],L[i][2],1)) == 0 ))
37    {
38      print(string(i)+"th component is embedded.");
39      j=s;
40    }
41  }
42}
43
44
45kill R,s,i,j;
46//======================  Exercise 4.2 =============================
47ring R = 0, (b,s,t,u,v,w,x,y,z), dp;
48ideal I =  wy-vz, vx-uy, tv-sw, su-bv, tuy-bvz;
49if(not(defined(normal))){ LIB "normal.lib"; }
50list nor = normal(I);
51for (int i=1; i<=size(nor); i++) { def R_nor(i) = nor[i]; }
52setring R_nor(1);
53R_nor(1); norid; normap;
54setring R_nor(2);
55R_nor(2); norid; normap;
56setring R_nor(3);
57R_nor(3); norid; normap;
58
59setring R;
60primdecGTZ(I);
61
62
63kill R,R_nor(1),R_nor(2),R_nor(3),i,nor;
64//======================  Exercise 4.3 =============================
65ring R = 32003, (x,y,z,a,b,c,d), dp;
66ideal I = a-xy, b-y2, c-yz, d-y3;
67I = I, z2-x5-y5;
68ideal J = eliminate(I,y);
69ring S = 32003, (x,z,a,b,c,d), dp;
70ideal J = imap(R,J);
71attrib(J,"isSB",1);
72J;
73dim(J);
74//-> 2
75if(not(defined(normal))){ LIB "normal.lib"; }
76list nor = normal(J);
77def R_nor = nor[1]; setring R_nor; 
78norid; normap;
79
80
81kill R,S,nor,R_nor;
82//======================  Exercise 4.4 =============================
83ring R = 0, x(1..3), dp;
84matrix D[3][3] = x(1), x(2), x(3)^2-1,
85                 x(2), x(3), x(1)*x(2)+x(3)+1,
86                 x(3)^2-1, x(1)*x(2)+x(3)+1, 0;
87ideal I = det(D), x(1)*x(3)-x(2)^2;
88matrix Df = jacob(I);
89I = std(I);
90int c = nvars(R) - dim(I); c;
91//-> 2
92if (not(defined(equidimMax))){ LIB "primdec.lib"; }
93ideal I_max = equidimMax(I);
94size(reduce(I_max,I));
95//-> 0
96
97ideal J = minor(Df,c), I;
98J = groebner(J);
99sat(I,J);
100//-> [1]:
101//->    _[1]=1
102//-> [2]:
103//->    2
104
105ideal J_max = equidimMax(J);
106size(reduce(J_max,J));
107//-> 0
108matrix DJ = jacob(J);
109ideal SLoc = minor(DJ,c), J;
110groebner(SLoc);
111//-> _[1]=1
112
113//---------Alternatively: check sqrt(J)=J----------
114ideal sqrtJ = radical(J);
115size(reduce(sqrtJ,J));
116//-> 0
117
118
119kill R,c;
120//======================  Exercise 4.5 =============================
121proc maximaldegree (ideal I)
122"USAGE:  maximaldegree(I);   I=ideal
123RETURN: integer; the maximum degree of the given
124                 generators for I
125NOTE:   return value is -1 if I=0
126"
127{
128  if (size(I)>0)
129  {
130    int i,dd;
131    int d = deg(I[1]);
132    for (i=2; i<=size(I); i++)
133    {
134      dd = deg(I[i]);
135      if (dd>d) { d = dd; }
136    }
137    return(d);
138  }
139  else
140  {
141    return(-1);
142  }
143}
144
145proc zeros (poly f)
146"USAGE:   zeros(f);   f poly
147ASSUME:  the active ring is univariate
148RETURN:  ring
149NOTE:    In the output ring, f decomposes into linear factors. The
150         list of solutions of f=0 may be accessed by typing, for
151         instance,
152            def RRR=zeros(f); setring RRR; ZEROS;
153"
154{
155  if (nvars(basering)!=1)
156    { ERROR("The active ring is not univariate"); }
157  if (deg(f)<=0){ ERROR("f is a constant polynomial"); }
158  def R_aux = basering;
159  ideal facts_f = factorize(f,1);
160  int d = maximaldegree(facts_f);
161  int counter=size(facts_f);
162  int i;
163  while (d>1)
164  {
165    if (defined(primitive)==0){ LIB "primitiv.lib"; }
166    while (deg(facts_f[counter])<=1) { counter--; }
167    if (defined(minpoly)) { poly g = minpoly; }
168      else { poly g=0; }
169    poly h = facts_f[counter];
170    if (g==0)
171    {
172      ring R_aux2 = (char(basering),a), x, dp;
173      map psi = R_aux,a;
174      minpoly = number(psi(h));
175      ideal facts_old = imap(R_aux,facts_f);
176      ideal facts_f;
177      for (i=1; i<=size(facts_old); i++)
178      {
179        facts_f = facts_f, factorize(facts_old[i],1);
180      } 
181    }
182    else
183    {
184      ring R_aux1 = char(basering), (a,x), dp;
185      poly g = imap(R_aux,g);
186      poly h = imap(R_aux,h);
187      ideal facts_f = imap(R_aux,facts_f);
188      ideal I = g,h;
189      ideal Prim = primitive(I);
190      poly MP_new = Prim[1];
191      poly a_new = Prim[2];
192      poly x_new = Prim[3];
193      ring R_aux2 = (char(basering),a), x, dp;
194      map psi = R_aux1,0,a;
195      minpoly = number(psi(MP_new));
196      poly a_new = psi(a_new);
197      poly x_new = psi(x_new);
198      map phi = R_aux1,a_new,x;
199      ideal facts_old = phi(facts_f);
200      ideal facts_f;
201      for (i=1; i<=size(facts_old); i++)
202      {
203        facts_f = facts_f, factorize(facts_old[i],1);
204      }
205      kill R_aux1;
206    }
207    facts_f = simplify(facts_f,2);
208    kill R_aux;
209    def R_aux = R_aux2;
210    setring R_aux;
211    kill R_aux2;
212    d = maximaldegree(facts_f);
213  }
214  ideal ZEROS = -subst(facts_f,x,0);
215  export(ZEROS);
216  return(R_aux);
217}
218
219ring R = 167, x, dp;
220def RRR = zeros(x100-13);
221setring RRR;
222ZEROS;
223poly g = 1;
224for (int i=1; i<=size(ZEROS); i++) { g = g*(x-ZEROS[i]); }
225g;
226//-> x100-13
227
228
229kill i,RRR,R;
230//======================  Exercise 4.6 =============================
231ring R = 0, (x,y,z), dp;
232ideal MI_2 = maxideal(2);
233ideal MI_3 = maxideal(3);
234ideal MI_4 = maxideal(4);
235ring R_ext = 0, (a(1..6),b(1..6),c(1..6),x,y,z), dp;
236ideal MI_2 = imap(R,MI_2);
237matrix A[6][1] = a(1..6);
238matrix B[6][1] = b(1..6);
239matrix C[6][1] = c(1..6);
240poly F(0) = (matrix(MI_2)*A)[1,1];
241poly F(1) = (matrix(MI_2)*B)[1,1];
242poly F(2) = (matrix(MI_2)*C)[1,1];
243ideal G = x2*F(0), xy*F(0), xz*F(0), yz*F(0), x2*F(1),
244          xy*F(1), xz*F(1), y2*F(1), yz*F(1), x2*F(2),
245          xy*F(2), xz*F(2), y2*F(2), yz*F(2), z2*F(2);
246ideal MI_4 = imap(R,MI_4);
247matrix M = coeffs(G,MI_4,xyz);
248poly D_0 = det(M);   
249deg(D_0); size(D_0);
250//-> 15
251//-> 37490
252
253ideal G_ext = x2*F(1), x2*F(2), y2*F(2);
254ideal MI_ext = x2y2, x2z2, y2z2;
255matrix M_ext = coeffs(G_ext,MI_ext,xyz);
256poly D_0_ext = det(M_ext);
257poly Res = D_0/D_0_ext;
258deg(Res); size(Res);
259//-> 12
260//-> 21894
261
262
263// session continued.....
264//======================  Exercise 4.7 =============================
265ring S = 0,(x,y,z,a(1..10)), (dp(2),dp(11));
266ideal MI_2 = imap(R,MI_2);
267ideal MI_3 = imap(R,MI_3);
268matrix A[10][1] = a(10),a(9),a(7),a(4),a(8),a(6),a(3),a(5),a(2),
269                  a(1);
270poly F = (matrix(MI_3)*A)[1,1];
271ideal J = diff(F,x), diff(F,y), diff(F,z);
272if (not(defined(flatten))){ LIB "matrix.lib"; }
273map phi = R_ext, flatten(coeffs(diff(F,x),MI_2,xyz)),
274                 flatten(coeffs(diff(F,y),MI_2,xyz)),
275                 flatten(coeffs(diff(F,z),MI_2,xyz));
276poly D = phi(Res);             
277deg(D); size(D);
278//-> 12
279//-> 2040   
280
281tst_status(1);$
282
Note: See TracBrowser for help on using the repository browser.