1 #!/usr/bin/python

     2 

     3 import pylab as pl

     4 from rowjoin import rowjoin

     5 from polyfuncs import polsize

     6 from t1calc import t1calc

     7 from seshft import seshft

     8 from colsplit import colsplit

     9 from clcoef import clcoef

    10 

    11 def left_prm(N, degN, D, degD, job=1, gap=None):

    12     """ """

    13     if gap == None:

    14         gap = 1.0e8

    15     

    16     F, degF = rowjoin(D, degD, N, degN)

    17     Frows, Fbcols = polsize(F, degF)

    18     Fcols = Fbcols * (degF+1)

    19     T1 = pl.array([])

    20     pr = pl.array([]) 

    21     degT1, T1rows, shft = 0, 0, 0

    22     S=F

    23     sel = pl.ones((Frows,1))

    24     T1bcols =1

    25     abar = pl.arange(Fbcols, Frows)

    26     if T1rows:

    27         T1 = pl.hstack((T1, pl.zeros((T1rows, Frows))))

    28     if T1.size == 0:

    29         T1rows = 0

    30     while T1.size == 0 or T1rows < (Frows - Fbcols):

    31         if T1.size == 0:

    32             T1rows = 0

    33         Srows = Frows*T1bcols     

    34         T1, T1rows, sel, pr = t1calc(S, Srows, T1, T1rows, sel, pr, Frows, Fbcols, abar, gap)

    35         T1rows, T1cols = pl.atleast_2d(T1).shape

    36         if T1.size == 0:

    37             T1rows = 0

    38         if T1rows < (Frows - Fbcols):

    39             if T1.size != 0:

    40                 T1 = pl.hstack((T1, pl.zeros((T1rows,Frows)).squeeze()))

    41             T1bcols += 1

    42             degT1 += 1

    43             shft += Fbcols

    44             S = seshft(S,F,shft);

    45             sel = pl.hstack((sel,sel[Srows-Frows:Srows+1]))

    46             rowvec = pl.arange((T1bcols-1)*Frows+(Fbcols), (T1bcols * Frows))

    47             abar = pl.hstack((abar, rowvec))

    48 

    49     B, degB, A, degA = colsplit(T1, degT1, Fbcols, Frows-Fbcols)

    50     B, degB = clcoef(B, degB)

    51     B = -B;

    52     A, degA = clcoef(A, degA)

    53     if job == 1:

    54         return B, degB, A, degA

    55 

    56     if job == 2:

    57         S = S[sel!=0,:]

    58         redSrows, Scols = S.shape

    59         C = pl.hstack(pl.eye(Fbcols), pl.zeros(Fbcols,Scols-Fbcols))

    60         T2 = pl.dot(C, pl.inv(S))

    61         T2 = makezero(T2,gap);

    62         T2 = move(T2, pl.find(sel), Srows);

    63         X, degX, Y, degY = colsplit(T2, degT1, Fbcols, Frows - Fbcols);

    64         X, degX = clcoef(X, degX);

    65         Y, degY = clcoef(Y, degY);

    66     elif job == 3:

    67         Y = S

    68         degY = sel

    69         X = degT1

    70         degX = Fbcols

    71 

    72     else:

    73         print 'Message from left_prm:no legal job number specified'

    74         exit()

    75 

    76     return B, degB, A, degA, Y, degY, X, degX

    77 

    78 

    79 if __name__== "__main__":

    80     D = pl.array([[1, 0, 0, 0, 0, 0],

    81                   [0, 1, 0, 1, 0, 0],

    82                   [0, 0, 1, 1, 1, 0]])

    83 

    84     N = pl.array([[1, 0, 0],

    85                   [0, 1, 0],

    86                   [0, 0, 1]])

    87 

    88     dD, dN = 1, 0

    89 

    90     result = left_prm(N, dN, D, dD)

    91     print result

    92     pass