def METHOD_NAME(self):
    mesh = self.mesh
    pde = self.pde
    NE = mesh.number_of_edges()
    NC = mesh.number_of_cells()
    itype = mesh.itype
    ftype = mesh.ftype
    A = self.get_left_matrix()
    b = self.get_right_vector()
    pc = mesh.entity_barycenter('cell')
    x = np.r_[self.uh, self.ph]#把self.uh,self.ph组合在一起
    b = b - A@x
    # Modify matrix
    bdIdx = np.zeros((A.shape[0],), dtype = itype)
    bdIdx[NE] = 1
    Tbd = spdiags(bdIdx, 0, A.shape[0], A.shape[1])
    T = spdiags(1-bdIdx, 0, A.shape[0], A.shape[1])
    AD = T@A@T + Tbd
    b[NE] = self.ph[0]
    A11 = AD[:NE,:NE]
    A11inv = inv(A11)
    A12 = AD[:NE,NE:NE+NC]
    A21 = AD[NE:NE+NC,:NE]
    # solve
    p1 = np.zeros((NC,),dtype = ftype)
    Anew = A21*A11inv*A12
    bnew = A21*A11inv*b[:NE] - b[NE:NE+NC]
    p1[1:NC] = spsolve(Anew[1:NC,1:NC],bnew[1:NC]) 
    p1[0] = pde.pressure(pc[0,:])
    u1 = A11inv*(b[:NE] - A12*p1)
    self.uh[:] = u1
    self.ph[:] = p1
    return x