diff -r 000000000000 -r 0efde00f9229 python/left_prm.py --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/python/left_prm.py Fri May 27 14:24:59 2011 +0530 @@ -0,0 +1,92 @@ +#!/usr/bin/python + +import pylab as pl +from rowjoin import rowjoin +from polyfuncs import polsize +from t1calc import t1calc +from seshft import seshft +from colsplit import colsplit +from clcoef import clcoef + +def left_prm(N, degN, D, degD, job=1, gap=None): + """ """ + if gap == None: + gap = 1.0e8 + + F, degF = rowjoin(D, degD, N, degN) + Frows, Fbcols = polsize(F, degF) + Fcols = Fbcols * (degF+1) + T1 = pl.array([]) + pr = pl.array([]) + degT1, T1rows, shft = 0, 0, 0 + S=F + sel = pl.ones((Frows,1)) + T1bcols =1 + abar = pl.arange(Fbcols, Frows) + if T1rows: + T1 = pl.hstack((T1, pl.zeros((T1rows, Frows)))) + if T1.size == 0: + T1rows = 0 + while T1.size == 0 or T1rows < (Frows - Fbcols): + if T1.size == 0: + T1rows = 0 + Srows = Frows*T1bcols + T1, T1rows, sel, pr = t1calc(S, Srows, T1, T1rows, sel, pr, Frows, Fbcols, abar, gap) + T1rows, T1cols = pl.atleast_2d(T1).shape + if T1.size == 0: + T1rows = 0 + if T1rows < (Frows - Fbcols): + if T1.size != 0: + T1 = pl.hstack((T1, pl.zeros((T1rows,Frows)).squeeze())) + T1bcols += 1 + degT1 += 1 + shft += Fbcols + S = seshft(S,F,shft); + sel = pl.hstack((sel,sel[Srows-Frows:Srows+1])) + rowvec = pl.arange((T1bcols-1)*Frows+(Fbcols), (T1bcols * Frows)) + abar = pl.hstack((abar, rowvec)) + + B, degB, A, degA = colsplit(T1, degT1, Fbcols, Frows-Fbcols) + B, degB = clcoef(B, degB) + B = -B; + A, degA = clcoef(A, degA) + if job == 1: + return B, degB, A, degA + + if job == 2: + S = S[sel!=0,:] + redSrows, Scols = S.shape + C = pl.hstack(pl.eye(Fbcols), pl.zeros(Fbcols,Scols-Fbcols)) + T2 = pl.dot(C, pl.inv(S)) + T2 = makezero(T2,gap); + T2 = move(T2, pl.find(sel), Srows); + X, degX, Y, degY = colsplit(T2, degT1, Fbcols, Frows - Fbcols); + X, degX = clcoef(X, degX); + Y, degY = clcoef(Y, degY); + elif job == 3: + Y = S + degY = sel + X = degT1 + degX = Fbcols + + else: + print 'Message from left_prm:no legal job number specified' + exit() + + return B, degB, A, degA, Y, degY, X, degX + + +if __name__== "__main__": + D = pl.array([[1, 0, 0, 0, 0, 0], + [0, 1, 0, 1, 0, 0], + [0, 0, 1, 1, 1, 0]]) + + N = pl.array([[1, 0, 0], + [0, 1, 0], + [0, 0, 1]]) + + dD, dN = 1, 0 + + result = left_prm(N, dN, D, dD) + print result + pass