Changeset 513673 in git

Timestamp:

Oct 15, 2013, 3:15:25 PM (11 years ago)

Author:

Andreas Steenpass <steenpass@…>

Branches:

(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')

Children:

3630419158783d9329a9e0d3d501208a60c71b10

Parents:

c21a6b6ea8d600f0dee5b13cbb0c042d2fe7e574

git-author:

Andreas Steenpass <steenpass@mathematik.uni-kl.de>2013-10-15 15:15:25+02:00

git-committer:

Andreas Steenpass <steenpass@mathematik.uni-kl.de>2013-10-15 15:19:17+02:00

Message:

chg: update realclassify.lib

The implementation of the splitting lemma has been changed.

Files:

• Singular/LIB/realclassify.lib (modified) (5 diffs)
• Tst/Short/realclassify_s.res.gz.uu (modified) (1 diff)
• Tst/Short/realclassify_s.stat (modified) (1 diff)

Singular/LIB/realclassify.lib

@@ -28 +28 @@
 ";
 
+LIB "linalg.lib";
 LIB "elim.lib";
 LIB "primdec.lib";
@@ -475 +476 @@
 USAGE:    realmorsesplit(f[, mu]); f poly, mu int
 RETURN:   a list consisting of the corank of f, the inertia index, an upper
-          bound for the determinacy, the residual form of f and
-          the transformation
+          bound for the determinacy, and the residual form of f
 NOTE:     The characteristic of the basering must be zero, the monomial order
           must be local, f must be contained in maxideal(2) and the Milnor
@@ -486 +486 @@
 EXAMPLE:  example morsesplit; shows an example"
 {
-  /* auxiliary variables */
-  int i, j;
+  int i;
 
   /* error check */
@@ -532 +531 @@
   }
 
-  /* preliminary stuff */
-  list S;
+  /* compute the determinacy */
   int k = determinacy(f, mu);
   f = jet(f, k);
-  def br = basering;
-  map Phi = br, maxideal(1);
-  map phi;
-  poly a, p, r;
-
-  /* treat the variables one by one */
+
+  /* get jet(f, 2) right */
+  matrix H = concat(jet(jacob(jacob(f)), 0)/2, unitmat(n));
+  H = sym_reduce(H);
+  intvec perm_zero;
+  intvec perm_neg;
+  intvec perm_pos;
+  int c;
+  int lambda;
   for(i = 1; i <= n; i++)
   {
-    if(jet(f, 2)/var(i) == 0)
-    {
-      S = insert(S, i);
-    }
-    else
-    {
-      f, a, p, r = rewriteformorsesplit(f, k, i);
-      if(jet(a, 0) == 0)
-      {
-        for(j = i+1; j <= n; j++)
-        {
-          if(jet(f, 2)/(var(i)*var(j)) != 0)
-          {
-            break;
-          }
-        }
-        phi = br, maxideal(1);
-        phi[j] = var(j)+var(i);
-        Phi = phi(Phi);
-        f = phi(f);
-      }
-      f, a, p, r = rewriteformorsesplit(f, k, i);
-      while(p != 0)
-      {
-        phi = br, maxideal(1);
-        phi[i] = var(i)-p/(2*jet(a, 0));
-        Phi = phi(Phi);
-        f = phi(f);
-        f, a, p, r = rewriteformorsesplit(f, k, i);
-      }
-    }
-  }
-
-  /* sort variables according to corank */
-  int cr = size(S);
-  phi = br, 0:n;
-  j = 1;
-  for(i = size(S); i > 0; i--)
-  {
-    phi[S[i]] = var(j);
-    j++;
-  }
+    if(H[i, i] == 0)
+    {
+      perm_zero = perm_zero, i;
+      c++;
+    }
+    if(H[i, i] < 0)
+    {
+      perm_neg = perm_neg, i;
+      lambda++;
+    }
+    if(H[i, i] > 0)
+    {
+      perm_pos = perm_pos, i;
+    }
+  }
+  intvec perm;
+  if(size(perm_zero) > 1)
+  {
+    perm = perm, perm_zero[2..size(perm_zero)];
+  }
+  if(size(perm_neg) > 1)
+  {
+    perm = perm, perm_neg[2..size(perm_neg)];
+  }
+  if(size(perm_pos) > 1)
+  {
+    perm = perm, perm_pos[2..size(perm_pos)];
+  }
+  perm = perm[2..size(perm)];
+  matrix T[n][n];
+  matrix D[1][n];
   for(i = 1; i <= n; i++)
   {
-    if(phi[i] == 0)
-    {
-      phi[i] = var(j);
-      j++;
-    }
-  }
-  Phi = phi(Phi);
+    T[1..n, i] = H[perm[i], (n+1)..(2*n)];
+    D[1, i] = H[perm[i], perm[i]];
+  }
+  map phi = basering, matrix(maxideal(1))*transpose(T);
   f = phi(f);
-
-  /* compute the inertia index lambda */
-  int lambda;
-  list negCoeff, posCoeff;
-  number ai;
-  poly f2 = jet(f, 2);
-  for(i = 1; i <= n; i++)
-  {
-    ai = number(f2/var(i)^2);
-    if(ai < 0)
-    {
-      lambda++;
-      negCoeff = insert(negCoeff, i);
-    }
-    if(ai > 0)
-    {
-      posCoeff = insert(posCoeff, i);
-    }
-  }
-
-  /* sort variables according to lambda */
-  phi = br, maxideal(1);
-  j = cr+1;
-  for(i = size(negCoeff); i > 0; i--)
-  {
-    phi[negCoeff[i]] = var(j);
-    j++;
-  }
-  for(i = size(posCoeff); i > 0; i--)
-  {
-    phi[posCoeff[i]] = var(j);
-    j++;
-  }
-  Phi = phi(Phi);
-  f = phi(f);
-
-  /* compute residual form */
-  phi = br, maxideal(1);
-  for(i = size(S)+1; i <= n; i++)
-  {
-    phi[i] = 0;
-  }
-  f = phi(f);
-
-  return(list(cr, lambda, k, f, Phi));
+  f = jet(f, k);
+
+  /* separate the variables */
+  phi = basering, maxideal(1);
+  map corank_part = basering, maxideal(1);
+  for(i = c+1; i <= n; i++)
+  {
+    corank_part[i] = 0;
+  }
+  poly h = f-jet(f, 2)-corank_part(f);
+  poly hi;
+  while(h != 0)
+  {
+    for(i = c+1; i <= n; i++)
+    {
+      hi = h/var(i);
+      phi[i] = var(i)-hi/(2*D[1, i]);
+      h = h-hi*var(i);
+    }
+    f = phi(f);
+    f = jet(f, k);
+    h = f-jet(f, 2)-corank_part(f);
+  }
+  poly g = cleardenom(f-jet(f, 2));
+
+  return(list(c, lambda, k, g));
 }
 example
@@ -649 +617 @@
   poly f = (x2+3y-2z)^2+xyz-(x-y3+x2z3)^3;
   realmorsesplit(f);
+}
+
+///////////////////////////////////////////////////////////////////////////////
+/*
+   symmetric Gauss algorithm
+
+   If A is not a square matrix, then the largest upper or left submatrix
+   is assumed to be symmetric.
+*/
+proc sym_reduce(matrix A)
+{
+  int r = nrows(A);
+  int c = ncols(A);
+  int n = r;
+  if(n > c)
+  {
+    n = c;
+  }
+  poly q;
+  int i, j;
+  for(i = 1; i <= n; i++)
+  {
+    for(j = i+1; j <= n; j++)
+    {
+      if(A[i, j] != 0)
+      {
+        while(A[i, i] == 0)
+        {
+          A[1..r, i] = A[1..r, i]+A[1..r, j];
+          A[i, 1..c] = A[i, 1..c]+A[j, 1..c];
+        }
+        q = A[i, j]/A[i, i];
+        A[1..r, j] = A[1..r, j]-q*A[1..r, i];
+        A[j, 1..c] = A[j, 1..c]-q*A[i, 1..c];
+      }
+    }
+  }
+  return(A);
 }
 

Tst/Short/realclassify_s.res.gz.uu

@@ -1 +1 @@
 begin 600 realclassify_s.res.gz
-M'XL("&H<3%(``W)E86QC;&%S<VEF>5]S+G)E<P!=4\%NHS`0O?LK1M$>0`$*
-M-DG39.&PFTND;2_M#:4132"R"B:RB=;FZ]?%B4W7!SSXS3S//!ZO;]O="P`D
-M.?S9_8)9+_JHH1^S#7J](3@'?7B@C/:>OT%?.^0Y\*ILCDTI!*W5042L^AN)
-MONQM&<G!QFEDR*<U_]VRF*0O(^"4G8%#!G$`G@Q4,/@!G(3+?\SATC4*:IWC
-M23PG*L2#_X[G4@VA)T-%YA(/Q'\GKF;UO6NO#F#&Z+&:Z;&*9+]&.F=;I/L0
-M%=B\#9H/*TV&"F).,-()8Y2B8F&B&!5+$^F\1Q.A-W6IH*M!Z$FN3<EIK]8W
-M]I>.MV4#M=[6<%_NJM\=+]FG0\:%T3-M6,>!7=N/BM_1%.U8Q7M:`F6G2MZ/
-M8[2M^HJWE)5'Y9A^9KI#*\>3DSR)HU&:MN.B$I=&?^IZ(@JV@L16"&*%""69
-MXP<B%;:*''1E)DV$]YE&U7P($YV$]3-9264PLL^4[2&9](-O%A#.`M\-D)")
-M`70#*G50FD,[BF6T&D=96G3BM$0[;20Y&Q*IG%L2[;"3$]$[:Y*%!5<3DJ=H
-4_$&^W'\57N)O?J!_US0V458#````
+M'XL("),_75(``W)E86QC;&%S<VEF>5]S+G)E<P!=4L%.PS`,O><KK(E#J[55
+MFW0#-N@!=ID$N\"M8E/9TBFB3:>D$TF_GM!N22&7.'[VB_WLM_?5>@,`208O
+MZR>8M+*-*O8Y6:*W"X(S,,X=XZSU_"7ZO2'+0-"BVE>%E*S4.QEQ^AW)MFAM
+M&LG`VFDTD(]S_OTR&X7/(Q",'T'`(\0!>"K00><'<)`N_C:#4U-I*$V,I_"4
+MZ!!W_A9/E>Y"3X6:3!7NB+\E+N?N;]5>&<"$LSV=F+;RY&.!3,PJ3S]"E./A
+MU1D^K`T9RLG@P<@$]%:*\ME@Q2B?#Y:)NQTL]*Y/%)H2I.GD7!6"M7IQ8=\T
+MHBXJ*,VU@.MQ7STWHN!?#ND/1J^LXHT`?JX_J;BB*5IS*EI6`.,'JJ[N&*UH
+M2T7->+'7CNGAT51HY;AWDB=QU$M3-T)2>:K,J,N1*-@*$ELAB!6"*!)BI;$E
+M2T;$^#)+Z6;Y=Y()&4U2D:E.'91F4/==#TWW-<TM.EJ9Q*Q,3W(<2)1V8T_,
+CJAR<&M[1D,PL>#<BN8_Z3?]=X[/T$G]Y@WX`&GCH!1\#````
 `
 end

Tst/Short/realclassify_s.stat

@@ -1 @@
-1 >> tst_memory_0 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:323576
-1 >> tst_memory_1 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:2322432
-1 >> tst_memory_2 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:2401344
-1 >> tst_timer_1 :: 1380719722:3160:3-1-6:ix86-Linux:mamawutz:16
+1 >> tst_memory_0 :: 1381842835:3160:3-1-6:ix86-Linux:mamawutz:323628
+1 >> tst_memory_1 :: 1381842835:3160:3-1-6:ix86-Linux:mamawutz:2322432
+1 >> tst_memory_2 :: 1381842835:3160:3-1-6:ix86-Linux:mamawutz:2401344
+1 >> tst_timer_1 :: 1381842835:3160:3-1-6:ix86-Linux:mamawutz:15