--* From postmaster%watson.vnet.ibm.com@yktvmv.watson.ibm.com Wed Sep 22 23:41:07 1993 --* Received: from yktvmv2.watson.ibm.com by radical.watson.ibm.com (AIX 3.2/UCB 5.64/900524) --* id AA17160; Wed, 22 Sep 1993 23:41:07 -0400 --* X-External-Networks: yes --* Received: from watson.vnet.ibm.com by yktvmv.watson.ibm.com (IBM VM SMTP V2R3) --* with BSMTP id 2293; Wed, 22 Sep 93 23:45:22 EDT --* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0" --* id --* for asbugs@watson; Wed, 22 Sep 93 23:45:21 -0400 --* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0" --* id 3565; Wed, 22 Sep 1993 23:45:20 EDT --* Received: from matthew.watson.ibm.com by yktvmv.watson.ibm.com --* (IBM VM SMTP V2R3) with TCP; Wed, 22 Sep 93 23:45:19 EDT --* Received: by matthew.watson.ibm.com (AIX 3.2/UCB 5.64/4.03) --* id AA15146; Wed, 22 Sep 1993 23:46:52 -0400 --* Date: Wed, 22 Sep 1993 23:46:52 -0400 --* From: gianni@matthew.watson.ibm.com --* Message-Id: <9309230346.AA15146@matthew.watson.ibm.com> --* To: asbugs@watson.ibm.com --* Subject: compiler dies running out of paging space [baco.as][30.5 (current)] --@ Fixed by: SSD Tue Sep 27 17:27:31 EDT 1994 --@ Tested by: none --@ Summary: In v0.37.0, a version of this file compiles successfully in about 60sec and 18.4Mb total. #if 0 PI: I don't get compiler out of paging space, but tinfer allocates a huge memory space. The msg is: *** Starting "tinfer" phase... Time 94.750 s Store 26619904 B : 25135056 B alloc - 1128600 B free - 0 gc = 25513228 "bug425.as", line 49: ERD ==> add -- ElementRecord Definition of its Category-functions ............^ [L49 C13] (Error) (After Macro Expansion) No suitable interpretation for the exp ression `add ...' The following exports were not provided: Missing sample: % (follow the expanded expression on several pages)_ #endif #include "aslib.as" +++ Record of two fields:Element and Position used as SparseRow Element entry ElementRecord (FL:Ring): ERC == ERD where ERC ==> BasicType with -- Element-Record-Category element : % -> FL ++ Selecting the Element entry in ElementRecord position : % -> SingleInteger ++ Selecting the position entry in ElementRecord create :(FL,SingleInteger) -> % ++ Create an ElementRecord from its fields elements -- coerce :% -> OutputForm -- ++ Coerce to an External Format = :(%,%) -> Boolean ++ Equality test ERD ==> add -- ElementRecord Definition of its Category-functions -- Internal representation of ElementRecord Rep ==> Record(el:FL,pos:SingleInteger) import from FL import from SingleInteger import from Rep -- Definition of Element-Selector element(s:%):FL == rep(s).el -- Definition of Position-Selector position(s:%):SingleInteger == rep(s).pos create(elle:FL,pp:SingleInteger):% == per [elle,pp] -- Definition of coerce functions -- coerce(s:%):OutputForm == -- out:OutputForm:= -- hconcat("("::OutputForm,hconcat(s.pos::OutputForm,")"::OutputForm)) -- out:=hconcat(s.el::OutputForm,out) -- out (s1:%) = (s2:%):Boolean == element(s1) = element(s2) and position(s1) = position(s2) apply(p: OutPort, x: %) : OutPort == import from String p("(pos: ")(position x)(" el: ")(element x)(")") Vector ==> Array BiModule(R:Ring,S:Ring) : Category == AbelianGroup with * : (R,%) -> % * : (%,S) -> % +++ Domain of Sparse Rows for the domain of Sparse Matrices SparseRow(FL:Ring): SRC == SparseRowDefinition where SRC ==> BiModule(FL,FL) with convert : List ElementRecord FL -> % ++ convert a list of records in a sparse row convert : % -> List ElementRecord FL ++ convert a sparse row in a list of Records component:(%,SingleInteger) -> ElementRecord FL ++ returning the i-th non zero component of the sparse row translate:(%,SingleInteger) -> % ++translate by m positions the row elements extract:(%,lRows : List SingleInteger)->% ++ extract the row elements with index in lRows nonZero? : (%,SingleInteger) -> Boolean ++ nonZero?(row,i) = true if the i-th element of row is not zero mesure : % -> SingleInteger ++ the number of non zero element *:(%,Vector FL) -> FL ++ scalar product sparse row by vector *:(Vector FL,%) -> FL ++ scalar product vector by sparse row *:(%,%) -> FL ++ scalar product of sparse rows sparse : Vector FL -> % ++ transforms a vector in a sparse row vectorise : (%,SingleInteger) -> Vector FL ++ transforms a sparse row in vector search:(%,SingleInteger) -> FL ++ search for the value of thye k-th entry set:(%,SingleInteger,FL) ->% ++ set(row,k,a) set the k-th entry equal to a -- if FL has FIELD then -- "/":(%,FL) -> % SparseRowDefinition ==> add import from FL import from List FL import from ElementRecord FL import from Integer import from Vector FL import from List SingleInteger import from SingleInteger Rep ==> List ElementRecord FL -- internal representation import from Rep (sr:%) * (r:FL) : % == per [create(element(e)*r,position(e)) for e in rep sr] (r:FL) * (sr: %):% == per [create(r*element(e),position(e)) for e in rep sr] -(l:%):% == per [create(-element(ll),position(ll)) for ll in rep l] ((sr:%) - (rs:%)):% == sr + (-rs) ((sr:%) = (rs:%)) : Boolean == lsr:List ElementRecord FL:= rep sr lrs:List ElementRecord FL := rep rs empty?(lsr) => empty?(lrs) empty?(lrs) => false import from List Boolean reduce(/\ ,[(s1 = s2 ) for s1 in lsr for s2 in lrs],true) ((sr:%) ~= (rs:%)) : Boolean == not(sr = rs) 0 == per [] (+ (sr:%)):% == sr ((sr:%) + (rs:%)): % == sr = 0 => rs rs = 0 => sr lres:List ElementRecord FL :=[] lsr:List ElementRecord FL:=rep sr lrs :List ElementRecord FL:=rep rs local xrs:ElementRecord FL local xsr:ElementRecord FL while (~empty?(lsr) or ~empty?(lrs)) repeat xsr:=first lsr xrs:=first lrs (pxsr:=position xsr) < (pxrs:=position xrs) => lres:=cons(xsr,lres) lsr:=rest lsr pxsr > pxrs => lres:=cons(xrs,lres) lrs:=rest lrs lres:=cons(create(element xrs + element xsr,pxrs),lres) lrs:=rest lrs lsr :=rest lsr empty?(lrs) => per concat(reverse lres,lsr) per concat(reverse lrs,lrs) -- if FL has Field then -- sr:% / el:FL == -- ++ division operation only if FL is a field -- el = 0 => error "you are dividing by zero" -- --[create((element(rc)::FL/el::FL),position(rc)) for rc in sr]@% -- inv(el)*sr convert(ls : List ElementRecord FL) : % == per ls convert(s:%) : List ElementRecord FL == rep s -- selecting the i-th element of the row component(row:%,i:SingleInteger): ElementRecord FL == for s in rep(row) repeat if position s = i then return s create(0,0) -- setting the j-th element set(row :%,j:SingleInteger,el:FL):% == row = 0 => el = 0 => row per [create(el,j)] lsr:List ElementRecord FL :=[] lrow:List ElementRecord FL :=rep row local xs : ElementRecord FL -- inizialize while lrow~=[] \/ position(xs:=first lrow)< j repeat lsr:=cons(xs,lsr) lrow:=rest lrow if position xs =j then lrow:=rest lrow if el ~= 0 then lrow:=cons(create(el,j),lrow) per concat (reverse lsr,lrow) translate(sr:%,j:SingleInteger):% == per [create(element(c),j+position(c)) for c in rep sr] -- extract the row elements with index in lInd extract(sr:%,lInd : List SingleInteger):% == per [esr for esr in rep sr | member?(position esr,lInd)] -- check if the i-th element is non zero nonZero?(row:%,i:SingleInteger) : Boolean == for s in rep(row) repeat if position s = i then return true false -- number of non zero elements in sparse row mesure(srow:%):SingleInteger == #(rep srow) -- scalar product between two sparse rows (row1: %) * (row2:%) : FL == row1=0 or row2=0 => 0 srow1:=rep row1 reduce(+, [element(s)*element(component(row2,position s)) for s in srow1], 0) -- scalar product of sparse row and a regular vector (srow: %) * (vect: Vector FL) : FL == reduce(+,[element(s)*vect(retract position s) for s in rep srow],0) -- scalar product of a vector and a sparse row (vect : Vector FL) * (srow : %) : FL == reduce(+,[vect(retract position s)*element(s) for s in rep srow],0) -- transforms a vector in a sparse row sparse(v:Vector FL):% == lrow :List ElementRecord FL := empty() for i in 1..#v repeat (el:=v(i)) ~= 0 => lrow:=cons(create(el,i),lrow) per reverse lrow -- transform a sparse row in a vector of dimension clm vectorise(row:%,clm:SingleInteger):Vector FL == v:Vector FL:=new(clm,0) for rc in rep row repeat v(retract position rc):=element rc v -- search for the value of thye k-th entry search(sr:%,k:SingleInteger) : FL == for c in rep sr repeat if (pc:=position c) = k then return element c if pc >k then return 0 0 -- integer times sparserow -- i * sr == [create((i*element(e)),position(e)) for e in sr]@% -- non negative integer times sparserow -- n * sr == [create((n*element(e)),position(e)) for e in sr]@% -- positive integer times sparserow -- p * sr == [create((p*element(e)),position(e)) for e in sr]@% 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 from 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) SparseMatrixOperations(FL:Field) : SMC == SMD where MS ==> SparseMatrix FL SMC ==> with rowEchelon : MS -> MS ++ Gauss-Jordan Form of the Sparse Matrix rank : MS -> SingleInteger ++ rank(m) is the rank of m nullity : MS -> SingleInteger ++ nullity(m) = number of (infinity of ) solutions of the homogeneous ++ linear system mX = 0 nullSpace : MS -> List Col ++ nullSpace(m) = a basis of the solutions of the homogeneous linear ++ system mX = 0 determinant : MS -> FL ++ the determinant function for square sparse matrix ++ obtained by Gaussian Elimination as ++ det(E1*..*En)*det(A)=1 => det(A)=inv( det(E1)*..*det(En) ) ++ i.e. the product of the pivot elements of the Gaussian Elimination -- inverse : MS -> Partial(MS) -- ++ invert the SparseMatrix via Gauss-Jordan Elimination SMD ==> add 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 SparseMatrix FL import from List FL import from SingleInteger rowAllzero? : (MS,SingleInteger) -> Boolean colAllzero? : (MS,SingleInteger) -> Boolean rowAllzero?(sm1:MS,i:SingleInteger): Boolean == (sMatrix sm1)(i) = 0 colAllzero?(m1:MS,i:SingleInteger) : Boolean == mm:=sMatrix m1 nr:= nRows m1 reduce(/\,[(search(mm(j),i))=0 for j in 1..nr],false) rowEchelon(sm1 : MS) : MS == m:=copy sm1 nr:=nRows m nc:=nCols m i:SingleInteger:=1 local c:ElementRecord FL n:SingleInteger:=0 for j in 1..nc repeat if i > nr then return m smt: SR :=(sMatrix m)(i) if smt = 0 then n:=0 else c:=component(smt,1) n:=position(c) if n ~= 0 then if i ~= n then Exchange_!(m,i,j) b:=element(c) multiplyRow_!(m,i,inv(b)) for k in 1..nr repeat el:=search((sMatrix m)(k),j) if k ~= i and el ~= 0 then addRows_!(m,k,i,-el) i:=i+1 m rank(sm1:MS) : SingleInteger == zero? sm1 => 0 sy:= (rk:=nRows sm1) >(rh:= nCols sm1) => rk:=rh transpose sm1 copy sm1 sy:=rowEchelon sy i:=nRows sy while rk >0 and rowAllzero?(sy,i) repeat i:=i-1 rk:=(rk-1) rk nullity(m : MS ) : SingleInteger == nCols m - rank m nullSpace(y : MS) : List Col == x:=rowEchelon y nrw:= nRows x ncl:= nCols x mm := sMatrix x basis:List Vector FL:= [] rk:=nrw rw:=nrw while rk >0 and rowAllzero?(x,rw) repeat rk:=(rk-1) rw:=(rw-1) ncl<=nrw and rk = ncl =>[new(ncl,0)] local j : SingleInteger rk = 0 => for j in 1..ncl repeat w: Vector FL:=new(ncl,0) w(j):= 1 basis:=cons(w,basis) basis v:Vector SingleInteger:=new(ncl,0) local c:ElementRecord FL for i in 1..rk repeat c:=component(mm(i),1) j:=position(c) v(j):= i j:=ncl l:=ncl while j >= 1 repeat w:Vector FL:=new(ncl,0) if v(j) = 0 then colAllzero?(x,j) => w(l) := 1 basis:=cons(w,basis) for k in 1 .. (j-1) for ll in 1 .. (l -1) repeat if v(k) ~=0 then w(ll) :=-search(mm(v(k)),j) w(l) := 1 basis:=cons(w,basis) j:=j-1 l:=l-1 basis determinant(sm1 : MS ) : FL == (nr:=nRows sm1)~=(nCols sm1) => error "SparseMatrix determinant:Matrix is not square" -- Gaussian Elimination nr = 1 => search((sMatrix sm1)(1),1) mx:= copy sm1 ans:FL:=1 local c : ElementRecord FL for i in 1..(nr-1) repeat -- Gaussian Elimination Loop j:=i while (position(c:=component((sMatrix mx)(j),1))~=i and j <= nr) repeat j:=j+1 -- searching for a non zero pivot if j > nr then return 0 if j ~= i then Exchange_!(mx,i,j) --exchange rows to have a non zero ans:=-ans --element as pivot for Gaussian Elim. el:=element(c) ans:=el*ans -- determinant induced by multiplyRow multiplyRow_!(mx,i,inv(el)) for k in (i+1)..(nr) repeat c:=component((sMatrix mx)(k),1) p:=position(c) el:=element(c) p ~= i => iterate -- the (k,i) element of matrix x is zero addRows_!(mx,k,i,-inv(el)) -- adds to k-th row the el mult. of j-th --ans:=1*ans the determinant of addRows Elementary Matrix ans:=ans*search((sMatrix mx)(nr),nr)