code: def factor_aug(z, DPhival, G, A):
      r
      M, N = G.shape
      P, N = A.shape
      """Multiplier for inequality constraints"""
      l = z[N+P:N+P+M]
  
      """Slacks"""
      s = z[N+P+M:]
  
      """Sigma matrix"""
      SIG = diags(l/s, 0)
      # SIG = diags(l*s, 0)
  
      """Convert A"""
      if not issparse(A):
          A = csr_matrix(A)
  
      """Convert G"""
      if not issparse(G):
          G = csr_matrix(G)
  
      """Since we expect symmetric DPhival, we need to change A"""
      sign = np.zeros(N)
      sign[0:N//2] = 1.0
      sign[N//2:] = -1.0
      T = diags(sign, 0)
  
      A_new = A.dot(T)
      
      W = AugmentedSystem(DPhival, G, SIG, A_new)
      return W
docs: r"""Set up augmented system and return.
  
      Parameters
      ----------
      z : (N+P+M+M,) ndarray
          Current iterate, z = (x, nu, l, s)
      DPhival : LinearOperator
          Jacobian of the variational inequality mapping
      G : (M, N) ndarray or sparse matrix
          Inequality constraints
      A : (P, N) ndarray or sparse matrix
          Equality constraints  
  
      Returns
      -------
      J : LinearOperator
          Augmented system