In [1]:

    

from sympy import *


    

In [7]:

    

n_u, n_b, m_p, m_r = symbols("n_u n_b m_p m_r")
F = MatrixSymbol("F", n_u, n_u)
M = MatrixSymbol("M", n_b, n_b)
C = MatrixSymbol("C", n_b, n_u)
B = MatrixSymbol("B", m_p, n_u)
D = MatrixSymbol("D", m_r, n_b)


    

In [12]:

    

B.T.shape


    

    
Out[12]:





(n_u, m_p)

In [33]:

    

A = BlockMatrix([[F,B.T,C.T,ZeroMatrix(n_u,m_r)],
                 [B,ZeroMatrix(m_p,m_p),ZeroMatrix(m_p,n_b),ZeroMatrix(m_p,m_r)],
                 [-C,ZeroMatrix(n_b,m_p),M,D.T],
                 [ZeroMatrix(m_r,n_u),ZeroMatrix(m_r,m_p),D,ZeroMatrix(m_r,m_r)]])
A


    

    
Out[33]:





[     F, B', C',  0]
[     B,  0,  0,  0]
[(-1)*C,  0,  M, D']
[     0,  0,  D,  0]

In [60]:

    

n_u, n_b, m_p, m_r = symbols("n_u n_b m_p m_r")
F = Symbol("F")
M = Symbol("M")
C = Symbol("C")
Ct = Symbol("Ct")
B = Symbol("B")
Bt = Symbol("Bt")
D = Symbol("D")
Dt = Symbol("Dt")


    

In [79]:

    

A = Matrix([[F,Bt,Ct,0],[B,0,0,0],[-C,0,M,Dt],[0,0,D,0]])
A


    

    
Out[79]:





[ F, Bt, Ct,  0]
[ B,  0,  0,  0]
[-C,  0,  M, Dt]
[ 0,  0,  D,  0]

In [80]:

    

A.LUdecomposition()


    

    
Out[80]:





([   1,    0,   0, 0]
[ B/F,    1,   0, 0]
[-C/F, -C/B,   1, 0]
[   0,    0, D/M, 1],
 [F,      Bt,      Ct,       0]
[0, -B*Bt/F, -B*Ct/F,       0]
[0,       0,       M,      Dt]
[0,       0,       0, -D*Dt/M],
 [])