--* From postmaster%watson.vnet.ibm.com@yktvmv.watson.ibm.com Thu Sep 23 11:15:33 1993 --* Received: from yktvmv2.watson.ibm.com by radical.watson.ibm.com (AIX 3.2/UCB 5.64/900524) --* id AA15478; Thu, 23 Sep 1993 11:15:33 -0400 --* X-External-Networks: yes --* Received: from watson.vnet.ibm.com by yktvmv.watson.ibm.com (IBM VM SMTP V2R3) --* with BSMTP id 1619; Thu, 23 Sep 93 11:19:49 EDT --* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0" --* id --* for asbugs@watson; Thu, 23 Sep 93 11:19:48 -0400 --* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0" --* id 2509; Thu, 23 Sep 1993 11:19:46 EDT --* Received: from matthew.watson.ibm.com by yktvmv.watson.ibm.com --* (IBM VM SMTP V2R3) with TCP; Thu, 23 Sep 93 11:19:44 EDT --* Received: by matthew.watson.ibm.com (AIX 3.2/UCB 5.64/4.03) --* id AA17472; Thu, 23 Sep 1993 11:21:20 -0400 --* Date: Thu, 23 Sep 1993 11:21:20 -0400 --* From: gianni@matthew.watson.ibm.com --* Message-Id: <9309231521.AA17472@matthew.watson.ibm.com> --* To: asbugs@watson.ibm.com --* Subject: Bug: Bad case 66 (line 1698 in file ../src/foam.c). [smatrixp.as][30.6] --@ Fixed by: SSD Tue Sep 27 17:44:28 EDT 1994 --@ Tested by: none --@ Summary: All three files (srow/elemrec1/smatrixp).as compile successfully under v0.37.0. #include "aslib.as" #library myrec "elemrec1.aso" import from myrec #library srow "srow.aso" import from srow Vector ==> Array SR ==> SparseRow FL SI ==> SingleInteger Col ==> Vector FL +++ sparse matrix domain defined as vector of sparse rows +++ having a BiModule structure for algebraic use. SparseMatrix(FL:Ring):SMC == SMD where SMC ==> BiModule(FL,FL) with transpose : % -> % ++ transpose a sparse matrix zero: (SingleInteger,SingleInteger) -> % ++Giving the zero matrix of dimension ri * ci zero? : % -> Boolean hlink: (%,%) -> % ++Horizontal linking of two matrices rm:=[sm1|sm2] to be used in ++Gauss Elimination Algorithm copy : % -> % ++make a copy of original matrix for distructive functions * : (SR,%) -> SR ++ sparserow by sparsematrix product * : (%,%) -> % ++ sparsematrix by sparsematrix product * : (%,Vector FL) -> Vector FL ++ sparse matrix by vector product setelt : (%,SingleInteger,SingleInteger,FL) -> FL ++Assign instruction : m(SingleInteger,j):=el-expr. elt:(%,SingleInteger,SingleInteger) ->FL ++ the (i,j) element of the matrix Exchange_! : (%,SingleInteger,SingleInteger) -> % ++ Exchange_!(m,i,j) exchange row i and row j addRows_! : (%,SingleInteger,SingleInteger ,FL) -> % ++ addRows_!(m,i,k,a) adds to the i-th row a*(k-th row) multiplyRow_! : (%,SingleInteger,FL) -> % ++ multiplyRow_!(m,i,a) multilpies by a the i-th row of m -- subMatrix : (%,List SingleInteger,List SingleInteger) -> % -- ++ subMatrix(m,lRows,lCols) extracts the submatrix with rows in lRows -- ++ and columns in lCols coerce : % -> Record(nrow:SingleInteger,ncol:SingleInteger,matrix:Vector SR) nRows : % -> SingleInteger nCols : % -> SingleInteger sMatrix : % -> Vector SR SMD ==> add Rep ==> Record(nrow:SingleInteger,ncol:SingleInteger,matrix:Vector SR) import from Integer import from FL import from SR import from List SR import from Vector SR import from List Boolean import from Vector FL import from Rep import from List FL import from SingleInteger 0 : % == per [0,0,new(0,0)] (m1:%)*(el:FL) : % == sm1:=rep m1 mm:=sm1.matrix nr:= sm1.nrow per[nr,sm1.ncol,[mm(i)*el for i in 1..nr]] (el:FL)*(m1:%) : % == sm1 := rep m1 mm:=sm1.matrix nr:=sm1.nrow per [nr,sm1.ncol,[el*mm(i) for i in 1..nr]] (m1:%) + (m2:%) : % == zero? m1 => m2 zero? m2 => m1 sm1:=rep m1 sm2:=rep m2 (nr:=sm1.nrow)~=sm2.nrow or (nc:=sm1.ncol)~=sm2.ncol => error "wrong dim" mm1:=sm1.matrix mm2:=sm2.matrix per [nr,nc,[(mm1(i)+ mm2(i)) for i in 1..nr]] +(sm1 : %): % == sm1 -(m1:% ) : % == sm1:=rep m1 mm:=sm1.matrix nr:=sm1.nrow per [nr,sm1.ncol,[(-mm(i)) for i in 1..nr]] (sm1: %) - (sm2 : %) : % == sm1 +(-sm2) (m1:%) = (m2:%) : Boolean == sm1:=rep m1 sm2 := rep m2 (nr:=sm1.nrow)~=sm2.nrow or (nc:=sm1.ncol)~=sm2.ncol => false mm1:=sm1.matrix mm2:=sm2.matrix reduce(/\,[(mm1(i) = mm2(i)) for i in 1..nr],false) transpose(m1 : %) : % == import from ElementRecord FL import from List ElementRecord FL sm1:=rep m1 mm:=sm1.matrix nr:=sm1.nrow nc:=sm1.ncol rm:Vector SR:=new(nc,0) for i in 1..nr repeat for c in convert mm(nr+1-i) repeat ind:=position(c) el:=element(c) rm(ind):=convert cons(create(el,nr+1-i),convert rm(ind)) per[nc,nr,rm] nRows(m:%) : SingleInteger == (rep m).nrow nCols(m:%) : SingleInteger == (rep m).ncol sMatrix(m:%) :Vector SR == (rep m).matrix zero?(sm1:%): Boolean == mm1:=rep sm1 (nr:=mm1.nrow) = 0 or mm1.ncol=0 => true mat:= mm1.matrix reduce(/\,[mat(i) = 0 for i in 1..nr],true) zero(ri:SingleInteger,ci:SingleInteger) : % == per [ri,ci,new(ri,0)] -- horizontal link of two matrices ris:=[sm1|sm2] hlink(sm1:%,sm2:% ):% == (nr1:=(rep sm1).nrow)~=(rep sm2).nrow => error "SparseMatrix hlink :incompatible matrices" mm1:=(rep sm1).matrix mm2:=(rep sm2).matrix nc1:=(rep sm1).ncol nc2:=(rep sm2).ncol per [nr1,nc1+(rep sm2).ncol, [(mm1(i)+translate(mm2(i),nc1)) for i in 1..nr1]] copy(m:%):% == x :=rep m mm:=x.matrix vm:Vector SR:=[ mm(i) for i in 1..x.nrow] per [x.nrow,x.ncol,vm] (sr:SR) * (sm1:%) : SR == sm2:=rep transpose copy sm1 mm:=sm2.matrix nr:=sm2.nrow sparse([sr*mm(i) for i in 1..nr]) (m1:%) * (m2:%) : % == sm1 := rep m1 sm2 := rep m2 sm1.ncol ~= sm2.nrow => error "SparseMatrix * :incompatible matrices" mm1:=sm1.matrix nr:=sm1.nrow nc:=sm2.ncol mm:Vector SR:=[(mm1(i)*m2) for i in 1..nr] per [nr,nc,mm] (m1:%) * (v:Vector FL) : Vector FL == sm1:=rep m1 mm:=sm1.matrix nr:=sm1.nrow [(mm(i)*v) for i in 1..nr] setelt(m:%,i:SingleInteger,j:SingleInteger,el:FL) : FL == mm:Vector SR:=(rep m).matrix mm(i):=set(mm(i),j,el) el elt(m1:%,i:SingleInteger,j:SingleInteger) : FL == sm1:= rep m1 (i< 1) or (i > sm1.nrow) => error "first index out of Range" (j < 1) or (j > sm1.ncol) => error "ssecond index out of range" search(sm1.matrix.i,j) coerce(sm1:%) : Record(nrow:SingleInteger,ncol:SingleInteger, matrix:Vector SR) == rep sm1 -- Exchanging the i-th and j-th rows Exchange_!(m1:%,i:SingleInteger,j:SingleInteger) : % == i = j => m1 sm1:=rep m1 mm:=sm1.matrix nr:=sm1.nrow nc:=sm1.ncol (i<=0 or i>nr)=> error "SparseMatrix Exchange :first Index Out of Range" (j<=0 or j>nr)=> error "SparseMatrix Exchange :second Index Out of Range" sr:SR:=mm(i) mm(i):=mm(j) mm(j):=sr per [nr,nc,mm] addRows_!(m1:%,k:SingleInteger,i:SingleInteger,a:FL) : % == sm1:=rep m1 nr:=sm1.nrow (k<=0 or k>nr) => error "SparseMatrix addRows :Out of Range Index" mm:=sm1.matrix nc:=sm1.ncol mm(k):=mm(k)+(a*mm(i)) per [nr,nc,mm] multiplyRow_!(m1:%,k:SingleInteger,a:FL) : % == sm1 := rep m1 mm:=sm1.matrix nr:=sm1.nrow nc:=sm1.ncol (k <= 0 or k>nr) =>error "SparseMatrix multiplyRow_! :Index Out of Range" mm(k):=a*mm(k) per [nr,nc,mm] -- subMatrix(sm1:%,lRows : List SingleInteger, lCols : List SingleInteger) : % == -- import List SingleInteger -- (minr<1 or minr > (nr:=sm1.nrow) or maxr < 1 or maxrnr)_ -- => error "row index out of range or incorrect" -- (minc<1 or minc>(nc:=sm1.ncol) or maxc<1 or maxcnc)_ -- => error "column index out of Range or incorrect" -- nr:=maxr-minr+1 -- nc:=maxc-minc+1 -- mm:=sm1.matrix -- mr:Vector SR:=new(nr,0) -- for i in minr .. Maxr repeat -- mr(i-minr+1):=translate(extract(mm(i),minc,maxc),-(minc)+1) -- [nr,nc,mr] --Functions temporarily commented -- coerce(mtr:M FL):$ == -- coercing a matrix in a sparse one -- nr:=nrows( mtr) -- nc:=ncols( mtr) -- mm:Vector SR:=new(nr,0) -- v:SR -- w:Vector FL -- for i in 1 .. nr repeat -- w:=row(mtr,i) -- v:=sparsing(w) -- mm(i):=v -- [nr,nc,mm] -- coerce(m:$):M FL == -- coerce a sparse matrix to a regular one -- mtr:M FL:=new(m.nrow,m.ncol,0) -- for i in 1 .. m.nrow repeat -- for c in list(m.matrix.i) repeat -- mtr(i,position(c)):=element(c) -- mtr -- i*sm1 == -- exponentation of SparseMatrix respect additive Group -- mm:=sm1.matrix -- nr:=sm1.nrow -- nc:=sm1.ncol -- [nr,nc,[(i*mm(j)) for j in 1 .. nr]] -- nni*sm1 == -- non negative exponantation of SparseMatrix respect additive Group -- mm:=sm1.matrix -- nr:=sm1.nrow -- nc:=sm1.ncol -- [nr,nc,[(nni*mm(i)) for i in 1..nr]] -- p*sm1 == -- positive exponantation of SparseMatrix respect additive Group -- mm:=sm1.matrix -- nr:=sm1.nrow -- nc:=sm1.ncol -- [nr,nc,[(p*mm(i)) for i in 1.. nr]] -- inverse(m1:%) : Partial(%) == -- sm1 := rep m1 -- (nr:= sm1.nrow)~= sm1.ncol => error "inverse: s.matrix not square" -- AB:%:=hlink(m1,scalarMatrix(sm1.nrow,1)) -- x:%:=rowEchelon AB -- elt(AB,nr,nr) = 0 =>"failed" -- subMatrix(AB,[1..nr],[nr+1..2*nr],(rep AB).ncol)