// $Id: random.lib,v 1.2 1997-04-28 19:27:23 obachman Exp $
//system("random",787422842);
//(GMG/BM, last modified 22.06.96)
///////////////////////////////////////////////////////////////////////////////

LIBRARY:  random.lib    PROCEDURES OF RANDOM MATRIX AND POLY OPERATIONS

 genericid(id[,p,b]);   generic sparse linear combinations of generators of id
 randomid(id,[k,b]);    random linear combinations of generators of id
 randommat(n,m[,id,b]); nxm matrix of random linear combinations of id
 sparseid(k,u[,o,p,b]); ideal of k random sparse poly's of degree d [u<=d<=o]
 sparsemat(n,m[,p,b]);  nxm sparse integer matrix with random coefficients
 sparsepoly(u[,o,p,b]); random sparse polynomial with terms of degree in [u,o]
 sparsetriag(n,m[..]);  nxm sparse lower-triag intmat with random coefficients
           (parameters in square brackets [] are optional)

LIB "inout.lib";
LIB "general.lib";
////////////////////////////////////////////////////////////////////////////////

proc genericid (id, list #)
USAGE:   genericid(id,[,p,b]);  id ideal/module, k,p,b integers
RETURN:  system of generators of id which are generic, sparse, triagonal linear
         combinations of given generators with coefficients in [1,b] and
         sparsety p percent, bigger p being sparser (default: p=75, b=30000)
NOTE:    For performance reasons try small bound b in characteristic 0
EXAMPLE: example genericid; shows an example
{
//----------------------------- set defaults ----------------------------------
   if( size(#)>=2 ) { int p=#[1]; int b=#[2]; }
   if( size(#)==1 ) { int p=#[1]; int b=30000}
   if( size(#)==0 ) { int p=75; int b=30000; }
//---------------- use sparsetriag for creation of genericid ------------------
   def i = simplify(id,10);
   i = i*sparsetriag(ncols(i),ncols(i),p,b);
   return(i);
}
example
{ "EXAMPLE:"; echo = 2;
   ring r=0,(t,x,y,z),ds;
   ideal i= x3+y4,z4+yx,t+x+y+z;
   genericid(i,0,10);
   module m=[x,0,0,0],[0,y2,0,0],[0,0,z3,0],[0,0,0,t4];
   print(genericid(m));
}
////////////////////////////////////////////////////////////////////////////////

proc randomid (id, list #)
USAGE:   randomid(id,[k,b]);  id ideal/module, b,k integers
RETURN:  ideal/module having k generators which are random linear combinations
         of generators of id with coefficients in the interval [-b,b]
         (default: b=30000, k=size(id))
NOTE:    For performance reasons try small bound b in characteristic 0
EXAMPLE: example randomid; shows an example
{
//----------------------------- set defaults ----------------------------------
   if( size(#)>=2 ) { int k=#[1]; int b=#[2]; }
   if( size(#)==1 ) { int k=#[1]; int b=30000; }
   if( size(#)==0 ) { int k=size(id); int b=30000; }
//--------------------------- create randomid ---------------------------------
   def i = id;
   i = matrix(id)*random(b,ncols(id),k);
   return(i);
}
example
{ "EXAMPLE:"; echo = 2;
   ring r=0,(x,y,z),dp;
   randomid(maxideal(2),2,9);
   module m=[x,0,1],[0,y2,0],[y,0,z3];
   show(randomid(m));
}
////////////////////////////////////////////////////////////////////////////////

proc randommat (int n, int m, list #)
USAGE:   randommat(n,m[,id,b]);  n,m,b integers, id ideal
RETURN:  nxm matrix, entries are random linear combinations of elements
         of id and coefficients in [-b,b]
         [default: (id,b) = (maxideal(1),30000)]
NOTE:    For performance reasons try small bound b in char 0
EXAMPLE:  example randommat; shows an example
{
//----------------------------- set defaults ----------------------------------
   if( size(#)>=2 ) { ideal id=#[1]; int b=#[2]; }
   if( size(#)==1 ) { ideal id=#[1]; int b=30000; }
   if( size(#)==0 ) { ideal id=maxideal(1); int b=30000; }
//--------------------------- create randommat --------------------------------
   id=simplify(id,2);
   int g=ncols(id);
   matrix rand[n][m]; matrix ra[1][m];
   for (int k=1; k<=n; k=k+1)
   {
      ra = id*random(b,g,m);
      rand[k,1..m]=ra[1,1..m];
   }
   return(rand);
}
example
{ "EXAMPLE:"; echo = 2;
   ring r=0,(x,y,z),dp;
   matrix A=randommat(3,3,maxideal(2),9);
   print(A);
   A=randommat(2,3);
   print(A);
}
///////////////////////////////////////////////////////////////////////////////

proc sparseid (int k, int u, list #)
USAGE:   sparseid(k,u[,o,p,b]);  k,u,o,p,b integers
RETURN:  ideal having k generators in each degree d, u<=d<=o, p percent of
         terms in degree d are 0, the remaining have random coefficients
         in the interval [1,b], (default: o=u=d, p=75, b=30000)
EXAMPLE: example sparseid; shows an example
{
//----------------------------- set defaults ----------------------------------
   if( size(#)>=3 ) { int o=#[1]; int p=#[2]; int b=#[3]; }
   if( size(#)==2 ) { int o=#[1]; int p=#[2]; int b=30000; }
   if( size(#)==1 ) { int o=#[1]; int p=75; int b=30000; }
   if( size(#)==0 ) { int o=u; int p=75; int b=30000; }
//------------------ use sparsemat for creation of sparseid -------------------
   int ii; ideal i; intmat m;
   for ( ii=u; ii<=o; ii=ii+1)
   {
       m = sparsemat(size(maxideal(ii)),k,p,b);
       i = i+ideal(matrix(maxideal(ii))*m);
   }
   return(i);
}
example
{ "EXAMPLE:"; echo = 2;
   ring r = 0,(a,b,c,d),ds;
   sparseid(3,4);"";
   sparseid(2,2,5,90,9);
}
///////////////////////////////////////////////////////////////////////////////

proc sparsemat (int n, int m, list #)
USAGE:   sparsemat(n,m[,p,b]);  n,m,p,b integers
RETURN:  nxm integer matrix, p percent of the entries are 0, the remaining
         are random coefficients >=1 and <= b; [defaults: (p,b) = (75,1)]
EXAMPLE: example sparsemat; shows an example
{
   int r,h,ii;
   int t = n*m;
   intmat v[1][t];
//----------------------------- set defaults ----------------------------------
   if( size(#)>=2 ) { int p=#[1]; int b=#[2]; }
   if( size(#)==1 ) { int p=#[1]; int b=1; }
   if( size(#)==0 ) { int p=75; int b=1; }
//------------------------- check trivial cases ------------------------------
   if( p<0 ) { p = 0; }
   if(p>100) { p=100; }
//--------------- this is faster for not very sparse matrices ----------------
   if( p<40 )
   {
      for( ii=1; ii<=t; ii=ii+1 )
      { r=( random(1,100)>p ); v[1,ii]=r*random(1,b); h=h+r; }
   }
  int bb = t*(100-p);
  if( 100*h > bb )
  {
     while( 100*h > bb )
     { r=random(1,t); h=h-( v[1,r]>0 ); v[1,r]=0; }
  }
  else
  {
//------------------- this is faster for sparse matrices ---------------------
     while ( 100*h < bb )
     { r=random(1,t); h=h+(v[1,r]==0); v[1,r]=random(1,b); }
  }
  intmat M[n][m] = v[1,1..t];
  return(M);
}
example
{ "EXAMPLE:"; echo = 2;
   sparsemat(5,5);"";
   sparsemat(5,5,95);"";
   sparsemat(5,5,5);"";
   sparsemat(5,5,50,100);
}
///////////////////////////////////////////////////////////////////////////////

proc sparsepoly (int u, list #)
USAGE:   sparsepoly(u[,o,p,b]);  u,o,p,b integers
RETURN:  poly having only terms in degree d, u<=d<=o, p percent of the terms
         in degree d are 0, the remaining have random coefficients in [1,b),
         (defaults: o=u=d, p=75, b=30000)
EXAMPLE:  example sparsepoly; shows an example
{
//----------------------------- set defaults ----------------------------------
   if( size(#)>=3 ) { int o=#[1]; int p=#[2]; int b=#[3]; }
   if( size(#)==2 ) { int o=#[1]; int p=#[2]; int b=30000; }
   if( size(#)==1 ) { int o=#[1]; int p=75; int b=30000; }
   if( size(#)==0 ) { int o=u; int p=75; int b=30000; }
   int ii; poly f;
//----------------- use sparseid for creation of sparsepoly -------------------
   for( ii=u; ii<=o; ii=ii+1 ) { f=f+sparseid(1,ii,ii,p,b)[1]; }
   return(f);
}
example
{ "EXAMPLE:"; echo = 2;
   ring r=0,(x,y,z),dp;
   sparsepoly(5);"";
   sparsepoly(3,5,90,9);
}
///////////////////////////////////////////////////////////////////////////////

proc sparsetriag (int n, int m, list #)
USAGE:   sparsetriag(n,m[,p,b]);  n,m,p,b integers
RETURN:  nxm lower triagonal integer matrix, diagonal entries equal to 1, about
         p percent of lower diagonal entries are 0, the remaining are random
         integers >=1 and <= b; [defaults: (p,b) = (75,1)]
EXAMPLE: example sparsetriag; shows an example
{
   int ii,min,l,r; intmat M[n][m];
   int t=(n*(n-1))/2;
//----------------------------- set defaults ----------------------------------
   if( size(#)>=2 ) { int p=#[1]; int b=#[2]; }
   if( size(#)==1 ) { int p=#[1]; int b=1; }
   if( size(#)==0 ) { int p=75; int b=1; }
//---------------- use sparsemat for creation of sparsetriag ------------------
   intmat v[1][t]=sparsemat(1,t,p,b);
   if( n<=m ) { min=n-1; M[n,n]=1; }
   else { min=m; }
   for( ii=1; ii<=min; ii=ii+1 )
   {
      l=r+1; r=r+n-ii;
      M[ii..n,ii]=1,v[1,l..r];
   }
   return(M);
}
example
{ "EXAMPLE:"; echo = 2;
   sparsetriag(5,7);"";
   sparsetriag(7,5,90);"";
   sparsetriag(5,5,0);"";
   sparsetriag(5,5,50,100);
}
///////////////////////////////////////////////////////////////////////////////