--- /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