Changeset 61549b in git for Singular/LIB/linalg.lib

Timestamp:

Nov 5, 2001, 10:17:05 AM (22 years ago)

Author:

Mathias Schulze <mschulze@…>

Branches:

(u'spielwiese', '4a9821a93ffdc22a6696668bd4f6b8c9de3e6c5f')

Children:

1f8111cd2acad6f5a631540b24272cf5aca44000

Parents:

5f403d5fffa49541688757a3ff48dc8d052361f2

Message:

*** empty log message ***


git-svn-id: file:///usr/local/Singular/svn/trunk@5663 2c84dea3-7e68-4137-9b89-c4e89433aadc

File:

• Singular/LIB/linalg.lib (modified) (13 diffs)

Singular/LIB/linalg.lib

@@ -1 +1 @@
 //GMG last modified: 04/25/2000
 //////////////////////////////////////////////////////////////////////////////
-version="$Id: linalg.lib,v 1.19 2001-08-25 14:56:41 mschulze Exp $";
+version="$Id: linalg.lib,v 1.20 2001-11-05 09:17:05 mschulze Exp $";
 category="Linear Algebra";
 info="
@@ -26 +26 @@
  pos_def(A,i);        test symmetric matrix for positive definiteness
  hessenberg(M);       transforms M to Hessenberg form
- spnormalize(a[,m]);  normalize spectrum
- eigenval(M);         eigenvalues of M with multiplicities
+ spnf(a[,m][,v|V]);   normalize spectrum
+ eigenvalues(M);      eigenvalues of M with multiplicities
  jordan(M);           Jordan data of constant square matrix M
  jordanbasis(M);      Jordan basis of constant square matrix M
@@ -1335 +1335 @@
 ///////////////////////////////////////////////////////////////////////////////
 
-static proc spappend(list l,number e,int m)
-{
-  if(size(l)==0)
-  {
-    l=list(ideal(e),intvec(m));
-  }
-  else
-  {
-    int n=ncols(l[1]);
-    l[1][n+1]=e;
-    l[2][n+1]=m;
-  }
-  return(l);
-}
-///////////////////////////////////////////////////////////////////////////////
-
-proc spnormalize(ideal e,list #)
-"USAGE:    spnormalize(a,w[,m]);
+proc spnf(ideal e,list #)
+"USAGE:   spnf(e[,m][,v|V]); ideal e, intvec m, module v, list V
 RETURN:
 @format
-list Sp: normalized spectrum (a,m)
+list Sp: normalized spectrum (e,m[,V])
   ideal Sp[1]: numbers in a in increasing order
   intvec Sp[2]:
-    int Sp[2][i]: multiplicity of number Sp[1][i] in a
+    int Sp[2][i]: multiplicity of number Sp[1][i] in e
+  list Spp[3]:
+    module Spp[3][i]: module associated to number Spp[1][i]
 @end format
-EXAMPLE:  example spnormalize; shows examples
+EXAMPLE: example spnf; shows examples
 "
 {
+  int n=ncols(e);
   intvec m;
+  module v;
+  list V;
   int i,j;
-  if(size(#)==0)
-  {
-    for(i=ncols(e);i>=1;i--)
+  while(i<size(#))
+  {
+    i++;
+    if(typeof(#[i])=="intvec")
+    {
+      m=#[i];
+    }
+    if(typeof(#[i])=="module")
+    {
+      v=#[i];
+      for(j=n;j>=1;j--)
+      {
+        V[j]=module(v[j]);
+      }
+    }
+    if(typeof(#[i])=="list")
+    {
+      V=#[i];
+    }
+  }
+  if(m==0)
+  {
+    for(i=n;i>=1;i--)
     {
       m[i]=1;
     }
   }
-  else
-  {
-    m=#[1];
-  }
-
-  list l;
-  number e0;
-  int m0;
-  for(i=ncols(e);i>=1;i--)
+
+  int k;
+  ideal e0;
+  intvec m0;
+  list V0;
+  number e1;
+  int m1;
+  for(i=n;i>=1;i--)
   {
     if(m[i]!=0)
@@ -1390 +1398 @@
           if(number(e[i])>number(e[j]))
           {
-            e0=number(e[i]);
+            e1=number(e[i]);
             e[i]=e[j];
-            e[j]=e0;
-            m0=m[i];
+            e[j]=e1;
+            m1=m[i];
             m[i]=m[j];
-            m[j]=m0;
+            m[j]=m1;
+            if(size(V)>0)
+            {
+              v=V[i];
+              V[i]=V[j];
+              V[j]=v;
+            }
           }
           if(number(e[i])==number(e[j]))
@@ -1401 +1415 @@
             m[i]=m[i]+m[j];
             m[j]=0;
+            if(size(V)>0)
+            {
+              V[i]=V[i]+V[j];
+            }
           }
         }
       }
-      l=spappend(l,number(e[i]),m[i]);
-    }
-  }
-
+      k++;
+      e0[k]=e[i];
+      m0[k]=m[i];
+      if(size(V)>0)
+      {
+        V0[k]=V[i];
+      }
+    }
+  }
+
+  if(size(V0)>0)
+  {
+    n=size(V0);
+    module U=std(V0[n]);
+    for(i=n-1;i>=1;i--)
+    {
+      V0[i]=simplify(reduce(V0[i],U),1);
+      if(i>=2)
+      {
+        U=std(U+V0[i]);
+      }
+    }
+  }
+
+  list l;
+  if(k>0)
+  {
+    l=e0,m0;
+    if(size(V0)>0)
+    {
+      l[3]=V0;
+    }
+  }
   return(l);
 }
@@ -1425 +1472 @@
     int l[2][i]: multiplicity of eigenvalue l[1][i]
 @end format
-EXAMPLE: example eigenval; shows examples
+EXAMPLE: example eigenvalues; shows examples
 "
 {
   M=jet(hessenberg(M),0);
   int n=ncols(M);
+  int k;
+  ideal e;
+  intvec m;
   number e0;
   intvec v;
-  list l,f;
-  int i,j,k;
+  list l;
+  int i,j;
   j=1;
   while(j<=n)
@@ -1454 +1504 @@
     if(size(v)==1)
     {
-      l=spappend(l,number(M[v,v]),1);
+      k++;
+      e[k]=M[v,v];
+      m[k]=1;
     }
     else
     {
-      f=factorize(det(submat(M,v,v)-var(1)));
-      for(i=size(f[1]);i>=1;i--)
+      l=factorize(det(submat(M,v,v)-var(1)));
+      for(i=size(l[1]);i>=1;i--)
       {
-        e0=number(jet(f[1][i]/var(1),0));
+        e0=number(jet(l[1][i]/var(1),0));
         if(e0!=0)
         {
-          l=spappend(l,number((e0*var(1)-f[1][i])/e0),f[2][i]);
+          k++;
+          e[k]=(e0*var(1)-l[1][i])/e0;
+          m[k]=l[2][i];
         }
       }
     }
   }
-  return(spnormalize(l[1],l[2]));
+  return(spnf(e,m));
 }
 example
@@ -1600 +1654 @@
   def e,m=#[1..2];
 
-  int i;
-  for(i=1;i<=ncols(e);i++)
+  for(int i=1;i<=ncols(e);i++)
   {
     if(deg(e[i]>0))
@@ -1610 +1663 @@
   }
 
-  int j,k,l;
+  int j,k,l,n;
   matrix N0,N1;
   module K0,K1;
@@ -1618 +1671 @@
   intvec w;
 
-  for(i=ncols(e);i>=1;i--)
+  for(i=1;i<=ncols(e);i++)
   {
     N0=M-e[i]*freemodule(ncols(M));
@@ -1648 +1701 @@
         for(j=l;j>=1;j--)
         {
-          U=module(u)+U;
-          w=2*j-l-1,w;
+          U=U+module(u);
+          n++;
+          w[n]=2*j-l-1;
           u=N0*u;
         }
@@ -1655 +1709 @@
     }
   }
-  w=w[1..size(w)-1];
+
   return(list(U,w));
 }
@@ -1696 +1750 @@
       for(i=s[k];i>=2;i--)
       {
-        J[j,j+1]=1;
+        J[j+1,j]=1;
         j++;
         J[j,j]=e[k];